与范围相关的（？）：创建的混合效果模型列表中的anova（）

debugcn 发表于 Dev

特拉班特耶

在使用进行混合效果建模的项目中lme，我试图比较具有不同相关结构和相等固定部分的模型。由于我将构建许多模型（针对不同的因变量），因此我尝试编写一个函数以生成具有不同相关结构的模型列表，如下面的示例所示（我实际上试图将其保持在最低限度）工作示例）。

如果我运行此列表中的元素，这个作品方差分析（），但只当fixedPart是我的全球环境。为什么会这样呢？有没有办法解决这个问题，以便我可以保留m并重新使用/删除fixedPart？

我认为这个问题与R中的（词法）作用域有关，但是我找不到真正解决它的方法。

提前致谢！

#Dependencies
library(multilevel)
library(multcomp)

#Generate sample data
nVals = 100
sData = rnorm(nVals, mean = 1, sd = 1)
dF <- data.frame(nSubject = 1:nVals, 
                 v1data = sData + rnorm(nVals, mean = 0, sd = 0.1),
                 v2data = sData + rnorm(nVals, mean = 0, sd = 0.1),
                 v3data = sData + rnorm(nVals, mean = 0, sd = 0.4))
dLongF = reshape(data=dF, varying=c("v1data","v2data","v3data"), v.names='data', direction="long", idvar="nSubject", times=1:3)

#Define function to assess different covariance structures
doAllCorrModels <- function(dataF, subjVarName, visitVarName, fixedPart){
    mList <- vector("list",2)
    mList[[1]] <- lme(fixedPart, #Random intercept, homogeneous variance
             random=as.formula(paste("~1|", subjVarName)),
             data=dataF,
             weights=NULL)
    mList[[2]] <- lme(fixedPart, #Random intercept, heterogeneous variance
             random=as.formula(paste("~1|", subjVarName)),
             data=dataF,
             weights=varIdent(form = as.formula(paste("~1|", visitVarName)))
    )
    mList
}

#Get different covariance structures
dataF <- dLongF
subjVarName <- "nSubject"
visitVarName <- "time"
fixedPart <- data ~ time
m <- doAllCorrModels(dataF, subjVarName, visitVarName, fixedPart)

#This works:
a1 <- anova(m[[1]], m[[2]])

#But this does not:
rm(fixedPart)
a2 <- anova(m[[1]], m[[2]])

罗兰

您可以使用do.call以下方法避免这种情况：

doAllCorrModels <- function(dataF, subjVarName, visitVarName, fixedPart){
  mList <- vector("list",2)
  mList[[1]] <- do.call(lme, list(fixed = fixedPart,
                                  random=as.formula(paste("~1|", subjVarName)),
                                  data=dataF,
                                  weights=NULL))

  mList[[2]] <- do.call(lme, list(fixed = fixedPart,
                                  random=as.formula(paste("~1|", subjVarName)),
                                  data=dataF,
                                  weights=varIdent(form = as.formula(paste("~1|", visitVarName)))))
  mList
}

本文收集自互联网，转载请注明来源。

如有侵权，请联系[email protected] 删除。