データ解析のための統計モデリング入門６章 / Handbook-of-statistical-modeling-for-data-analysis-section6

Slide 1

Slide 1 text

統計モデリング入門６章担当：masso データラーニングギルド輪読会

Slide 2

Slide 2 text

本章の目的さまざまなGLMを紹介する章 - ポアソン回帰以外でどんなGLMがあるのか？ - それらはどのような特徴やありがたみがあるのか？などを理解すること

Slide 3

Slide 3 text

目次 1. 様々な種類のデータで応用できるGLM（6.1節） 2. 二項分布で表現する「あり・なし」カウントデータ（6.3節） 3. ロジスティック回帰とロジットリンク関数（6.4節） 4. 交互作用項の入った線形予測子（6.5節） 5. 割り算値の統計モデリングはやめよう（6.6節） 6. 正規分布とその尤度／ガンマ分布のGLM（6.7節／6.8節）

Slide 4

Slide 4 text

6.1節

Slide 5

Slide 5 text

様々な種類のデータで応用できるGLM（6.1節）まとめ GLMでは応答変数のばらつきを表現する確率分布は正規分布だけでなく、ポアソン分布・二項分布・ガンマ分布などが選択できる

Slide 6

Slide 6 text

GLMとは「確率分布」「リンク関数」「線形予測子」の組み合わせで表現する。つまり、これらの組み合わせの数だけ、GLMは存在するということ。とはいえ、「よく使われるパターン」があるので紹介する。

Slide 7

Slide 7 text

その前にリンク関数と線形予測子おさらい線形予測子とは、説明変数xiの（べき乗も含む）線形結合のこと *P47 （{βi}をベクトル、{xi}を係数とした線形結合）線形予測子 = β1 + β2 xi + β3 xi^2 リンク関数fとは、応答変数の平均値と線形予測子をつなぐ関数 f (応答変数の平均値) = β1 + β2 xi + β3 xi^2 ←自信ない…どっちがベクトル？

Slide 8

Slide 8 text

GLM よくある組み合わせ

Slide 9

Slide 9 text

直線回帰とリンク関数を用いた手法は違う GLMなど使わなくても、変数変換して、直線回帰をすればいいと考える人もいるかもしれないが、このような方法はリンク関数とは全く別のものである。（P114/115）「変数変換」というのが何をさしているかよくわからないが、直線回帰ということは、パラメータが「切片」「傾き」しかない。状態これでは、もし説明変数が２個以上あった場合でも、その自由度を無視することになる。（ということを言っているのかなと解釈した）

Slide 10

Slide 10 text

6.2節

Slide 11

Slide 11 text

6.2節例題：上限のあるカウントデータ上限のあるカウントデータとは？（応答変数y∈{0,1...N}） → 「N個体の内y個が○、N-y個が□」 -------------------------------- 本例題の説明　個体：i∈{1, 2, … 100} 個体iに対して　観察種子数 8個（Ni）　生存種子数 yi個 ∈{0...8} 　生存確率 qi 個体iの1個の種子が生きている確率 i ∈{1...50}：施肥処理ありﬁ=C i ∈{51...100}：施肥処理なしﬁ=T

Slide 12

Slide 12 text

6.2節例題：上限のあるカウントデータ summary(d)の出力結果データ可視化施肥処理と生存確率の関係体サイズと生存確率の関係

Slide 13

Slide 13 text

6.3節

Slide 14

Slide 14 text

二項分布で表現する「あり・なし」カウントデータ（6.3節）まとめ「N個の観察対象のうちk個で反応がみられた」というタイプのデータにみられるばらつきを表すために二項分布が使えるというより、「上限のあるデータではポアソン分布は適さない」

Slide 15

Slide 15 text

ちょっと脱線：二項分布とは

Slide 16

Slide 16 text

6.2節の例題を二項分布で表現

Slide 17

Slide 17 text

6.4節

Slide 18

Slide 18 text

ロジスティック回帰とロジットリンク関数（6.4節）まとめ生起確率と線形予測子を結びつけるロジットリンク関数を使ったGLMのあてはめは、ロジスティック回帰とよばれる

Slide 19

Slide 19 text

ロジット関数はロジスティック関数の逆関数 https://mathwords.net/logitkansu

Slide 20

Slide 20 text

ロジットリンク関数ロジット関数をリンク関数としたもの。ロジット関数をリンク関数として利用すると・線型予測子と[0, 1]の数（個体iの生存確率qi）を繋げることができる

Slide 21

Slide 21 text

最も当てはまりの良いパラメータ {βj } を推定尤度関数：L(y | θ) = Π{p(y | θ)} 対数尤度関数：log {L(y | θ)} = log [Π{p(y | θ)}]

Slide 22

Slide 22 text

最も当てはまりの良いパラメータ {βj} を推定 Rでの命令文 > glm( cbind( y, N-y ) ~ x + f, data = d, family = binomial ) 応答変数の指定：cbind( 生存数, 死亡数 ) 二項分布を指定：family = binomial

Slide 23

Slide 23 text

オッズとは？ logit関数のlogの中身をオッズという。（qi / (1 - qi)の部分）したがって、ロジットリンク関数と線型予測子の関係性より、確率qiが0.5（五分五分）のときオッズは1倍、qiが0.8ならオッズは4倍

Slide 24

Slide 24 text

ロジットリンク関数とオッズロジットリンク関数のおかげで、「オッズ ∝ 要因の掛け算」で表現できるため、パラメータを最尤推定などで決めて、固定した後に、「説明変数が1単位変動した時に、オッズがどれぐらい変動するか」といった解釈が容易になる。 ※要因＝パラメータβ×説明変数結果として、オッズ比 ≒ リスクの概算がしやすくなる

Slide 25

Slide 25 text

オッズ比 ≒ リスクの概算生活習慣Ｘによってナントカ病の発病リスクが７倍になります発病確率をロジスティック回帰で求め、最尤推定値βs=1.95となったとする。この時のリスク≒オッズ比は、 ≒ 7.028

Slide 26

Slide 26 text

6.5節

Slide 27

Slide 27 text

交互作用項の入った線形予測子（6.5節）まとめ線形予測子の構成要素として、複数の説明変数の積の効果をみる交互作用項が使える

Slide 28

Slide 28 text

交互作用項によってより複雑な状態を表す体サイズxiと施肥処理有無ﬁが組み合わさった場合の影響度を表す Rの命令文では、モデル式が変わる > glm( cbind( y, N - y ) ~ x * f, family = binomial, data = d ) x * f は、x + f + x:f の省略形

Slide 29

Slide 29 text

交互作用項の導入は必ずしも良い結果につながらない今回の例題では、交互作用項をいれても生存種子数のモデルはあまり変わらない。（下図） AICを計算すると、むしろ悪化する。（272 -> 274）

Slide 30

Slide 30 text

交互作用項を取り扱う際の注意点 ● むやみに追加しない。なぜなら、 ○ 説明変数が多い場合、「組合せ論的爆発」で増加していく ○ それが何を表しているのか解釈ができなくなることがある ● 現実問題では、交互作用項を多く含むモデルのAICが最良になることがよくある。しかし、 ○ 交互作用の効果を過大推定してる可能性あり（つじつま合わせ） ○ 現実では、説明変数では説明できない「個体差」「場所差」が発生するが、それらを考慮しないGLMをあてはめると過度に複雑なモデルが最良になる傾向がある。

Slide 31

Slide 31 text

6.6節

Slide 32

Slide 32 text

割り算値の統計モデリングはやめよう（6.6節）まとめデータ解析でしばしばみられる観測値同士の割り算値作成や応答変数の変数変換の問題点があるが、ロジスティック回帰やオフセット項の工夫をすれば、情報消失の原因となる「データ加工」は不要になる

Slide 33

Slide 33 text

観測値をこねくり回して指標を創作しないようによくある創作・割算値：観測データ / 観測データ・変数変換：応答変数 = log( 観測データ ) とか = avg( 観測データ ) 「N個のうちy個で事象が生じる確率」を明示的に扱う二項分布を使うことによって、「y / N」などといった観測データ同士の割算を避けられる。

Slide 34

Slide 34 text

観測データ同士の割算値がもたらす悲劇 ● 情報の消失 ○ 1000打数300安打の3割打者と10打数3安打の3割打者は、どちらも同じ程度に確からしい「3割打者」ではない。確からしさが失われる。 ● 変換された値の分布が不明 ○ 分子分母それぞれに誤差が入った数量同士の割算値はどんな確率分布か？（一般に難しい問題。分子分母が独立でない場合は、さらにややこしくなる。）

Slide 35

Slide 35 text

割算値いらずのオフセット項わざ例題 ● 森林のあちこちに100箇所で調査 i ∈ {1, 2, … 100} ● 調査値 i における面積 Ai （ほんとは固定にすべきだけど…あえて） ● 調査地 i の「明るさ」xi ● 調査地 i における植物個体数 yi を記録＝応答変数 ● O：調査地 i における植物個体の人口密度が明るさxiにどう影響されるかを知りたい yi / Ai という割算値を作る必要はない！

Slide 36

Slide 36 text

ポアソン分布＋オフセット項で人口密度をモデル化人口は正の量なので、指数関数と明るさ xiの依存性を組み合わせて Aiを指数関数の中にいれてやれば、線型予測子：β1+β2xi+logAi リンク関数：対数リンク関数 logAiは、係数がない＝オフセット項

Slide 37

Slide 37 text

オフセット項わざの使い所 ● GLM（とそれを発展させた統計モデル）で応用可能 ● 「単位○○あたりのカウントデータ」に使える ○ ○○（例えば面積）の対数を線型予測子に追加する ● もっと一般化すると、「連続値 / 連続値」となる比率・密度などに使える、と言える。 ● 分子分母共に誤差を含む場合は、ベイズ統計モデルで工夫すれば観測値どうしの割算を回避できる（本書の対象外）

Slide 38

Slide 38 text

6.7節 / 6.8節

Slide 39

Slide 39 text

正規分布とその尤度／ガンマ分布のGLM（6.7節／6.8節）まとめ連続値の確率変数のばらつきを表現する確率分布としては、正規分布・ガンマ分布などがあり、これらを統計モデルの部品として扱うときには、離散値の確率分布（ポアソン分布や二項分布）との違いに注意しなければならない

Slide 40

Slide 40 text

正規分布とその尤度 ● 平均μ、標準偏差σをパラメータとする連続値の確率分布を表す f(x) は、条件付き確率の形式で記述すると、p( y | μ, σ ) ● これをRで図示するならば > y <- seq( -5, 5, 0.1 ) > plot( y, dnorm( y, mean = 0, sd = 1 ), type = “l” )

Slide 41

Slide 41 text

（連続変数）確率密度関数の尤度計算法確率＝確率密度 × ⊿y として、尤度を計算するけど、最終的に⊿yの項は無視する省略する

Slide 42

Slide 42 text

（連続変数）確率密度関数の尤度計算法結果として、「対数尤度が正になり得る」「AICや逸脱度が負の数になり得る」といった離散型確率分布とは異なる性質を持っている。 ※幅⊿y = 0.5としている

Slide 43

Slide 43 text

正規分布 : 最小二乗法＝最尤推定 σがμと無関係かつ定数であるとした場合「対数尤度関数の最大化」 = 「Σ(yi - μ)^2 の最小化」になる最尤推定最小二乗法確率分布：正規分布、リンク関数：恒等関数、線型予測子：β1+β2 xi

Slide 44

Slide 44 text

ガンマ分布 ● 0以上の連続値の確率分布 ● shapeパラメータs、rateパラメータrから成る ● 平均：s / r　分散：s / r^2　※s=1のときは指数分布 ● Γ(s)はガンマ関数　（Wikipedia 「ガンマ関数」）

Slide 45

Slide 45 text

ガンマ分布 Rの命令文 > dgamma( y = ○, shape = ○, rate = ○ )

Slide 46

Slide 46 text

ガンマ分布のGLM ● 確率分布関数：ガンマ分布 p(y | s, r) ● リンク関数：対数かな？ ● 線型予測子：問題次第

Slide 47

Slide 47 text

例題：ガンマ分布でGLMをする本例題の説明　個体：i ∈ {1, 2, … 50} 各個体 i に対して　葉の重量：xi　←説明変数　花の重量：yi　←応答変数

Slide 48

Slide 48 text

例題：ガンマ分布でGLMをするなんらかの生物学的根拠により　μi = Axi^b と表せるとする対数リンク関数を用いて、線型予測子は、 log( μi ) = a + b log(xi) ※A = exp(a) と置いた Rの命令文 > glm( y ~ log(x) , family = Gamma(link = “log”), data = d )

Slide 49

Slide 49 text

例題：ガンマ分布でGLMをするなんらかの生物学的根拠により　μi = Axi^b と表せるとする対数リンク関数を用いて、線型予測子は、 log( μi ) = a + b log(xi) ※A = exp(a) と置いた Rの命令文 > glm( y ~ log(x) , family = Gamma(link = “log”), data = d ) glm()による推定では、平均μiを決める線形予測子とリンク関数だけを指定すれば良い平均・分散をshape, rateパラメータとどう対応付けるか考えなくて良い

Slide 50

Slide 50 text

参考文献

Slide 51

Slide 51 text

参考文献 - 30分だけでは決してよくわからないとてもとても難しい一般化線形モデル with R - 線形結合モデルを科学的説明たりうるか - Wikipedia 「線型結合」 - Wikipedia 「線型性」 - 確率密度関数と確率質量関数 - Wikipedia 「二項分布」