Пользовательская функция для установки модели lme: ошибка

Question

Я начинаю писать функцию, которая строит линейные смешанные модели с 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) 

Answer 1

Это не так (как вы и я оба ожидали) проблему с оценкой в ​​разных кадрах; скорее, это проблема согласованности между именами переменных между формулой и данными. R чувствителен к регистру, поэтому важно, используется ли или valueValue, id или ID, и т.д .. Кроме того, формула интерпретация использует нестандартную оценку (NSE), поэтому если вы имеете переменную value равен символу Value, value ~ 1не волшебным образом преобразуется в Value ~ 1. То, что я изложил ниже, работает, передавая функции переменных ответа, времени и идентификатора функции, поскольку это самый простой подход. Это немного более элегантно для конечного пользователя, если вы используете нестандартную оценку, но это немного сложнее программировать (и, следовательно, понимать, отлаживать и т. Д.).

Ниже легкий/boneheaded подход, я также обсудить, как реализовать 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)