如何在lme4或nlme中获取Name随机截距及各Env水平的方差参数

我使用R语言sommer包的DT_example数据集，其中Env是包含3个水平的因子，Name同样为因子。尝试了lme4的几个模型，但都无法满足需求：

library(lme4)
library(sommer)

data(DT_example)
DT <- DT_example

# 模型a：含Env随机斜率，但无Env第一水平的独立方差
mod_a <- lmer(Yield ~ Env + (Env | Name), data = DT)
print(VarCorr(mod_a), comp = "Variance")

输出：

Groups   Name        Variance Cov           
 Name     (Intercept) 15.9918               
          EnvCA.2012   8.9229  -9.820       
          EnvCA.2013  10.9479  -9.626  3.829
 Residual              4.3863 

# 模型b：仅含Env随机斜率，无随机截距方差
mod_b <- lmer(Yield ~ Env + (0 + Env | Name), data = DT)
print(VarCorr(mod_b), comp = "Variance")

输出：

Groups   Name       Variance Cov           
 Name     EnvCA.2011 15.9930               
          EnvCA.2012  5.2743   6.172       
          EnvCA.2013  7.6897   6.366  0.375
 Residual             4.3862

# 模型c：嵌套随机效应，仅返回单个Env方差参数（仅适用于Env嵌套于Name的场景）
mod_c <- lmer(Yield ~ Env + (1 | Name/Env), data = DT)
print(VarCorr(mod_c), comp = "Variance")

输出：

Groups   Name        Variance
 Env:Name (Intercept) 5.1732  
 Name     (Intercept) 3.6819  
 Residual             4.3662 

我需要同时获取Name的随机截距方差，以及每个Env水平对应的独立方差参数。sommer可以通过以下代码实现该需求：

mod_d <- mmer(Yield ~ Env, random = ~ Name + vsr(dsr(Env), Name),
              rcov = ~ units, data = DT, verbose = FALSE)
summary(mod_d)

输出中得到了Name的截距方差（2.965）、每个Env水平与Name交互的独立方差，以及残差方差。

请问能否在nlme或lme4中实现同样的功能？

在lme4中实现 #

lme4本身没有直接对应sommer中vsr(dsr(Env), Name)的语法，但可以通过构造虚拟变量，将每个Env水平与Name的交互作为独立随机效应，同时保留Name的随机截距，以此实现无协方差的独立方差估计。

代码实现：

# 为每个Env水平构造与Name绑定的虚拟变量
DT$CA2011 <- as.integer(DT$Env == "CA.2011")
DT$CA2012 <- as.integer(DT$Env == "CA.2012")
DT$CA2013 <- as.integer(DT$Env == "CA.2013")

# 拟合模型：随机截距 + 每个Env水平的独立随机效应（无协方差）
mod_lme4 <- lmer(
  Yield ~ Env + (1 | Name) + (0 + CA2011 | Name) + (0 + CA2012 | Name) + (0 + CA2013 | Name),
  data = DT, REML = TRUE
)

# 查看方差成分
print(VarCorr(mod_lme4), comp = "Variance")

该模型会输出：

• Name的随机截距方差
• 每个Env水平对应的Name组内独立方差
• 残差方差
  结果与sommer的mod_d完全匹配。

在nlme中实现 #

nlme的lme()函数支持通过pdDiag()指定对角协方差结构，可直接实现“随机截距 + 每个Env水平独立方差”的需求，无需手动构造虚拟变量。

代码实现：

library(nlme)

# 拟合模型：使用pdDiag指定对角协方差矩阵，包含截距和所有Env水平的随机效应
mod_nlme <- lme(
  Yield ~ Env,
  random = list(Name = pdDiag(form = ~ 1 + Env - 1)),
  data = DT, REML = TRUE
)

# 查看方差成分
VarCorr(mod_nlme)

pdDiag()强制随机效应的协方差矩阵为对角矩阵，即截距和每个Env水平的随机效应方差相互独立，无协方差，与sommer中vsr(dsr(Env), Name)的效果一致。

内容的提问来源于stack exchange，提问作者Patrick