用户自定义函数用于lme模型拟合：错误

Question

用户自定义函数用于lme模型拟合：错误

3

我正在开始编写一个使用nlme构建线性混合模型的函数。我遇到了一个错误：Error in eval(expr, envir, enclos) : object 'value' not found，我认为这是由于R不知道在哪里找到数据框变量（例如value）。如果这确实是导致错误的原因，那么如何告诉函数value和timepoint属于下面（可重现）代码中的Dat变量？

require(nlme)
Dat <- data.frame(
    id = sample(10:19),
    Time = sample(c("one", "two"), 10, replace = T),
    Value = sample(1:10)
)
nlme_rct_lmm <- function (data, value, timepoint,
                      ID) {

    #base_level intercept only model
    bl_int_only <- gls(value ~ 1, 
                       data = data, 
                       method = "ML", 
                       na.action="na.omit")        
    #vary intercept across participants
    randomIntercept <- lme(value ~ 1, 
                           data = data, 
                           random = ~1|ID, 
                           method = "ML", 
                           na.action = "na.omit")       
    #add timepoint as a fixed effect
    timeFE <- lme(value ~ timepoint, 
                  data = data, 
                  random = ~1|ID, 
                  method = "ML", 
                  na.action = "na.omit")
}
nlme_rct_lmm(Dat, Value, Time, id)

- Santi Allende

1个回答

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

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

- Ben Bolker · Accepted Answer

这并不是我们预期的在不同框架内进行评估的问题，而是公式和数据之间变量名称一致性的问题。R是区分大小写的，因此您使用value或Value，id或ID等都很重要。此外，公式解释使用非标准评估（NSE），因此如果您有一个变量value等于符号Value，value ~ 1并不会神奇地变成Value ~ 1。下面我所概述的方法通过将响应、时间和ID变量的名称传递给函数来实现，因为这是最简单的方法。如果您使用非标准评估，对最终用户来说可能更加优雅，但这需要编程技巧（因此理解、调试等也更加困难）。

在简单/愚蠢的方法下面，我还讨论了如何实现NSE方法（向下滚动到底部...）

请注意，您的示例没有返回任何内容；在R中，这意味着完成函数后所有结果都将被丢弃。您可能希望将结果作为列表返回（或者您的实际函数将执行一系列模型测试等其他操作，并将这些答案作为结果返回...）。

require(nlme)

Dat <- data.frame(
    ID = sample(10:19),
    Time = sample(c("one", "two"), 10, replace = T),
    Value = sample(1:10)
)

nlme_rct_lmm <- function (data, value, timepoint,
                      ID) {

    nullmodel <- reformulate("1",response=value)
    fullmodel <- reformulate(c("1",timepoint),response=value)
    remodel <- reformulate(paste("1",ID,sep="|"))

    #base_level intercept only model
    bl_int_only <- gls(nullmodel,
                       data = data, 
                       method = "ML", 
                       na.action="na.omit")

    #vary intercept across participants
    randomIntercept <- lme(nullmodel,
                           data = data, 
                           random = remodel,
                           method = "ML", 
                           na.action = "na.omit")

    #add timepoint as a fixed effect
    timeFE <- lme(fullmodel,
                  data = data, 
                  random = remodel,
                  method = "ML", 
                  na.action = "na.omit")
}

nlme_rct_lmm(Dat, "Value", "Time", "ID")

如果你想要更加优雅（但内部较为难懂）的东西，你可以用以下代码替换定义模型的行。内部的substitute()调用检索作为参数传递给函数的符号；外部的substitute()调用将这些符号插入到公式中。

nullmodel <- formula(substitute(v~1,list(v=substitute(value))))
fullmodel <- formula(substitute(v~t,list(v=substitute(value),
                                 t=substitute(timepoint))))
remodel <- formula(substitute(~1|i,list(i=substitute(ID))))

现在，不需要将变量指定为字符串，就可以按照您的期望运行以下代码：nlme_rct_lmm(Dat, Value, Time, ID)。