計算置信區間模型的平均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