计算置信区间模型的平均R中

Question

我试图使用基于谢弗物流曝光方法来运行巢生存模型中使用的收缩数据，2004年我有一系列的参数，并希望比较所有可能的模式，然后用Burnham和Anderson在2002年估算的模型平均参数。然而，我很难找出如何估计收缩调整参数的置信区间。

是否有可能估计置信区间使用收缩估计模型平均参数？我可以很容易地使用model.average $ coef.shrinkage提取模型平均参数的平均估计值，但我不清楚如何得到相应的置信区间。

任何帮助表示衷心感谢。我目前正在使用MuMIn软件包，因为我在链接功能方面遇到了AICcmodavg错误。

下面是我使用的代码的简化版本：

library(MuMIn) 

# Logistical Exposure Link Function 
# See Shaffer, T. 2004. A unifying approach to analyzing nest success. 
# Auk 121(2): 526-540. 

logexp <- function(days = 1) 
{ 
    require(MASS) 
    linkfun <- function(mu) qlogis(mu^(1/days)) 
    linkinv <- function(eta) plogis(eta)^days 
    mu.eta <- function(eta) days * plogis(eta)^(days-1) * 
    .Call("logit_mu_eta", eta, PACKAGE = "stats") 
    valideta <- function(eta) TRUE 
    link <- paste("logexp(", days, ")", sep="") 
    structure(list(linkfun = linkfun, linkinv = linkinv, 
      mu.eta = mu.eta, valideta = valideta, name = link), 
     class = "link-glm") 
} 

# randomly generate data 
nest.data <- data.frame(egg=rep(1,100), chick=runif(100), exposure=trunc(rnorm(100,113,10)), density=rnorm(100,0,1), height=rnorm(100,0,1)) 
    nest.data$chick[nest.data$chick<=0.5] <- 0 
    nest.data$chick[nest.data$chick!=0] <- 1 

# run global logistic exposure model 
glm.logexp <- glm(chick/egg ~ density * height, family=binomial(logexp(days=nest.data$exposure)), data=nest.data) 

# evaluate all possible models 
model.set <- dredge(glm.logexp) 

# model average 95% confidence set and estimate parameters using shrinkage 
mod.avg <- model.avg(model.set, beta=TRUE) 
(mod.avg$coef.shrinkage) 

如何提取任何想法/生成相应的置信区间？

由于 艾米

Answer 1

琢磨一下这一段时间后，我已经提出了在卢卡奇，P. M.，伯纳姆，K. P.，&安德森，D.R。（2009）基于等式5以下的溶液。模型选择偏差和弗里德曼的悖论。 年鉴统计数学，62研究所（1），117-125的。任何意见的有效性将不胜感激。

# MuMIn generated shrinkage estimate 
    shrinkage.coef <- mod.avg$coef.shrinkage 

# beta parameters for each variable/model combination 
coef.array <- mod.avg$coefArray 
    coef.array <- replace(coef.array, is.na(coef.array), 0) # replace NAs with zeros 

# generate empty dataframe for estimates 
shrinkage.estimates <- data.frame(shrinkage.coef,variance=NA) 

# calculate shrinkage-adjusted variance based on Lukacs et al, 2009 
for(i in 1:dim(coef.array)[3]){ 
    input <- data.frame(coef.array[,,i],weight=model.set$weight) 

    variance <- rep(NA,dim(input)[2]) 
    for (j in 1:dim(input)[2]){ 
    variance[j] <- input$weight[j] * (input$Std..Err[j]^2 + (input$Estimate[j] - shrinkage.estimates$shrinkage.coef[i])^2) 
    } 
    shrinkage.estimates$variance[i] <- sum(variance) 
} 

# calculate confidence intervals 
shrinkage.estimates$lci <- shrinkage.estimates$shrinkage.coef - 1.96*shrinkage.estimates$variance 
shrinkage.estimates$uci <- shrinkage.estimates$shrinkage.coef + 1.96*shrinkage.estimates$variance