如何解决尝试在没有来自 bam() 输出的随机影响的情况下进行预测时出错
我有一个数据集,我想在 mgcv 包中与 bam() 配合使用。该模型有一个二元结果,我需要为每个动物 ID 指定随机截距。下面是数据的一个子集(我的实际数据更大,协变量更多):
dat2 <- read.csv('https://github.com/silasbergen/example_data/raw/main/dat2.csv')
dat2$Animal_id <- factor(dat2$Animal_id)
> head(dat2)
Animal_id DEM_IA Anyrisk
1 105 279.94 0
2 105 278.68 0
3 106 329.13 0
4 106 329.93 0
5 106 332.25 0
6 106 333.52 0
> summary(dat2)
Animal_id DEM_IA Anyrisk
105: 2 Min. :156.3 Min. :0.0000
106: 83252 1st Qu.:246.8 1st Qu.:0.0000
107: 22657 Median :290.1 Median :0.0000
108:104873 Mean :284.8 Mean :0.3619
109:142897 3rd Qu.:318.0 3rd Qu.:1.0000
110: 53967 Max. :411.8 Max. :1.0000
我想拟合模型并预测没有随机效应的新数据:
library(mgcv)
mod <- bam(Anyrisk ~s(Animal_id,bs="re") + s(DEM_IA),data = dat2,family = "binomial",discrete=TRUE)
topred <- data.frame(DEM_IA = c(280,320))
predict(mod,newdata = topred,exclude="s(Animal_id)",newdata.guaranteed = TRUE)
但这会引发错误:
Error in eval(predvars,data,env) : object 'Animal_id' not found
当我特别告诉它从预测中排除这个词时,它为什么需要 Animal_id
?这也特别奇怪,因为我可以运行 ?random.effects
mgcv
帮助文件中的类似示例,没问题,即使我修改这些示例以使用 bam() 而不是 gam()!任何帮助将不胜感激!
编辑
我可能找到了解决办法;显然,如果在 discrete=TRUE
模型中使用 bam()
,那么 predict.bam()
也使用 discrete=TRUE
,这将不在缺少随机效应的情况下工作,但这有效:
mod<- bam(Anyrisk ~s(Animal_id,topred,newdata.guaranteed = TRUE,discrete=FALSE)
输出:
1 2
-0.4451066 -0.0285989
解决方法
tl;dr 通过将 something 放入 Animal_id
来解决此问题,您指定的值(不是 NA
虽然...)
为什么?如果不深入研究代码就不能确定,但是……使用 model.frame(formula,newdata)
作为计算所需模型矩阵的步骤通常很方便。 (例如,可以先构建整个模型矩阵,然后将要忽略的列归零……)弄清楚可以从公式中删除哪些项可能是一个单独的、更困难的步骤。 (我不知道为什么它在 bam
和 gam
中的工作方式不同......)
这似乎工作正常:
topred <- data.frame(DEM_IA = c(280,320),Animal_id=dat2$Animal_id[1])
predict(mod,newdata = topred,exclude="s(Animal_id)",newdata.guaranteed = TRUE)
检查您为 Animal_id
指定的内容是否真的无关紧要:
res <- lapply(levels(dat2$Animal_id),function(i) {
dd <- transform(topred,Animal_id=i)
predict(mod,newdata = dd,newdata.guaranteed = TRUE)
})
do.call(rbind,res)
结果:
1 2
[1,] -0.4451066 -0.0285989
[2,] -0.4451066 -0.0285989
[3,] -0.4451066 -0.0285989
[4,] -0.4451066 -0.0285989
[5,] -0.4451066 -0.0285989
[6,] -0.4451066 -0.0285989
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。