Несколько случайных эффектов без взаимодействия в nlme :: lme - PullRequest
Новая
       9

Несколько случайных эффектов без взаимодействия в nlme :: lme

0 голосов
/ 07 сентября 2018

Рассмотрим эту простую модель, состоящую из двух независимых случайных эффектов: 

$$y_{ijk} = \mu + \delta1_i + \delta2_j + \epsilon_{i,j,k}$$

, где \delta1_i и \delta2_j - независимые случайные эффекты (т. Е. \delta1_i \sim N(0,\sigma_1^2) iid, \delta2_i \sim N(0,\sigma_1^2) iid, оба независимые).

Определить и вывести компоненты дисперсии такой модели просто, используя lme4 :: lmer (извините, приведенный ниже пример, скорее всего, не имеет физического значения): 

library(lme4)
library(nlme)
data(Orthodont,package="nlme")   

lmer(data=Orthodont,distance~1+(1|Sex)+(1|Subject))
#Linear mixed model fit by REML ['lmerMod']
#Formula: distance ~ 1 + (1 | Sex) + (1 | Subject)
#   Data: Orthodont
#REML criterion at convergence: 510.3937
#Random effects:
# Groups   Name        Std.Dev.
# Subject  (Intercept) 1.596   
# Sex      (Intercept) 1.550   
# Residual             2.220   
#Number of obs: 108, groups:  Subject, 27; Sex, 2
#Fixed Effects:
#(Intercept)  
#      23.83

Как можноэто легко смоделировать с помощью nlme :: lme?Указание случайного эффекта с использованием списка случайных эффектов неявно вводит взаимодействие в зависимости от порядка случайных эффектов, что можно проверить здесь (изменение порядка спецификации, изменение результатов): 

VarCorr(lme(data=Orthodont,distance~1,random=list(~1|Subject,~1|Sex)))
#            Variance     StdDev  
#Subject =   pdLogChol(1)         
#(Intercept) 1.875987     1.369667
#Sex =       pdLogChol(1)         
#(Intercept) 1.875987     1.369667
#Residual    4.929784     2.220312

VarCorr(lme(data=Orthodont,distance~1,random=list(~1|Sex,~1|Subject)))
#            Variance     StdDev  
#Sex =       pdLogChol(1)         
#(Intercept) 2.403699     1.550387
#Subject =   pdLogChol(1)         
#(Intercept) 2.546702     1.595839
#Residual    4.929784     2.220312
...