データの構造と私が作成した結果をできるだけ多くお見せしたいと思います。
データの構造は次のとおりです。
GroupID Person Factor2 Factor1 Rating
<int> <int> <fctr> <fctr> <int>
1 2 109 2 0 1
2 2 109 2 1 -2
3 2 104 1 0 4
4 2 236 1 1 1
5 2 279 1 1 2
6 2 179 2 1 0
Person は参加者 ID、GroupID は評価された刺激の種類、因子 1 (レベル 0 と 1) と因子 2 (レベル 1 と 2) は固定因子、評価は結果変数です。
線形混合効果モデルにおける重要な交互作用のプロットを印刷しようとしています。パッケージ lme4 と lmerTest を使用してデータを分析しました。
これが実行したモデルです。
> model_interaction <- lmer(Rating ~ Factor1 * Factor2 + ( 1 | Person) +
(1 | GroupID), data)
> model_interaction
Linear mixed model fit by REML ['merModLmerTest']
Formula: Rating ~ Factor1 * Factor2 + (1 | Person) + (1 | GroupID)
Data: data
REML criterion at convergence: 207223.9
Random effects:
Groups Name Std.Dev.
Person (Intercept) 1.036
GroupID (Intercept) 1.786
Residual 1.880
Number of obs: 50240, groups: Person, 157; GroupID, 80
Fixed Effects:
(Intercept) Factor11 Factor22 Factor11:Factor22
-0.43823 0.01313 0.08568 0.12440
summary() 関数を使用すると、R は次の出力を返します。
> summary(model_interaction)
Linear mixed model fit by REML
t-tests use Satterthwaite approximations to degrees of freedom
['lmerMod']
Formula: Rating ~ Factor1 * Factor2 + (1 | Person) + (1 | GroupID)
Data: data
REML criterion at convergence: 207223.9
Scaled residuals:
Min 1Q Median 3Q Max
-4.8476 -0.6546 -0.0213 0.6516 4.2284
Random effects:
Groups Name Variance Std.Dev.
Person (Intercept) 1.074 1.036
GroupID (Intercept) 3.191 1.786
Residual 3.533 1.880
Number of obs: 50240, groups: Person, 157; GroupID, 80
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) -4.382e-01 2.185e-01 1.110e+02 -2.006 0.047336 *
Factor11 1.313e-02 2.332e-02 5.004e+04 0.563 0.573419
Factor22 8.568e-02 6.275e-02 9.793e+03 1.365 0.172138
Factor11:Factor22 1.244e-01 3.385e-02 5.002e+04 3.675 0.000238 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) Fctr11 Fctr22
Factor11 -0.047
Factor22 -0.135 0.141
Fctr11:Fc22 0.034 -0.694 -0.249
線形混合効果モデルの p 値を解釈できないことはわかっています。そこで、因子 1 と因子 2 の主効果のみを含むモデルと相互作用モデルを比較して、追加の anova を実行しました。
> model_Factor1_Factor2 = lmer(Rating ~ Factor1 + Factor2 +
( 1 | Person) + (1 | GroupID), data)
> anova(model_Factor1_Factor2, model_interaction)
Data: data
Models:
object: Rating ~ Factor1 + Factor2 + (1 | Person) + (1 | GroupID)
..1: Rating ~ Factor1 * Factor2 + (1 | Person) + (1 | GroupID)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
object 6 207233 207286 -103611 207221
..1 7 207222 207283 -103604 207208 13.502 1 0.0002384 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
私はこのアウトプットを次のように解釈しました。Factor1 と Factor2 の交互作用は、Factor1 と Factor2 の主効果のみを使用したモデルと比較して、結果測定値の追加の分散を説明しています。
線形混合効果モデルの出力を解釈するのは難しいので、Factor1 と Factor2 の相互作用を示すグラフを印刷したいと思います。私はlsmeansパッケージを使用してそうしました（最初にplot（allEffects）を使用しましたが、混合効果モデルで係数とその信頼区間を取得する方法の質問を読んだ後、線形混合効果モデルのグラフを印刷する正しい方法ではないことに気付きました）。
これが私がやったことです（このWebサイトhttp://rcompanion.org/handbook/G_06.htmlに従って）
> leastsquare = lsmeans(model_interaction, pairwise ~ Factor2:Factor1,
adjust="bon")
> CLD = cld(leastsquare, alpha=0.05, Letters=letters, adjust="bon")
> CLD$.group=gsub(" ", "", CLD$.group)
> CLD
Factor2 Factor1 lsmean SE df lower.CL upper.CL .group
1 0 -0.4382331 0.2185106 111.05 -0.9930408 0.1165746 a
1 1 -0.4251015 0.2186664 111.36 -0.9803048 0.1301018 a
2 0 -0.3525561 0.2190264 112.09 -0.9086735 0.2035612 a
2 1 -0.2150234 0.2189592 111.95 -0.7709700 0.3409233 b
Degrees-of-freedom method: satterthwaite
Confidence level used: 0.95
Conf-level adjustment: bonferroni method for 4 estimates
P value adjustment: bonferroni method for 6 tests
significance level used: alpha = 0.05
これは私が使用したプロット機能です
> ggplot(CLD, aes(`Factor1`, y = lsmean, ymax = upper.CL,
ymin = lower.CL, colour = `Factor2`, group = `Factor2`)) +
geom_pointrange(stat = "identity",
position = position_dodge(width = 0.1)) +
geom_line(position = position_dodge(width = 0.1))
プロットは、このリンクを使用して見つけることができます (まだ画像を投稿することは許可されていません。回避策をお許しください)。
ここで私の質問は次のとおりです。有意な相互作用と、この相互作用によって説明される分散がかなりあるのに、プロットの信頼区間が重なっているのはなぜですか? 信頼区間で何か間違ったことをしたと思いますか？それとも、線形混合効果モデルの有意指数を解釈できないからでしょうか?