如何解决如何使用参数自举计算协方差矩阵?
使用以下数据集“数据”:load(url("https://www.math.ntnu.no/emner/TMA4315/2020h/hoge-veluwe.Rdata"))
我已经安装了Poisson GLM model = glm(y~t + I(t^2),family = poisson,data)
。我现在想通过使用带有1000个模拟的参数自举估计从GLM回归获得的β系数的协方差矩阵。到目前为止,我的代码是:
ysim = simulate(mod_quad,1)
betahat = matrix(0,nrow = 1000,ncol = 3)
for (i in 1:1000){
sim_data = cbind(ysim,data$t)
betahat[i,]= glm(ysim ~ data$t + I(data$t^2),data = sim_data )$coefficients
ysim = simulate(glm(ysim ~ data$t + I(data$t^2)family = poisson,data = sim_data ),1)
}
var(betahat[1000,])
每次都将betahat矩阵变成零矩阵,所以我不确定我的方法中缺少什么?
解决方法
这实际上是编码而不是统计问题。精简,我会这样做:
## set vals to NA (not 0) to make detection of problems easier
betahat <- matrix(NA,nrow = 1000,ncol = 3)
for (i in 1:1000) {
## replace response with parametric simulation
sim_data <- transform(data,y=simulate(mod_quad,1)[[1]])
## refit model with new data
newfit <- update(mod_quad,data=sim_data)
## store new coefficients
betahat[i,] <- coef(newfit)
}
## compute variance
var(betahat)
版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。