關(guān)于“模型選擇（model selection）”，具體信息請查閱 Burnham 和 Anderson（2002）的原著，也可參考?Wolfgang Viechtbauer 和 斯幸峰 的博文，網(wǎng)址如下

Model selection and multimodel inference: a practical information-theoretic approach_Burnham &?Anderson_2002

https://caestuaries.opennrm.org/assets/06942155460a79991fdf1b57f641b1b4/application/pdf/burnham_anderson2002.pdf

Model Selection using the glmulti and MuMIn Packages_Viechtbauer（演示得很詳細(xì)）

https://www.metafor-project.org/doku.php/tips:model_selection_with_glmulti_and_mumin

模型選擇與多模型推斷_斯幸峰（需要科學(xué)上網(wǎng)?？茖W(xué)網(wǎng)也有作者的轉(zhuǎn)載）

http://sixf.org/cn/2014/03/model-selection-multimodel-inference/

1? 起因

今早收到“土壤微生物組”公眾號的推送，一篇 2020 年發(fā)在 GCB 上的 Meta 分析文章，配圖非常漂亮

Ji Chen et al., Soil carbon loss with warming: New evidence from carbon-degrading enzymes, Global Change Biology, 2020.?https://doi.org/10.1111/gcb.14986

圖（b）在其他 Meta 分析中也偶有出現(xiàn)。以“Sum of Akaike weights”這一指標(biāo)的大小來判斷不同變量的相對重要性雖然直觀，但也讓第一次接觸這一標(biāo)準(zhǔn)的人感到有些迷茫：這是個啥玩意兒？

2? 經(jīng)過

2.1? 模仿

參考 Viechtbauer 的博文，示例的部分 R 代碼如下

2.1.1? 數(shù)據(jù)準(zhǔn)備

? library(metafor)

? help(dat.bangertdrowns2004)

? dat <- dat.bangertdrowns2004? #? 導(dǎo)入數(shù)據(jù)

? rbind(head(dat, 10), tail(dat, 10))

? dat <- dat[!apply(dat[, c(“l(fā)ength”, “wic”, “feedback”, “info”, “pers”, “imag”, “meta”)], 1, anyNA), ]? #? 示例所納入的變量

2.1.2? “模型分析”

? library(glmulti)

? rma.glmulti <- function(formula, data, ...)

? ? rma(formula, vi, data = data, method = “ML”, ...)

? #? 這和一般的自定義函數(shù)寫法不同，沒搞清楚

? res <- glmulti(yi ~ length + wic + feedback + info + pers + imag + meta, data = dat, level=? 1, fitfunction = rma.glmulti, crit = "aicc", confsetsize = 128)

? #? 這一步即為“模型分析”。調(diào)用其結(jié)果即可得到變量的“相對重要性”值

? plot(res, type = “s”)? #? 以圖的形式表現(xiàn)“相對重要性”

2.1.3? 變量的相對重要性

按照 Viechtbauer 的演示，可以通過如下代碼整理各變量的相對重要性

? mmi <- as.data.frame(coef(res, varweighting = "Johnson"))

? mmi <- data.frame(Estimate = mmi$Est, SE = sqrt(mmi$Uncond),?

? ? ? ? ? ? ? ??? ? Importance = mmi$Importance, row.names = row.names(mmi))

? mmi$z <- mmi$Estimate / mmi$SE

? mmi$p <- 2*pnorm(abs(mmi$z), lower.tail = FALSE)

? names(mmi) <- c("Estimate", "Std. Error", "Importance", "z value", "Pr(>|z|)")

? mmi$ci.lb <- mmi[[1]] - qnorm(.975) * mmi[[2]]

? mmi$ci.ub <- mmi[[1]] + qnorm(.975) * mmi[[2]]

? mmi <- mmi[order(mmi$Importance, decreasing = TRUE), c(1, 2, 4:7, 3)]

? round(mmi, 4)?#? 結(jié)果如下

通過文獻(xiàn)的方法部分，容易得知“相對重要性”與“Akaike 權(quán)重/模型權(quán)重”有關(guān)。然而，“權(quán)重”是如何計(jì)算的？它是否可以被更直觀地體驗(yàn)到？

Chen 等（2020）解釋道，“簡而言之，特定變量的相對重要性被計(jì)算為包含該變量的所有模型的 Akaike 權(quán)重之和；該值被視為所有可能模型中各變量的總體支持度。”

這一解釋似乎還是有欠直觀，下面以變量 imag 為例。

2.2? 實(shí)踐

summary(res)$modelweights 命令會給出各模型的權(quán)重，注意到，各模型的權(quán)重之和為 1

> sum(summary(res)$modelweights)?

結(jié)果為

[1] 1

示例納入 7 個變量，理論上模型應(yīng)該有?2^7=128 (個)。不妨猜測，imag 的“相對重要性”、或說“Akaike 權(quán)重之和”即是包括 imag 變量在內(nèi)的模型的“模型解釋度”之和。借助中間字符匹配函數(shù) grepl( )

> res_modelweights <- weightable(res)

> sum(res_modelweights[which(grepl("imag", res_modelweights$model) == T), "weights"])

結(jié)果為

[1] 0.8478242

這與前述表格所列一致！變量的相對重要性/權(quán)重被計(jì)算為相應(yīng)的模型權(quán)重之和，這是一種奇妙的轉(zhuǎn)換。

3? 結(jié)果/理論

3.1? 斯幸峰?的解釋

斯幸峰 的解釋是，“（計(jì)算模型權(quán)重）得到各個模型的 AIC 值后，按照 AIC 從小到大排列，然后每個模型的 AIC 值與最小的AIC值相減，得到 ΔAIC。通過得到的 ΔAIC，計(jì)算各個模型的模型權(quán)重，即 Akaika weight(wi)。

……模型權(quán)重越大，表示該模型是真實(shí)模型的可能性就越大。比如第二個模型的 w2 為 0.31，則表示這個模型為真實(shí)模型（best possible model）的可能性為 31%。通過模型權(quán)重還可以計(jì)算各個參數(shù)的重要值（importance）。方法很簡單，比如參數(shù) 1，則挑出含參數(shù) 1 的所有模型，然后把這些模型的權(quán)重相加，即是該參數(shù)的權(quán)重。各個參數(shù)的權(quán)重值一比，就知道哪個參數(shù)最重要了?！?/p>

3.2??Burnham 和 Anderson 的推導(dǎo)

在 Burnham 和 Anderson（2002）的原著 Model selection and multimodel inference: a practical information-theoretic approach 中，有關(guān)計(jì)算如下

式中，i 表示相應(yīng)的“第 i 個模型”，R 是計(jì)算 Akaike 權(quán)重時納入的模型總數(shù)。