为什么引入随机斜率效应会扩大斜率的 SE?

机器算法验证 r 混合模式 lme4-nlme 随机效应模型
2022-03-30 14:21:54

我正在尝试分析 Year 对特定群体(我有 3 个群体)的变量 logInd 的影响。最简单的模型:

> fix1 = lm(logInd ~ 0 + Group + Year:Group, data = mydata)
> summary(fix1)

Call:
lm(formula = logInd ~ 0 + Group + Year:Group, data = mydata)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.5835 -0.3543 -0.0024  0.3944  4.7294 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
Group1       4.6395740  0.0466217  99.515  < 2e-16 ***
Group2       4.8094268  0.0534118  90.044  < 2e-16 ***
Group3       4.5607287  0.0561066  81.287  < 2e-16 ***
Group1:Year -0.0084165  0.0027144  -3.101  0.00195 ** 
Group2:Year  0.0032369  0.0031098   1.041  0.29802    
Group3:Year  0.0006081  0.0032666   0.186  0.85235    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.7926 on 2981 degrees of freedom
Multiple R-squared: 0.9717,     Adjusted R-squared: 0.9716 
F-statistic: 1.705e+04 on 6 and 2981 DF,  p-value: < 2.2e-16 

我们可以看到第 1 组显着下降,第 2 组和第 3 组增加但不显着。

显然个体应该是随机效应,所以我为每个个体引入随机截距效应:

> mix1a = lmer(logInd ~ 0 + Group + Year:Group + (1|Individual), data = mydata)
> summary(mix1a)
Linear mixed model fit by REML 
Formula: logInd ~ 0 + Group + Year:Group + (1 | Individual) 
   Data: mydata 
  AIC  BIC logLik deviance REMLdev
 4727 4775  -2356     4671    4711
Random effects:
 Groups     Name        Variance Std.Dev.
 Individual (Intercept) 0.39357  0.62735 
 Residual               0.24532  0.49530 
Number of obs: 2987, groups: Individual, 103

Fixed effects:
              Estimate Std. Error t value
Group1       4.6395740  0.1010868   45.90
Group2       4.8094268  0.1158095   41.53
Group3       4.5607287  0.1216522   37.49
Group1:Year -0.0084165  0.0016963   -4.96
Group2:Year  0.0032369  0.0019433    1.67
Group3:Year  0.0006081  0.0020414    0.30

Correlation of Fixed Effects:
            Group1 Group2 Group3 Grp1:Y Grp2:Y
Group2       0.000                            
Group3       0.000  0.000                     
Group1:Year -0.252  0.000  0.000              
Group2:Year  0.000 -0.252  0.000  0.000       
Group3:Year  0.000  0.000 -0.252  0.000  0.000

它具有预期的效果 - 斜率的 SE(系数 Group1-3:Year)现在较低,剩余 SE 也较低。

个体的斜率也不同,所以我还介绍了随机斜率效应:

> mix1c = lmer(logInd ~ 0 + Group + Year:Group + (1 + Year|Individual), data = mydata)
> summary(mix1c)
Linear mixed model fit by REML 
Formula: logInd ~ 0 + Group + Year:Group + (1 + Year | Individual) 
   Data: mydata 
  AIC  BIC logLik deviance REMLdev
 2941 3001  -1461     2885    2921
Random effects:
 Groups     Name        Variance  Std.Dev. Corr   
 Individual (Intercept) 0.1054790 0.324775        
            Year        0.0017447 0.041769 -0.246 
 Residual               0.1223920 0.349846        
Number of obs: 2987, groups: Individual, 103

Fixed effects:
              Estimate Std. Error t value
Group1       4.6395740  0.0541746   85.64
Group2       4.8094268  0.0620648   77.49
Group3       4.5607287  0.0651960   69.95
Group1:Year -0.0084165  0.0065557   -1.28
Group2:Year  0.0032369  0.0075105    0.43
Group3:Year  0.0006081  0.0078894    0.08

Correlation of Fixed Effects:
            Group1 Group2 Group3 Grp1:Y Grp2:Y
Group2       0.000                            
Group3       0.000  0.000                     
Group1:Year -0.285  0.000  0.000              
Group2:Year  0.000 -0.285  0.000  0.000       
Group3:Year  0.000  0.000 -0.285  0.000  0.000

但是现在,出乎意料的是,斜率的 SE(系数 Group1-3:Year)现在要高得多,甚至比完全没有随机效应的情况下还要高!

这怎么可能?我希望随机效应会“吃掉”无法解释的可变性并增加估计的“确定性”!

然而,残差 SE 的行为与预期一样——它低于随机截距模型。

如果需要,这是数据。

编辑

现在我意识到了惊人的事实。如果我分别对每个人进行线性回归,然后在所得斜率上运行 ANOVA,我得到的结果与随机斜率模型完全相同!你知道为什么吗?

indivSlope = c()
for (indiv in 1:103) {
    mod1 = lm(logInd ~ Year, data = mydata[mydata$Individual == indiv,])
    indivSlope[indiv] = coef(mod1)['Year']
}

indivGroup = unique(mydata[,c("Individual", "Group")])[,"Group"]


anova1 = lm(indivSlope ~ 0 + indivGroup)
summary(anova1)

Call:
lm(formula = indivSlope ~ 0 + indivGroup)

Residuals:
      Min        1Q    Median        3Q       Max 
-0.176288 -0.016502  0.004692  0.020316  0.153086 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)
indivGroup1 -0.0084165  0.0065555  -1.284    0.202
indivGroup2  0.0032369  0.0075103   0.431    0.667
indivGroup3  0.0006081  0.0078892   0.077    0.939

Residual standard error: 0.04248 on 100 degrees of freedom
Multiple R-squared: 0.01807,    Adjusted R-squared: -0.01139 
F-statistic: 0.6133 on 3 and 100 DF,  p-value: 0.6079 

如果需要,这是数据。

1个回答

我认为问题出在您的期望上:) 请注意,当您为每个人添加随机截距时,截距的标准误差会增加由于每个人都可以有自己的截距,因此组平均值不太确定。随机斜率也发生了同样的事情:您不再估计一个常见的(组内)斜率,而是估计不同斜率的平均值。

编辑:为什么没有更好的模型给出更精确的估计?

让我们反过来想一想:为什么初始模型低估了标准误差?它假设不独立的观察是独立的。第二个模型放宽了该假设(以影响截距的方式),第三个模型进一步放宽了它。

编辑 2:与许多患者特定模型的关系

您的观察是一个已知属性(如果您只有两年,那么随机效应模型将等同于配对 t 检验)。我不认为我可以管理一个真正的证明,但也许写出这两个模型会使关系更清晰。让我们忽略分组变量,因为它只会使符号复杂化。我将使用希腊字母表示随机效果,使用拉丁字母表示固定效果。

随机效应模型为 ( - 主体, - 在主体内复制): 其中ij

Yij=a+αi+(b+βi)xij+ϵij,
(αi,βi)N(0,Σ)ϵijN(0,σ2)

当您为每个主题拟合单独的模型时,则 其中

Yij=ai+bixij+ϵij,
ϵijN(0,σi2)

[注:以下真的只是挥手而已:]

您可以看到这两个模型之间有很多相似之处,对应于对应于的平均值对应于,因为随机效应平均为 0。随机截距和斜率的无约束相关性导致模型可以单独拟合。我不确定单个假设如何与特定于主题的啮合,但我会假设会有所不同。aia+αibib+βibibσσiαi