时间戳

首页 > 解决方案 > 使用重新构造（）后更新 nlme 模型

r - 使用重新构造（）后更新 nlme 模型

问题描述

nlme在公式参数中使用重新制定后更新模型有一个小问题lme()

这是一些数据

set.seed(345)
A0 <- rnorm(4,2,.5)
B0 <- rnorm(4,2+3,.5)
A1 <- rnorm(4,6,.5)
B1 <- rnorm(4,6+2,.5)
A2 <- rnorm(4,10,.5)
B2 <- rnorm(4,10+1,.5)
A3 <- rnorm(4,14,.5)
B3 <- rnorm(4,14+0,.5)
score <- c(A0,B0,A1,B1,A2,B2,A3,B3)
id <- rep(1:8,times = 4, length = 32)
time <- factor(rep(0:3, each = 8, length = 32))
group <- factor(rep(c("A","B"), times =2, each = 4, length = 32))
df <- data.frame(id = id, group = group, time = time,  score = score)

现在说我想将变量指定为lme函数外部的对象......

t <- "time"
g <- "group"
dv <- "score"

...然后重新制定它们...

mod1 <- lme(fixed = reformulate(t, response = "score"),
            random = ~1|id, 
            data = df)

summary(mod1)

Linear mixed-effects model fit by REML
 Data: df 
       AIC      BIC    logLik
  101.1173 109.1105 -44.55864

Random effects:
 Formula: ~1 | id
        (Intercept)  Residual
StdDev:   0.5574872 0.9138857

Fixed effects: reformulate(t, response = "score") 
                Value Std.Error DF   t-value p-value
(Intercept)  3.410345 0.3784804 21  9.010626       0
time1        3.771009 0.4569429 21  8.252693       0
time2        6.990972 0.4569429 21 15.299445       0
time3       10.469034 0.4569429 21 22.911036       0
 Correlation: 
      (Intr) time1  time2 
time1 -0.604              
time2 -0.604  0.500       
time3 -0.604  0.500  0.500

Standardized Within-Group Residuals:
       Min         Q1        Med         Q3        Max 
-1.6284111 -0.5463271  0.1020036  0.5387158  2.1784156 

Number of Observations: 32
Number of Groups: 8

到现在为止还挺好。但是，如果我们想使用 向模型的固定效应部分添加项update()怎么办？

mod2 <- update(mod1, reformulate(paste(g,"*",t), response = "score"))

我们收到错误消息

Error in reformulate(t, response = "score") : 
  'termlabels' must be a character vector of length at least one

显然我可以在不使用的情况下再次写出模型，update()但我只是想知道是否有办法使更新工作。

我认为问题在于lme使用reformulate.

非常感谢任何解决方案。

标签： rnlme

解决方案

问题是，当您不在对 的调用中输入公式文字时lme，某些类型的函数将不起作用。特别是，错误来自的地方是

formula(mod1)
#  Error in reformulate(t, response = "score") : 
#  'termlabels' must be a character vector of length at least one

nlme:::formula.lme尝试在错误的环境中评估参数。构建第一个模型的另一种方法是

mod1 <- do.call("lme", list(
  fixed = reformulate(t, response = "score"),
  random = ~1|id, 
  data = quote(df)))

当你这样做时，这会将公式注入到调用中

formula(mod1)
# score ~ time

这将允许update函数更改公式。

推荐阅读