glmer 预测与 allow.new.levels=TRUE答案

【问题标题】：glmer prediction with allow.new.levels=TRUEglmer 预测与 allow.new.levels=TRUE
【发布时间】：2015-10-07 04:48:28
【问题描述】：

我查看了Prediction with lme4 on new levels，其中引用了 R 文档的 allow.new.levels=TRUE。

我不明白“如果 allow.new.level=TRUE，那么预测将对具有先前未观察到的水平（或 NA）的数据使用无条件（人口水平）值”的实际含义。

就我而言，我有一个名为 exdata 的 data.frame

 subjects y        x1         x2
   1 0 1.6179339  0.9517194
   1 0 1.4128789  1.0248514
   1 0 0.9127448  1.8073684
   1 0 1.5729219  2.1003925
   1 0 1.6254359  1.2471660
   2 0 1.6626074  8.5102559
   2 0 1.3903638  6.0425018
   2 0 1.1438239  2.5654422
   2 0 1.1393088  2.9982242
   2 0 1.1564141  2.8395960
   3 0 1.1688192 13.9791461
   3 0 0.9255715 18.8544778
   3 0 1.2369097  4.2376671
   3 0 1.3021943  9.6894289
   3 0 1.2296961 12.4789910
   4 0 1.0978131  2.0577688
   4 0 1.1405409  1.4339044
   4 0 1.0355546  1.9496732
   4 0 1.1370849  1.7402332
   4 0 1.1942591  1.3509880
   5 0 0.4141535  2.1723957
   5 0 0.9129311  0.8274350
   5 0 0.9658796  1.2754419
   5 0 0.8370701  2.1998756
   5 0 0.5509546  2.3590774
   6 0 1.2827411  1.5474088
   6 0 1.1636606  0.7746669
   6 0 1.1782936  1.1566909
   6 0 1.1630238  1.7486415
   6 0 1.1565711  0.6984409
   7 0 0.8600331  0.1382253
   7 0 0.8303510  0.1927431
   7 0 0.8087967  0.6065926
   7 0 0.7815187  0.9464185
   7 0 0.7532042  0.9771646
   8 0 1.1638190  1.3456340
   8 0 0.5867126  1.4862727
   8 0 0.6523964  0.5138441
   8 0 0.9513971  2.3932337
   8 0 0.9278743  2.3273670
   9 1 1.0978606  1.2585635
   9 1 1.0414897  1.2946008
   9 0 0.6215353  0.2907148
   9 0 1.0267046  1.0173432
   9 0 1.1470992  0.7014759
  10 0 0.9505266  0.4247866
  10 0 0.8624758  0.2276577
  10 0 0.8279061  0.2314898
  10 0 0.7856832  0.3143003
  10 0 0.7569739  0.7880622

有 10 个受试者，y 作为响应 (0-1) 和两个解释变量。我想将其建模为 logit glm 随机效应模型，所以我使用了

model <- glmer(y~x1+x2+ (1|subjects), family=binomial(link="logit"), data=exdata[exdata$subjects!=5,], nAGQ=1)

我只使用前 4 个主题，因为我想使用最后一个进行预测。模型的结果是随机效应对 82.55 的方差和估计值

(Intercept) -14.8377
x1            4.0366
x2            0.1056

对新主题的预测是

prediction <- predict(model, newdata=exdata[exdata$subjects==5,], type="response", allow.new.levels=TRUE)

给予

2.408299e-06 1.564704e-05 2.031392e-05 1.331586e-05 4.266690e-06

这些值是如何预测的？我的希望是计算 subject=5 的 ranef() ，并使用模型中的系数来计算概率。我可以以任何方式打印 subject=5 的 ranef() 吗？由于我不了解 allow.new.levels 的 R 文档，我无法理解预测是如何做出的。

【问题讨论】：

标签： r predict lme4

【解决方案1】：

在这个特定示例中，仅使用模型的固定部分（因为这是具有单个随机效应的模型的总体级别模型）：

y1 <- as.matrix(cbind(1, exdata[exdata$subjects==5, c("x1", "x2")])) %*% fixef(model) 
c(exp(y1)/(1+exp(y1)))
#[1] 2.408298e-06 1.564704e-05 2.031392e-05 1.331586e-05 4.266689e-06

【讨论】：

对，这意味着如果她想为主题 5 生成包含该主题随机效应的预测，她需要使用现在扩展的样本重新估计原始模型，其中包括对主题的观察5，然后使用该模型运行predict()。
另一种可能性（相当难以实施！）是对主题 5 进行预测，条件是 (1) 固定效应参数的所有； (2) 与其他 9 名受试者的估计受试者间方差； (3) 主题 5 的观测数据——即收缩估计/预测。我不确定它是否真的有意义，但你可以想象得到。与@ulfelder 建议的不同之处在于，在包含主题 5 时，您不会重新估计主题间的差异......

【解决方案2】：

所以这意味着 lme4（或另一个 R 包）没有一个函数可以预测新级别，但仍然估计新级别的随机效应..？ @ulfelder 的建议是最简单的，但是每次添加新主题时它都会改变我的估计，并且主题 5 的预测是基于一个模型，该模型使用来自主题 5 的数据。我想避免这种情况，因为然后，该模型被设计为预测其自身的值，可以这么说。 @BenBolker 避免了这种情况，但我不知道如何实现。

您如何看待结合这两种方法。

使用包含主题 5 观察结果的数据运行新模型 (model1)。
为主题 5 找到 ranef()。
像@Roland 那样计算预测结果（仍然使用初始模型的估计值等）
添加随机效果。

那是

model1 <- glmer(y~x1+x2+ (1|subjects), family=binomial(link="logit"), data=exdata, nAGQ=1)
y1 <- as.matrix(cbind(1, exdata[exdata$subjects==5, c("x1", "x2")])) %*% fixef(model) 
y2 <- y1 + unlist(ranef(model1))[5] 
c(exp(y2)/(1+exp(y2)))

2.398097e-06 1.558076e-05 2.022787e-05 1.325946e-05 4.248617e-06

不幸的是，与模型相比，模型 1 中的估计值（以及随机效应）发生了变化，但我试图避免在预测主题 5 时使用这些变化。

【讨论】：