若第 i 个申请者得到研究所入学许可的机率 pi 与其 gpa 与 gre 有关,
模式为: log(pi/1-pi) = -6 + gpa + 0.005gre
知gpa ~ N (3.1,0.3) ,gre ~ N(580,80)。假设样本数为 30,模拟重复次数(number of repetitions)设为 1000 次,种子数指定为 1,2,...,1000。以蒙地卡罗模拟 来看回归系数 MLE 的偏差
提示:
先产生n个人的gpa及gre
给定 gpa 及 gre,我们可求出每个人得到入学许可的机率pi
以 sample(c(0,1),1,c(1-pi ,pi ),replace=F)指令来决定第 i 个申请者是否得到入学许可
产生完资料后,以牛顿法来解 MLE,请自行确认你的结果glm(admit~gpa+gre,family=binomial) 符合。
betaco <- c(-6,1,0.005)
n <- 30
Y <- c()
no.rep <- 1000
MLE <- matrix(NA,no.rep,3)
MLElm <- matrix(NA,no.rep,3)
for(i in 1:no.rep){
set.seed(i)
gpa <- rnorm(n,3.1,0.3)
gre <- rnorm(n,580,80)
random.error <- rnorm(n,0,1)
pi <- exp(X %% betaco)/(1+exp(X %% betaco))
admit <- sample(c(0,1),1,c(1-pi[n],pi[n]),replace=F)
.....目前到这卡住了,请高手帮忙,谢谢!