R：在 nlme 函数中使用因子变量答案

【问题标题】：R: using factor variables in nlme functionR：在 nlme 函数中使用因子变量
【发布时间】：2017-10-17 16:55:30
【问题描述】：

library(nlme)
model <- nlme(height ~ (R0) + 1,
              data = Loblolly,
              fixed = list(R0 ~ 1),
              random = list(Seed = pdDiag(list(R0 ~ 1))),
              start = list(fixed = c(R0 = -8.5)))

这是一个只有 1 个固定效应参数的简单模型。这个模型很合适，但是当我想引入因子水平协变量（即年龄）时，我遇到了以下错误。

Loblolly$age2 <- as.factor(ifelse(Loblolly$age < 12.5, 0, 1))
model2 <- nlme(height ~ (R0 + age2) + 1,
              data = Loblolly,
              fixed = list(R0 ~ 1 + (age2)),
              random = list(Seed = pdDiag(list(R0 ~ 1))),
              start = list(fixed = c(R0 = -8.5, age2 = 1)))

Error in chol.default((value + t(value))/2) : 
  the leading minor of order 1 is not positive definite
In addition: Warning messages:
1: In Ops.factor(R0, age2) : ‘+’ not meaningful for factors
2: In Ops.factor(R0, age2) : ‘+’ not meaningful for factors
3: In Ops.factor(R0, age2) : ‘+’ not meaningful for factors

这似乎是一个语法错误，但我不确定如何修复它。

【问题讨论】：

您的示例中是否缺少某些内容？当我运行你的第一个块时，我得到：“eval 中的错误（predvars，data，env）：找不到对象'组'”
我不这么认为？我尝试在新的 R 会话中运行，它工作正常。
很确定group 是Loblolly 数据集中的一个变量
我得到同样的错误。
我检查了?Loblolly。它只有 3 个变量：身高、年龄和种子。当我运行你的代码时，我得到了同样的错误，找不到对象“组”。

标签： r model nlme

【解决方案1】：

首先，您的模型规范不正确：当您在fixed = list(R0 ~ 1 + (age2)) 中将固定效应定义为RO 时，您必须在模型定义中使用此定义。

模型拟合指令变为：

model2 <- nlme(height ~ (R0) + 1,
          data = Loblolly,
          fixed = list(R0 ~ 1 + (age2)),
          random = list(Seed = pdDiag(list(R0 ~ 1))),
          start = list(fixed = c(R0 = -8.5, age2 = 1)))

现在这会导致一条新的错误消息：

Error in nlme.formula(height ~ (R0) + 1, data = Loblolly, fixed = list(R0 ~  : 
  step halving factor reduced below minimum in PNLS step

请注意，nlme 有一个 verbose 参数（在我们的例子中信息量不大）。

但似乎在没有收敛时会发生此错误。在这种情况下，这是由于您的起始值不再适用于此模型规范。

我只是尝试了一组不同的值，例如：

model2 <- nlme(height ~ (R0) + 1,
               data = Loblolly,
               fixed = list(R0 ~ 1 + (age2)),
               random = list(Seed = pdDiag(list(R0 ~ 1))),
               start = list(fixed = c(R0 = 0, age2 = 30)), verbose=TRUE)

一个收敛并提供一个模型

> model2
Nonlinear mixed-effects model fit by maximum likelihood
  Model: height ~ (R0) + 1 
  Data: Loblolly 
  Log-likelihood: -305.1093
  Fixed: list(R0 ~ 1 + (age2)) 
R0.(Intercept)       R0.age21 
      12.96167       36.80548 

Random effects:
 Formula: R0 ~ 1 | Seed
        R0.(Intercept) Residual
StdDev:   0.0002761926 9.145988

Number of Observations: 84
Number of Groups: 14

【讨论】：

谢谢。为诸如age2之类的二元协变量选择起始值的合理方法是什么？您是如何得出这些具体的起始值的？
如前所述，它真的很接近随机值——我可能很幸运！为什么接近随机：只需查看应该具有正相关关系的变量名称（身高和年龄）。一种好方法是在没有随机效应的情况下拟合回归并从估计开始。