lmer() 和 lme() 的结果完全不同

机器算法验证 r 混合模式 lme4-nlme
2022-04-11 18:12:24

我从 lmer() 和 lme() 得到了完全不同的结果!只需查看系数的 std.errors。两种情况完全不同。为什么会这样,哪个模型是正确的?

> mix1c = lmer(logInd ~ 0 + crit_i + Year:crit_i + (1 + Year|Taxon), data = datai)
> mix1d = lme(logInd ~ 0 + crit_i + Year:crit_i, random = ~ 1 + Year|Taxon, data = datai)
> 
> summary(mix1d)
Linear mixed-effects model fit by REML
 Data: datai 
       AIC      BIC    logLik
  4727.606 4799.598 -2351.803

Random effects:
 Formula: ~1 + Year | Taxon
 Structure: General positive-definite, Log-Cholesky parametrization
            StdDev       Corr  
(Intercept) 9.829727e-08 (Intr)
Year        3.248182e-04 0.619 
Residual    4.933979e-01       

Fixed effects: logInd ~ 0 + crit_i + Year:crit_i 
                 Value Std.Error   DF   t-value p-value
crit_iA      29.053940  4.660176   99  6.234515  0.0000
crit_iF       0.184840  3.188341   99  0.057974  0.9539
crit_iU      12.340580  5.464541   99  2.258301  0.0261
crit_iW       5.324854  5.152019   99  1.033547  0.3039
crit_iA:Year -0.012272  0.002336 2881 -5.253846  0.0000
crit_iF:Year  0.002237  0.001598 2881  1.399542  0.1618
crit_iU:Year -0.003870  0.002739 2881 -1.412988  0.1578
crit_iW:Year -0.000305  0.002582 2881 -0.118278  0.9059
 Correlation: 
             crit_A crit_F crit_U crit_W cr_A:Y cr_F:Y cr_U:Y
crit_iF       0                                              
crit_iU       0      0                                       
crit_iW       0      0      0                                
crit_iA:Year -1      0      0      0                         
crit_iF:Year  0     -1      0      0      0                  
crit_iU:Year  0      0     -1      0      0      0           
crit_iW:Year  0      0      0     -1      0      0      0    

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-6.98370498 -0.39653580  0.02349353  0.43356564  5.15742550 

Number of Observations: 2987
Number of Groups: 103 
> summary(mix1c)
Linear mixed model fit by REML 
Formula: logInd ~ 0 + crit_i + Year:crit_i + (1 + Year | Taxon) 
   Data: datai 
  AIC  BIC logLik deviance REMLdev
 2961 3033  -1469     2893    2937
Random effects:
 Groups   Name        Variance   Std.Dev. Corr   
 Taxon    (Intercept) 6.9112e+03 83.13360        
          Year        1.7582e-03  0.04193 -1.000 
 Residual             1.2239e-01  0.34985        
Number of obs: 2987, groups: Taxon, 103

Fixed effects:
               Estimate Std. Error t value
crit_iA      29.0539403 18.0295239   1.611
crit_iF       0.1848404 12.3352135   0.015
crit_iU      12.3405800 21.1414908   0.584
crit_iW       5.3248537 19.9323887   0.267
crit_iA:Year -0.0122717  0.0090916  -1.350
crit_iF:Year  0.0022365  0.0062202   0.360
crit_iU:Year -0.0038701  0.0106608  -0.363
crit_iW:Year -0.0003054  0.0100511  -0.030

Correlation of Fixed Effects:
            crit_A crit_F crit_U crit_W cr_A:Y cr_F:Y cr_U:Y
crit_iF      0.000                                          
crit_iU      0.000  0.000                                   
crit_iW      0.000  0.000  0.000                            
crit_iA:Yer -1.000  0.000  0.000  0.000                     
crit_iF:Yer  0.000 -1.000  0.000  0.000  0.000              
crit_iU:Yer  0.000  0.000 -1.000  0.000  0.000  0.000       
crit_iW:Yer  0.000  0.000  0.000 -1.000  0.000  0.000  0.000
> 
2个回答

lmer使用拉普拉斯近似,当随机效应的整个正态分布以其模式近似时。众所周知,这种近似会产生下偏的方差分量的估计值。lme通过高斯正交近似使用更彻底的近似,但我既不知道积分点的默认数量也不知道操纵这个数字的方法。

另请注意,这lmer产生了 -1 的随机效应(截距与年份)的相关性。这非常糟糕,并且表明存在数字问题。它在它无法克服的可能性的近似值中碰到了某种山脊(难怪,鉴于这种近似值很差)。lme做得更好一些,并得出了 0.619 的相关性。如果你真的有一年,比如 2010 年、2011 年、2012 年,这对于混合模型来说是一个非常糟糕的主意,因为因子之间的高度多重共线性可能意味着数值稳定性差和收敛时间长。time如果将其转换为接近零的中心,而不是接近 2010 年,你会好得多。

结果的对数似然显着更高lmer(),这表明它更接近真正的最优值。

由于lmer()拟合是在参数空间的边界上,它lme()重新参数化到无穷大(“Log-Cholesky”参数化),因此lme没有找到更好的解决方案是合理的。

这并不能完全回答哪个答案更有用的问题,但它确实反对那个答案lme