将SAS中的混合模型公式转换为R语言

Question

将SAS中的混合模型公式转换为R语言

5

我希望使用R中的nlme包拟合混合模型，该模型等价于以下SAS代码：

proc mixed data = one;
class var1 var2  year loc rep;
model yld = var1 * var2;
random loc year(loc) rep*year(loc);

编辑：实验内容解释

同样的 var1 和 var2 组合进行了多次重复测试（rep- 重复测试编号为 1:3）。重复测试（rep）是随机的。这组实验在不同地点（loc）和年份（year）中重复。虽然为方便起见，在每个地点和年份中，重复测试都被编号为 1:3，但一个地点和一年内的重复测试 1 与其他地点和年份内的重复测试 1 没有相关性。

我尝试了以下代码：

 require(nlme) 
    fm1 <- lme(yld ~ var1*var2, data = one, random = loc + year / loc + rep * year / loc)

我的代码是否正确？

编辑：根据建议修改了数据和模型。您可以从以下链接下载示例数据文件： https://sites.google.com/site/johndatastuff/mydata1.csv

data$var1 <- as.factor(data$var1)
data$var2 <- as.factor(data$var2)
data$year <- as.factor(data$year)
data$loc <- as.factor(data$loc)
data$rep <- as.factor(data$rep)

following suggestions from the comments below:
fm1 <- lme(yld ~ var1*var2, data = data, random = ~ loc + year / loc + rep * year / loc)

Error in getGroups.data.frame(dataMix, groups) : 
  Invalid formula for groups

根据SAS输出预期

Type 3 tests of fixed effects 
var1*var2         14         238       F value 16.12 Pr >F = < 0.0001

Covariance parameters:
loc = 0, year(loc) = 922161, year*rep(loc) = 2077492, residual = 1109238

我试过下面这个模型，但还是遇到了一些错误：

Edits: Just for information I tried the following model
require(lme4)  
 fm1 <- lmer(yld ~ var1*var2 + (1|loc) +  (1|year / loc) + (1|rep : (year / loc)),  
            data = data)  
Error in rep:`:` : NA/NaN argument 
In addition: Warning message: 
In rep:`:` : numerical expression has 270 elements: only the first used

- John

2

当您尝试它时发生了什么？您能否提供一个小的可重现示例？ - Aaron left Stack Overflow

你尝试在公式中加入 ~ 吗？fm1 <- lme(yld ~ var1*var2, data = one, random = ~ loc + year / loc + rep * year / loc) - Max Gordon

嗯...我没有使用过lme/lmer，但是任何比lmer(yld ~ var1*var2 +(1|year), data=data)更高级的东西都会失败。 - Max Gordon

2

@John，请不要在SE网站之间重复发布帖子。作为一个经验法则，纯R语言问题应该在SO上发布，统计学问题则应该在这里发布；否则，请选择一个网站并标记您的问题以引起管理员的注意。我将把这个问题移动到SO上，并与您的其他问题合并，以保留评论历史记录。（感谢Aaron的注意！） - chl

非常抱歉，我试图转移问题，因为我无法得到解决方案...但是不知道...然后只是将其平移...谢谢您处理这种情况。 - jon

显示剩余3条评论

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Aaron left Stack Overflow · Accepted Answer

感谢提供更详细的信息。我将数据存储在d中，以避免与data函数和参数混淆；命令可以使用任何一种方式，但是避免使用data通常被认为是一个好习惯。

请注意，交互很难适应，因为var和var2之间缺乏平衡；以下是交叉表的参考：

> xtabs(~var1 + var2, data=d)
    var2
var1  1  2  3  4  5
   1 18 18 18 18 18
   2  0 18 18 18 18
   3  0  0 18 18 18
   4  0  0  0 18 18
   5  0  0  0  0 18

通常为了仅适应交互作用（而不包括主效应），需要使用 : 而不是 *，但在这种情况下，最好将其设置为单个因素，如下所示：

d$var12 <- factor(paste(d$var1, d$var2, sep=""))

然后使用nlme，尝试：

fm1 <- lme(yld ~ var12, random = ~ 1 | loc/year/rep, data = d)
anova(fm1)

使用lme4，尝试以下操作：

fm1 <- lmer(yld ~ var12 + (1 | loc/year/rep), data=d)
anova(fm1)

注意，因为 nlme 和 lme4 的函数名称有重叠，所以在 R 会话中一次只能加载一个；如果需要切换，需要关闭 R 并重新启动。（还有其他方法，但这是最简单的解释方式。）