将新数据传递给模拟.glmmTMB答案

【问题标题】：Passing New Data to simulate.glmmTMB将新数据传递给模拟.glmmTMB
【发布时间】：2019-03-18 21:14:38
【问题描述】：

在上一个问题 (Generate a predicted count distribution from a ZINB model of class glmmTMB) 中，我询问了如何为“glmmTMB”类的零膨胀负二项式模型生成预测计数分布。此后我发现的一个解决方案是函数模拟.glmmTMB (https://www.rdocumentation.org/packages/glmmTMB/versions/0.2.3/topics/simulate.glmmTMB)。但是，我想对测试数据进行模拟以验证模型的预测能力，我只了解如何对用于拟合模型的相同数据进行模拟。

在下面的示例中，我如何模拟 newdata 数据框的结果？

library(glmmTMB)
data("bioChemists", package = "pscl")
zinb <- glmmTMB(art ~ fem + mar + kid5 + phd + ment, ziformula = ~ ., data = 
bioChemists, family = nbinom2(link = "log"))
sim_1 <- simulate(zinb) #works as expected

#make new dataframe
newdata = unique(bioChemists[,c("fem","mar","kid5","phd","ment")])
sim_2 <- simulate(zinb, newdata = newdata) #ignores newdata

【问题讨论】：

这可能必须手动完成......
感谢@Ben Bolker 的快速回复。这会涉及到这个答案（stackoverflow.com/a/14967981/6365720）的内容吗？我试图将该解决方案外推到这种情况下，但在将其推广到更复杂的模型结构时遇到了麻烦。我很高兴花时间手动完成这项工作，但坦率地说，我不知道从哪里开始。再次感谢。

标签： r regression

【解决方案1】：

我认为这可行（可以封装在函数等中）：

n <- nrow(newdata)
## construct model matrices for conditional and Z-I
##  components (.[-2] drops the LHS of the formula)
X <- model.matrix(formula(zinb)[-2],newdata)
X_zi <- model.matrix(formula(zinb,component="cond")[-2],newdata)
## extract coefficients
beta <- fixef(zinb)$cond
beta_zi <- fixef(zinb)$zi
## draw random values for each component
cond <- rnbinom(n, mu=exp(X %*% beta), size=sigma(zinb))
zi <- rbinom(1, prob=plogis(X_zi %*% beta_zi), size=1)
cond*zi

最后一步有点太聪明了：ifelse(zi==0,0,cond) 可能更清楚，或者对于最后三个步骤，您可以使用 emdbook 包中的 rzinbinom 函数...

一般来说，我认为simulate() 方法应该允许newdata 和newparams - 它为参数引导、后验模拟等提供了广泛的可能性。等等。

更紧凑（可能更健壮）：

cond <- predict(zinb, newdata=newdata, type="conditional")
zi <- predict(zinb, newdata=newdata, type="zprob")
emdbook::rzinbinom(nrow(newdata),
                   mu=cond, size=sigma(zinb),
                   zprob=zi)

【讨论】：

感谢您抽出宝贵的时间来指导我完成这件事 - 我真的很感激。听到您同意 simulate() 应该允许 newdata 也是正确的。
@benbolker 您是否知道一种方法，但是当您对色散参数也有预测时？所以原始公式可能类似于glmmTMB(art ~ fem + mar + kid5 + phd + ment, dispformula = ~ fem, ...)。我找不到从生成的模型中提取对色散参数的预测的方法。
我知道做出这些预测的唯一方法是艰难的方法，即f <- formula(model, component="disp"); X <- model.matrix(f, data=newdata); beta <- fixef(model)[["disp"]]; eta <- X %*% beta; sd_pred <- exp(eta)
@BenBolker 您如何为在 glmmTMB 中估计的 nbinom1 模型实施此解决方案？我假设simulate.glmmTMB() 为 NB-1 模型提供计数这一事实表明某处存在某种类似于rnbinom1() 的函数，但所有模拟都被转移到该包中的 TMB 和 C 函数，我可以'不知道在哪里。