No.21042 dredgeとmodel.selのAICが異なる  【gibaco 】 2014/05/19(Mon) 12:15

済みませんが,教えて下さい。
以下のような,オフセット項,ランダム効果を加えたGLMM(glmmadmb)をつくり,AIC基準でモデル選択をしようとしています。

ana <- glmmadmb(monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)) + (1 | site), data =d, family= "nbinom")

モ デル選択をするためにdredgeとmodel.selをそれぞれ使って,モデル間の比較を行いました。ところが,この二つのやり方ではlogLik, dfの値が異なり,当然ながらAICも違っています。また,各モデルのAIC順も異なっていました。個々のモデルについて,summaryを出すと, model.selの結果とは一致し,dredgeとは一致しません。どうしてこんなことが起きるのか,どちらを使うべきか,教えていただければと思いま す。

 色々試してみたので以下に結果の例をのせておきます。うまくファイル添付できなかったので,長くなってしまいました。申し訳ありま せん。ここではmodel.selにはフルモデルしか入れていないので,dredgeの結果についてもフルモデル(8番のモデル)のみ注目して頂ければと 思います。
どうも,dredgeを使うと,ランダム効果がないものとしてAIC等を算出しているように見えます。ただし,dredgeの入っているMuMInパッケージ説明ではglmmadmbがつかえるようなことが書いてあります(http://cran.r-project.org/web/packages/MuMIn/MuMIn.pdf)。どう考えればいいのでしょうか。dredgeは使うべきではないのでしょうか?よろしくお願いします。

結果の例
その1 glmmadmbを使い,dredgeした場合

> ana <- glmmadmb(monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)) + (1 | site), data =d, family= "nbinom")
> dredge(ana, fixed = "offset(log(Rdcamt))", rank="AIC")
Global model call: glmmadmb(formula = monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)) +
(1 | site), data = d, family = "nbinom")
---
Model selection table
(Int) b15 mnD p15 off(log(Rdc)) df logLik AIC delta weight
4 -2.894 0.2754 -3.760 + 4 -850.338 1708.7 0.00 0.597
8 -2.965 0.2838 -4.018 0.030100 + 5 -850.138 1710.3 1.60 0.268
2 -2.991 0.2669 + 3 -853.137 1712.3 3.60 0.099
6 -2.995 0.2672 0.001188 + 4 -853.136 1714.3 5.60 0.036
5 -1.546 -0.155200 + 3 -878.674 1763.3 54.67 0.000
7 -1.488 -2.011 -0.144700 + 4 -878.002 1764.0 55.33 0.000
3 -1.644 -3.030 + 3 -882.215 1770.4 61.75 0.000
1 -1.752 + 2 -883.755 1771.5 62.83 0.000
Random terms (all models):
‘1 | site’

その2 glmmadmbを使い(その1と同じモデル),model.selを使って,フルモデルの統計量を算出。その1と異なる値になる。なお,summaryで算出されるAICや係数はこれ(その2)と同じ値でした。

> ana <- glmmadmb(monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)) + (1 | site), data =d, family= "nbinom")
> model.sel(ana, rank = AIC)
Model selection table
(Int) b15 mnD p15 off(log(Rdc)) df logLik AIC delta weight
ana   -3.114 0.2455 -3 0.08234 + 6 -834.01 1680 0 1
Random terms (all models):
‘1 | site’

その3 glmmadmb(その1)からランダム効果項を削除して計算し,dredgeした場合=>その1の結果と一致

> ana <- glmmadmb(monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)), data =d, family= "nbinom")
> dredge(ana, fixed = "offset(log(Rdcamt))", rank="AIC")
Global model call: glmmadmb(formula = monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)),
data = d, family = "nbinom")
---
Model selection table
(Int) b15 mnD p15 off(log(Rdc)) df logLik AIC delta weight
4 -2.894 0.2754 -3.760 + 4 -850.338 1708.7 0.00 0.597
8 -2.965 0.2838 -4.018 0.030100 + 5 -850.138 1710.3 1.60 0.268
2 -2.991 0.2669 + 3 -853.137 1712.3 3.60 0.099
6 -2.995 0.2672 0.001188 + 4 -853.136 1714.3 5.60 0.036
5 -1.546 -0.155200 + 3 -878.674 1763.3 54.67 0.000
7 -1.488 -2.011 -0.144700 + 4 -878.002 1764.0 55.33 0.000
3 -1.644 -3.030 + 3 -882.215 1770.4 61.75 0.000
1 -1.752 + 2 -883.755 1771.5 62.83 0.000
> model.sel(ana_deer, rank = AIC)

その4 glmmadmb(その1)からランダム効果項をむりやり削除して計算し,model.selを使って,フルモデルの統計量を確認=>その1とその3と結果一致

> ana <- glmmadmb(monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)) + (1 | site), data =d, family= "nbinom")
> model.sel(ana, rank = AIC)
Model selection table
(Int) b15 mnD p15 off(log(Rdc)) df logLik AIC delta weight
ana   -2.965 0.2838 -4.018 0.0301 + 5 -850.138 1710.3 0 1

混合モデルではなく,GLMで試してみました。
その5 glm.nbを使い(当然ランダム効果は入っていません),dredgeした場合=>その1,3,4と一致

> ana <- glm.nb(monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)), data =d)
> dredge(ana, fixed = "offset(log(Rdcamt))", rank="AIC") #
Global model call: glm.nb(formula = monkN ~ b150 + p150 + meanmnkhuntD + offset(log(Rdcamt)),
data = d, init.theta = 0.6566762117, link = log)
---
Model selection table
(Int) b15 mnD p15 off(log(Rdc)) df logLik AIC delta weight
4 -2.894 0.2754 -3.760 + 4 -850.338 1708.7 0.00 0.597
8 -2.965 0.2838 -4.018 0.030100 + 5 -850.138 1710.3 1.60 0.268
2 -2.991 0.2669 + 3 -853.137 1712.3 3.60 0.099
6 -2.995 0.2672 0.001188 + 4 -853.136 1714.3 5.60 0.036
5 -1.546 -0.155200 + 3 -878.674 1763.3 54.67 0.000
7 -1.488 -2.011 -0.144700 + 4 -878.002 1764.0 55.33 0.000
3 -1.644 -3.030 + 3 -882.215 1770.4 61.75 0.000
1 -1.752 + 2 -883.755 1771.5 62.83 0.000

その6 glm.nbを使い(その5と同じ),model.selを使って,フルモデルの統計量を算出=>その1,3,4,5と一致

> model.sel(ana, rank=AIC)
Model selection table
(Int) b15 mnD p15 off(log(Rdc)) df logLik AIC delta weight
ana   -2.965 0.2838 -4.018 0.0301 + 5 -850.138 1710.3 0 1

● 「統計学関連なんでもあり」の過去ログ--- 046 の目次へジャンプ
● 「統計学関連なんでもあり」の目次へジャンプ
● 直前のページへ戻る