収縮を使用してモデル平均データの信頼区間を計算するR

Question

Shaffer、2004に基づくロジスティック曝露法を使用して巣生存モデルを実行しようとしています。さまざまなパラメータがあり、バーナムとアンダーソン（2002年）のように収縮を使用してモデル平均のパラメータを推定する。しかし、収縮調整されたパラメータの信頼区間をどのように見積もるかはわかりません。

収縮を使用して推定されたモデル平均パラメータの信頼区間を推定することは可能ですか？ model.average $ coef.shrinkageを使用して、モデル平均化されたパラメータの平均推定値を縮小して簡単に抽出できますが、対応する信頼区間を取得する方法は不明です。

ご協力いただきありがとうございます。私は現在、リンク機能に関してAICcmodavgでエラーが発生するので、MuMInパッケージを使用しています。

以下は私が使用しているコードの簡易版は、次のとおりです。

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