Вложенные случайные эффекты в `lme {nlme}`

Question

Мне нужно подогнать LMM со случайным эффектом взаимодействия, но без маргинального случайного эффекта, используя команду lme. То есть, я хочу подгонять модель в oats.lmer (см. Ниже), но используя функцию lme из пакета nlme.

Код

require("nlme") 
require("lme4") 

oats.lmer <- lmer(yield~nitro + (1|Block:Variety), data = Oats) 

summary(oats.lmer) 
#Linear mixed model fit by REML ['lmerMod'] 
#Formula: yield ~ nitro + (1 | Block:Variety) 
# Data: Oats 
# 
#REML criterion at convergence: 598.1 
# 
#Scaled residuals: 
#  Min  1Q Median  3Q  Max 
#-1.66482 -0.72807 -0.00079 0.56416 1.85467 
# 
#Random effects: 
# Groups  Name  Variance Std.Dev. 
# Block:Variety (Intercept) 306.8 17.51 
# Residual     165.6 12.87 
#Number of obs: 72, groups: Block:Variety, 18 
# 
#Fixed effects: 
#   Estimate Std. Error t value 
#(Intercept) 81.872  4.846 16.90 
#nitro   73.667  6.781 10.86 
# 
#Correlation of Fixed Effects: 
#  (Intr) 
#nitro -0.420 

Я начал играть с этим

oats.lme <- lme(yield~nitro, data = Oats, random = (~1|Block/Variety)) 
summary(oats.lme) 
#Linear mixed-effects model fit by REML 
# Data: Oats 
#  AIC  BIC logLik 
# 603.0418 614.2842 -296.5209 
# 
#Random effects: 
# Formula: ~1 | Block 
#  (Intercept) 
#StdDev: 14.50596 
# 
# Formula: ~1 | Variety %in% Block 
#  (Intercept) Residual 
#StdDev: 11.00468 12.86696 
# 
#Fixed effects: yield ~ nitro 
#    Value Std.Error DF t-value p-value 
#(Intercept) 81.87222 6.945273 53 11.78819  0 
#nitro  73.66667 6.781483 53 10.86291  0 
# Correlation: 
#  (Intr) 
#nitro -0.293 
# 
#Standardized Within-Group Residuals: 
#  Min   Q1   Med   Q3   Max 
#-1.74380770 -0.66475227 0.01710423 0.54298809 1.80298890 
# 
#Number of Observations: 72 
#Number of Groups: 
#    Block Variety %in% Block 
#     6     18 

, но проблема в том, что она ставит также маргинальное случайный эффект для Variety, который я хочу, чтобы опустить.

Вопрос заключается в следующем: как определить случайные эффекты в oats.lme так, чтобы oats.lme был идентичным (не менее структурным) до oats.lmer?

Answer 1

Это может быть простым, как следующие:

library(nlme) 
data(Oats) 

## construct an auxiliary factor `f` for interaction/nesting effect 
Oats$f <- with(Oats, Block:Variety) 
## use `random = ~ 1 | f` 
lme(yield ~ nitro, data = Oats, random = ~ 1 | f) 

#Linear mixed-effects model fit by REML 
# Data: Oats 
# Log-restricted-likelihood: -299.0328 
# Fixed: yield ~ nitro 
#(Intercept)  nitro 
# 81.87222 73.66667 
# 
#Random effects: 
# Formula: ~1 | f 
#  (Intercept) Residual 
#StdDev: 17.51489 12.86695 
# 
#Number of Observations: 72 
#Number of Groups: 18