205 年 4 月 28 日 @ 統計モデリング 統計モデリング 第三回配布資料 文献 : A. J. Dobso ad A. G. Barett: A Itroducto to Geeralzed Lear Models. 3rd ed., CRC Press. J. J. Faraway: Etedg the Lear Model wth R. CRC Press. 配布資料の PDF は以下からも DL できます. 短縮 URL http://tyurl.com/lb7kb8 担当 : 田中冬彦
Chap. 3: データは主に Faraway から Google map から転載
Google map から転載 Locato
生データ 数値と記号のられつ ) Gve Data /2) > galapagos NS ES Area Aear Dst DstSC Elevato EM Baltra 58 23 25.09.84 0.6 0.6 NA 0 Bartolome 3 2.24 572.33 0.6 26.3 09 Caldwell 3 3 0.2 0.78 2.8 58.7 4 Champo 25 9 0.0 0.8.9 47.4 46 Coamao 2 0.05 903.82.9.9 NA 0 Daphe_Maor 8 0.34.84 8.0 8.0 NA 0 Darw 0 7 2.33 2.85 34. 290.2 68 Ede 8 4 0.03 7.95 0.4 0.4 NA 0 Ederby 2 2 0.8 0.0 2.6 50.2 2 Espaola 97 26 58.27 0.57. 88.3 98 Ferada 93 35 634.49 4669.32 4.3 95.3 494 Garder 58 7 0.57 58.27. 93. 49 Garder2 5 4 0.78 0.2 4.6 62.2 227 Geovesa 40 9 7.35 29.49 47.4 92.2 76 Isabela 347 89 4669.32 634.49 0.7 28. 707 Marchea 5 23 29.49 59.56 29. 85.9 343 Oslow 2 2 0.0 0.0 3.3 45.9 25 Pta 04 37 59.56 29.49 29. 9.6 777 Pzo 08 33 7.95 0.03 0.7 0.7 458 Las_Plazas 2 9 0.23 25.09 0.5 0.6 NA 0 Rabda 70 30 4.89 572.33 4.4 24.4 367 Sa_Crstobal 280 65 55.62 0.57 45.2 66.5 76 Sa_Salvador 237 8 572.33 4.89 0.2 9.8 906 Sata_Cruz 444 95 903.82 0.52 0.6 0.0 864 Sata_Fe 62 28 24.08 0.52 6.5 6.5 259 Sata_Mara 285 73 70.92 0.0 2.6 49.2 640 Seymour 44 6.84 25.09 0.6 9.6 NA 0 Tortuga 6 8.24 7.95 6.8 50.9 86 Wolf 2 2 2.85 2.33 34. 254.7 253 CRAN alr3 パッケージ内, galapagos データセット alr3 は, S. Wesberg による.
Gve Data 2/2) データの説明 galapagos 各島ごとの亀の種類 サンプルサイズ 29) 島の地理情報 NS=Number of speces, 島で観測された亀の種類 Area= 島の面積 [hr], Aear = 一番近い島の面積 [hr], Dst= 一番近い島との距離 [km], DstSC=Sata Cruz 島からの距離 [km], Elevato = 島の高度 [m] > galapagos NS Area Aear Dst DstSC Elevato Baltra 58 25. 09.84 0.6 0.6 NA Bartolome 3.24 572.33 0.6 26.3 09 Caldwell 3 0.2 0.78 2.8 58.7 4 データの図示 重要!) 見てわかること 定量的な確認
散布図 0 2000 4000 0 500 500 散布図 > parsgala.cor); NS 0 00 300 はずれ値? 散布図をみると Area でひとつでかいのが目立つ Isabera slad 0 2000 4000 Area Dst 0 0 20 30 40 0 500 500 Eleva 0 00 300 0 0 20 30 40
Chap. 3: データは主に Faraway から 面積 Area) の はずれ値! Google map から転載
Processed Data データの説明 gala.test = galapagos データの一部を削除 説明の都合 ) NS=Number of speces, 島で観測された亀の種類 Area= 島の面積 [hr], Aear = 一番近い島の面積 [hr], Dst= 一番近い島との距離 [km], DstSC=Sata Cruz 島からの距離 [km], Elevato = 島の高度 [m] > gala.test NS Area Aear Dst DstSC Elevato Bartolome 3.24 572.33 0.6 26.3 09 Caldwell 3 0.2 0.78 2.8 58.7 4 Champo 25 0.0 0.8.9 47.4 46... Sata_Mara 285 70.92 0.0 2.6 49.2 640 Tortuga 6.24 7.95 6.8 50.9 86 Wolf 2 2.85 2.33 34. 254.7 253 分析の課題 島ごとの亀の種類を説明する統計モデルを考える まずは線形モデルで )
線形モデルでむりやり解析 線形回帰 > gala.lm.res <- lmns~., data=gala.test); 線形回帰モデル 5 Y = α + β + ε = ε ~ 2 N0, σ ) あてはめた値 Ftted Value) 5 yˆ = ˆ α + = ˆ β α = 29.6, ˆ β = 0.285, ˆ β = 0.0667, ˆ β ˆ 2 3 = ˆ β = 0.33, ˆ β 4 5 = 0.44 0.290,
線形モデルでの残差プロット 回帰式は 6 次元なので図示できない 残差プロットでモデルのよさを検討 y yˆ Galapagos Tortose あてはめた値 Ftted Value) yˆ = ˆ α + 5 = ˆ β Resduals -50 0 50 00 残差 Resduals) y ˆ y 0 00 200 300 400 Ftted ŷ > plotpredctgala.lm.res), resdualsgala.lm.res), lab="ftted", ylab="resduals", ma="galapagos Tortose" ); > ableh=0, col="red"); モデルが正しいなら, 残差は平均 0 分散一定の正規分布から発生しているようにみえるはずだが
一般化線形モデルの導入 問題点 NS=Number of speces, 島で観測された亀の種類 離散値 Categorcal Data) 残差分析の結果 分散が一定とはいいがたい モデルが不適切であることを示唆 ) 分析の課題 線形モデルよりもよいモデルを考える Remark * 一般にAICや検定統計量の計算以前に残差はチェック 相関, 等分散 etc.)
ポアソン回帰モデル /2) 一般のポアソン回帰モデル Y Po µ ) ~ E[ Y ] = µ =,, log µ = p = β 説明変数の意味 µ µ = 0, = ) = = 0) 男 女, 喫煙 喫煙なし, etc.) e β 他の条件が同じ場合 説明変数が連続量でも同様に β を解釈できる.
ポアソン回帰モデル 2/2) ポアソン回帰モデルでの解析例 > gala.glm.res <- glmns~., famly=posso, data=gala.test); ポアソン回帰モデル Y Po µ ) ~ E[ Y ] = µ log µ = α + 5 = β あてはめた値 Ftted Value) 5 log ˆ µ = ˆ α + = ˆ β ˆ 4 ˆ 4 α = 3.48, β 5.26 0, 5.26 0, ˆ = β2 = β3 ˆ 3 ˆ 3 β4 = 4.55 0, β = 2.4 0 3 ˆ 5 = 7.00 0,
ポアソン回帰の残差プロット /3) 残差 * いくつかある ) Devace resduals スケーリングされた残差の一種 ) r = sg y ˆ µ ) DR d = d y, ˆ µ ) = 2 y log y / ˆ µ y + ˆ µ ) d Devace Resduals -5 0 5 0 Galapagos Tortose 00 200 300 400 ^ µˆ > plotpredctgala.glm.res, type="respose"), resdualsgala.glm.res), lab=epressohatmu)), ylab="devace Resduals", ma="galapagos Tortose"); > そのままのスケールだと, 左によってしまっている log mu) を横軸にとって残差プロットする
ポアソン回帰の残差プロット 2/3) Galapagos Tortose 残差 * いくつかある ) Devace resduals r = sg y ˆ µ ) DR d = d y, ˆ µ ) = 2 y log y / ˆ µ y + ˆ µ ) d Devace Resduals -5 0 5 0 3.0 3.5 4.0 4.5 5.0 5.5 6.0 ^ log > plotpredctgala.glm.res, type= lk"), resdualsgala.glm.res), lab=epressohatlogmu))), ylab="devace Resduals", ma="galapagos Tortose"); > 特徴的なパターンは見つからない
ポアソン回帰の残差プロット 3/3) Galapagos Tortose 残差 * いくつかある ) Respose resduals r = µˆ RR y ポアソン分布の分散 Var Y 5 ) = µ = + ep α β = Respose Resduals -00-50 0 50 00 50 3.0 3.5 4.0 4.5 5.0 5.5 6.0 ^ log >plotpredctgala.glm.res, type="lk"), resdualsgala.glm.res, type="respose"), lab=epressohatlogmu))), ylab="respose Resduals", ma="galapagos Tortose"); 実はもっと良いモデルが作れる 詳細は Faraway Chap.3, Chap.6 を参照 )
Offset ポアソン回帰の注意点 N 例 : 各地区ごとに人口 ) が違う場合の患者数 p p µ = N ep β = ep log N + β = = Y というモデル化をする offset term kow costat) Overdsperso 本来, ポアソン分布は期待値と分散は一致 Var Y ) = E[ ] = µ Y しかし, 実際のデータは分散が大きくみえる 負の二項分布でモデル化できる * 詳細は文献を参照
一般化線形モデル
Epoetal Famly /2) Y 確率変数 θ θ,, θ ) = パラメータ p Y の確率密度関数 もしくは確率分布 ) が以下の形でかける時 変数の ) 指数型分布族 Epoetal Famly) という. p y θ ) = s y) t θ )ep a y) b θ )) = ep a y) b θ ) + c θ ) + d y)) s y) > 0, t θ ) > 0 Caocal form p z θ ) = ~ s z) t θ )ep z b θ )) Natural Parameter η : = b θ ) 多くの確率分布が指数型分布族になっている.
)) ) )ep ) ) θ θ θ b y a t y s y p = )) ) ) ) ep y d c b y a + + = θ θ 例 :,) N m 2 ) 2 2 ) m e m p = π 平均 m, 分散の正規分布 ガウス分布 ) σ 2 = = ) log2 2 2 2 ep 2 2 π µ m 二項分布 回の試行, は固定 ),,,2, = 0 q q q p = ) ) = q q q log ep ) Epoetal Famly 2/2)
独立な確率変数 Geeralzed Lear Models /2) Y,,Y 各々, 同じ指数型分布族に従うとき 変数の ) 一般化線形モデル Geeralzed Lear Models; GLM) という. p y θ ) = ep yb θ ) + c θ ) + d y )) = p y θ ) ep yb θ ) + c θ ) + = = ここで興味あるパラメータは = E[ Y ] = µ GLM では共変量 説明変数 ) で以下を仮定 g µ ) = β p = = d y gµ) lk fucto リンク関数 ) )
ポアソン回帰の例 p y 以下を仮定 Geeralzed Lear Models 2/2) µ ) = µ ) y! ここで興味あるパラメータは log µ ) = p = y e β µ E[ Y ] = µ 以下の対数尤度関数をについて数値最大化 最尤推定 ) = = log p y µ ) = µ + y log µ β,, =,2,, = GLM ではプログラムによる数値解法が前提になっている β p = log y!
統計モデリング ~ 次回予告 Schstosoma apocum の卵 *) Schstosoma apocum の成虫 イメージイラスト ) 成虫 Schstosoma apocum) の写真は日本住血吸虫で検索してください 閲覧注意 *http://www.dpd.cdc.gov/dpd/html/schstosomass.htm