r语言 - 计算glm系数矩阵、偏差和均方误差



>我有以下使用 R 的代码:

x0 <- matrix(rnorm(100,1))
x <- as.matrix(cbind("Intercept"=1, x0))
n <- dim(x0)[[1]]
z <- cbind(rep(1,n),x0)
p <- dim(x0)[[2]]+1
for(i in 1:n) {
  gstart <- glm(y~x0,family=binomial)$coef
}

我想在n样本中计算先前广义线性模型的估计值,并为n实例创建一个估计矩阵,然后计算biasmean square error,其中参数矩阵由以下代码给出:

n=100 #is the number of samples
parameter.mat<-cbind(rep(2,n),rep(2,n))  

我认为 您想检查glm返回的系数与平均非参数bootstrap系数之间的差异。下面的示例,首先给出了一种使用 boot package 的方法,然后使用循环(类似于您的问题)

# some example data - set seed for reproducibility
set.seed(1)
dat <- data.frame(y = rbinom(100, 1, 0.5),  x = rnorm(100))
# samples
n <- 1000
# glm estimates
mod <- glm(y ~ x, family="binomial", data=dat)$coef

# alternative method using boot package -----------------------------------
library(boot)
# function to extract model coefficients
f <- function(dat, i) glm(y ~ x, family="binomial", data=dat[i, ])$coef
# run bootstrap
set.seed(1)
boot(dat, f, R=n)  
# manual bootstrap  - sample with replacement -----------------------------
out <- vector("list", length=n)
for(i in 1:n) {
     newdat <- dat[sample(1:nrow(dat), , T), ]
     out[[i]] <- glm(y ~ x, family="binomial", data=newdat)$coef
     }
# matrix of bootstrap coefficients
bc <- do.call("rbind", out)
# bootstrap means
bc.mn <- colMeans(bc)
bias <- mod - bc.mn

最新更新