缺少R级交互作用的lme4异常/处理对比编码

Question

我有一个混合效应模型(与lme4)，有一个双向相互作用项，每个术语有多个层次(每个4个)，我想参照它们的大平均值来研究它们的影响。我在这里从car数据集中介绍了这个示例，并省略了错误项，因为这个示例不需要它：

## shorten data frame for simplicity
df=Cars93[c(1:15),]
df=Cars93[is.element(Cars93$Make,c('Acura Integra', 'Audi 90','BMW 535i','Subaru Legacy')),]
df$Make=drop.levels(df$Make)
df$Model=drop.levels(df$Model)

## define contrasts (every factor has 4 levels)
contrasts(df$Make) = contr.treatment(4)
contrasts(df$Model) = contr.treatment(4)

## model
m1 <- lm(Price ~ Model*Make,data=df)
summary(m1)

正如您所看到的，在交互项中省略了第一个级别。我希望输出中的所有4个级别，引用大平均值(通常指异常编码)。这些是我看过的来源：https://marissabarlaz.github.io/portfolio/contrastcoding/#coding-schemes和How to change contrasts to compare with mean of all levels rather than reference level (R, lmer)?。不过，最后一个引用并不报告交互。

Answer 1

简单的答案是你想要的东西是不可能直接实现的。你必须使用一种稍微不同的方法。

在具有交互作用的模型中，您希望使用均值为零而不是特定级别的对比。否则，低阶效应(即主效应)不是主要效应，而是简单效应(当其他因素水平处于参考水平时评估)。这在我关于混合模型的章节中有更详细的解释：http://singmann.org/download/publications/singmann_kellen-introduction-mixed-models.pdf

为了得到你想要的，你必须以一种合理的方式来拟合这个模型，然后把它传递给emmeans来与截距(即未加权的大平均值)进行比较。这也适用于交互，如下所示(由于您的代码不起作用，我使用warpbreaks)。

afex::set_sum_contrasts() ## uses contr.sum globally
library("emmeans")

## model
m1 <- lm(breaks ~ wool * tension,data=warpbreaks)
car::Anova(m1, type = 3)

coef(m1)[1]
# (Intercept) 
#    28.14815 

## both CIs include grand mean:
emmeans(m1, "wool") 
#  wool emmean   SE df lower.CL upper.CL
#  A      31.0 2.11 48     26.8     35.3
#  B      25.3 2.11 48     21.0     29.5
# 
# Results are averaged over the levels of: tension 
# Confidence level used: 0.95 

## same using test
emmeans(m1, "wool", null = coef(m1)[1], infer = TRUE) 
#   wool emmean   SE df lower.CL upper.CL null t.ratio p.value
#  A      31.0 2.11 48     26.8     35.3 28.1  1.372  0.1764 
#  B      25.3 2.11 48     21.0     29.5 28.1 -1.372  0.1764 
# 
# Results are averaged over the levels of: tension 
# Confidence level used: 0.95 

emmeans(m1, "tension", null = coef(m1)[1], infer = TRUE) 
#  tension emmean   SE df lower.CL upper.CL null t.ratio p.value
#  L         36.4 2.58 48     31.2     41.6 28.1  3.196  0.0025 
#  M         26.4 2.58 48     21.2     31.6 28.1 -0.682  0.4984 
#  H         21.7 2.58 48     16.5     26.9 28.1 -2.514  0.0154 
# 
# Results are averaged over the levels of: wool 
# Confidence level used: 0.95 

emmeans(m1, c("tension", "wool"), null = coef(m1)[1], infer = TRUE) 
#  tension wool emmean   SE df lower.CL upper.CL null t.ratio p.value
#  L       A      44.6 3.65 48     37.2     51.9 28.1  4.499  <.0001 
#  M       A      24.0 3.65 48     16.7     31.3 28.1 -1.137  0.2610 
#  H       A      24.6 3.65 48     17.2     31.9 28.1 -0.985  0.3295 
#  L       B      28.2 3.65 48     20.9     35.6 28.1  0.020  0.9839 
#  M       B      28.8 3.65 48     21.4     36.1 28.1  0.173  0.8636 
#  H       B      18.8 3.65 48     11.4     26.1 28.1 -2.570  0.0133 
# 
# Confidence level used: 0.95 

请注意，对于coef()，您可能希望对lme4模型使用fixef()。