いろいろするblog

プログラミングとか統計とかえーあい的なあれこれとか

データ解析のための統計モデリング入門をpythonでやる（3）-一般化線形モデル-

【あらすじ】
種子データを一般化線形モデルであらわしてみた。

【本文】
説明変数を組み込んだ統計モデルを扱う。
つまり、個々のデータが独自の説明変数の値を持ち、それによって従うモデルが変わるようになる。
このようなモデルを一般化線形モデルという。

今回はある植物のデータとして
種子数:y_i
体サイズ:x_i
施肥処理:C or T
であるデータを扱う。

前回のポアソン分布

\frac{\lambda^{k} exp(-\lambda)}{k!}

において、\lambdaが体サイズx_iによって変わることを考える。
今回は以下の式に従うと考える。

\lambda_i = exp(\beta_1 + \beta_x x_i)

ここで、logを取ると

\log \lambda_i = \beta_1 + \beta_x x_i

となる。
右辺のように\betaの線形結合になるものを線形予測子と呼ぶ。
このように
(\lambda_i の関数) = (線形予測子)
となる場合、左辺の\lambda_i の関数をリンク関数と呼ぶ。
今回は対数を取っているので対数リンク関数である。

平均がこの場合の対数尤度を求める。
つまり、

\log_10 L(\beta_1 , \beta_2) = \sum{(y_i \log \lambda - \lambda - \sum_{k}^{y_i} \log k)}

が最大となる \beta_1 , \beta_2 を求める。
pythonでざっくり当てはめるにはstatsmodelを使う。

[python のコード]

model = sm.formula.api.glm(formula="y ~ x", data=data, family=sm.api.families.Poisson(link=sm.api.families.links.log))

formula="y ~ x" は yの説明変数はx であることを示している。
もしy の説明変数をfと考える場合（つまり、体サイズx は関係なく、施肥の有無によってのみ種子数がかわる）は"y ~ f"とする。
また、y の説明変数をxとfの両方と考える場合（体サイズと施肥の両方が関係する）は"y ~ x + f"とする

yに対する説明変数がxのみの場合、一般化線形モデルではxの値の位置に依存した異なる分布でyの値が表現される。

データの特徴を反映した確率分布とリンク関数を使いましょう。