如何修改lme4中的槽>1.0

3

我一直使用下面的代码成功地修改了使用 lme4 版本 <1.0 的模型中的 'Zt'、'L' 和 'A' 插槽。今天我刚更新到了 lme4 1.0-4,发现模型对象不同了。有人能提供关于如何在新的 lmer 模型对象中修改这些插槽的见解/指导吗?

dat<-structure(list(pop1 = c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 10L, 
                    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 
                    4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 7L, 7L, 
                    7L, 8L, 8L, 9L), pop2 = c(2L, 3L, 4L, 5L, 6L, 7L, 8L, 9L, 1L, 
                    3L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 4L, 5L, 6L, 7L, 8L, 9L, 10L, 
                    5L, 6L, 7L, 8L, 9L, 10L, 6L, 7L, 8L, 9L, 10L, 7L, 8L, 9L, 10L, 
                    8L, 9L, 10L, 9L, 10L, 10L), X = c(0.49136, 0.75587, 0.93952, 
                    0.61278, 0.79934, 1.07918, 1.13354, 1.15836, 1.2014, 0.43136, 
                    0.77815, 0.716, 0.93952, 1.13672, 1.16137, 1.18184, 1.21748, 
                    0.65321, 0.86332, 1.04922, 1.19866, 1.20412, 1.22272, 1.24797, 
                    0.89763, 1.08991, 1.19033, 1.15836, 1.17319, 1.18752, 0.64345, 
                    0.93952, 0.98227, 1.01703, 1.07188, 0.78533, 0.94939, 0.99564, 
                    1.06819, 0.64345, 0.716, 0.85126, -0.04576, 0.4624, 0.30103), 
           Y = c(0.491694, 0.394703, 0.113303, 0.156597, 0.450924, 0.487845, 
                 0.21821, 0.129027, -0.131522, 0.35156, -0.116826, 0.18941, 
                 0.306608, 0.258401, 0.008552, -0.024369, -0.305258, -0.013628, 
                 0.215715, 0.13783, 0.467272, 0.088882, 0.084295, -0.172337, 
                 -0.206725, -0.084339, -0.191651, -0.001586, -0.079501, -0.195094, 
                 0.232045, 0.17102, 0.003742, -0.023688, -0.26085, 0.205326, 
                 0.172809, 0.133219, -0.159054, 0.082231, 0.011025, -0.238611, 
                 0.732679, 0.478058, 0.325698)), .Names = c("pop1", "pop2", 
                  "X", "Y"), class = "data.frame", row.names = c(NA, -45L))

library(lme4) # lme4 versions >1.0 have different model output    

# Specify the model formula 
    lmer_mod <- as.formula("Y ~ X + (1|pop1)") 

# Create the Zl and ZZ matrices 
    Zl <- lapply(c("pop1","pop2"), function(nm) Matrix:::fac2sparse(dat[[nm]], "d", drop=FALSE)) 
    ZZ <- Reduce("+", Zl[-1], Zl[[1]]) 

# Fit lmer model to the data 
    mod <- lmer(lmer_mod, data = dat, REML = TRUE) 

# Replace the following slots in the fitted model
# These slots don't exist in this form in the new lmerMod objects
    mod@Zt <- ZZ 
    mod@A <- ZZ 
    mod@L <- Cholesky(tcrossprod(ZZ), LDL=FALSE, Imult=1) 

# Refit the model to the same response data 
    Final.mod <- refit(mod, dat[,Y])

任何关于如何修改这些插槽的帮助或见解将不胜感激。同时,我猜我会继续使用旧版本的lme4来建立这些模型。

请先阅读新的lme4包中的?modular页面以开始... - Ben Bolker
1个回答

2

这个是否满足您的需求?(这个跟?modular很相似…)

创建Zl和ZZ矩阵:

Zl <- lapply(c("pop1","pop2"),
         function(nm) Matrix:::fac2sparse(dat[[nm]], "d", drop=FALSE)) 
ZZ <- Reduce("+", Zl[-1], Zl[[1]]) 

构建随机效应数据结构:
lf <- lFormula(Y ~ X + (1|pop1), data=dat)

修改它们:

lf$reTrms$Zt <- ZZ

继续完成剩余的模型构建和拟合步骤:

dfun <- do.call(mkLmerDevfun,lf)   ## create dev fun from modified lf
opt <- optimizeLmer(dfun)          ## fit the model
## make the results into a 'merMod' object
fit <- mkMerMod(environment(dfun), opt, lf$reTrms,
         fr = lf$fr)

这看起来很有潜力,但是偏差函数产生了这个错误:Error in Lambdat %*% Ut: Cholmod error 'A and B inner dimensions must match' at file ../MatrixOps/cholmod_ssmult.c, line 82。 - Bill Peterman
很有趣,对我有效。你能给出sessionInfo()吗?这是使用lme4 1.1-0、Matrix 1.0-14进行的...你可能需要确保最近重新安装了Matrix和RcppEigen——编译基础库必须匹配... - Ben Bolker
抱歉,Ben。你的解决方案有效。我在原始代码中没有正确指定我的Zl矩阵(应该是“dat”,而不是“data”)。有了正确的Zl和ZZ矩阵,你的代码可以工作,并且结果与预期输出相匹配。谢谢!(我已经编辑了我的原始代码,以便可以重现这个问题) - Bill Peterman

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接