Blog

Pythonで「データ解析のための統計モデリング入門」:第3章

Pythonで緑本の演習をしていきます。前記事はこちら

第3章、説明変数をくみこんだ統計モデルについてです。

 

第3章 一般化線形モデル(GLM)

第2章で使用したサンプルデータは、架空の植物50個体からなる集団を調査して得られた各個体の種子数を数えたものでした。このデータで、個体ごとに異なる説明変数(個体の属性)によって平均種子数が変化する統計モデルを観測データにあてはめることを、一般化線形モデル(GLM)と言います。

第3章では、一般化線形モデルの中のポアソン回帰を実施します。

 

3.1 例題:個体ごとに平均種子数が異なる場合

第3章の例題では、架空の植物100個体を調査して得られた、個体ごとの種子数のデータを用います。植物個体iの種子数は \( y_{i} \) 個であり、個体の属性のひとつである体サイズ \( x_{i} \) が観測されているとします。また、全個体のうち50個体は何も処理をしておらず、残り50個体には肥料を加える処理(処理T)をほどこしているとします。

 

3.2 観測されたデータの概要を調べる

サンプルデータ(csvファイル)はこちらからダウンロードできます。pandasのread_csvを使ってDataFrame形式で読み込みます。

dの列ごとにデータを表示させてみます。

pandasのinfo()メソッドで各列のデータの型について知ることができます。

describe()メソッドでデータの概要を示します。

value_counts()メソッドで \( f )\ の要素の数を得られます。

 

3.3 統計モデリングの前にデータを図示する

まずは観測データに余計な手を加えないで、データ全体を見てみます。

横軸に \( x \) 列、縦軸に \( y \) 列をとった散布図をプロットします。

pandasのプロット機能を使っても良いかと思ったのですが(簡単に書ける場合が多い)、後々複雑なプロットをしたいと思った時を考えて、はじめからmatplotlibでいきます。

scatter

次に、横軸を施肥処理、縦軸に \( y \) 列をとった箱ひげ図をプロットします。

box_plot

 

3.4 ポアソン回帰の統計モデル

個体ごとの平均種子数 \( \lambda_{i} \) が体サイズ \( x \) や施肥処理 \( f \) に影響されるようなモデルを設計します。

まず、個体 \( i \) の体サイズ \( x_{i} \) だけに依存する統計モデルについて考えます。切片 \( \beta_{1} \) 、 傾き \( \beta_{2} \) の関係を図示すると次のようになります。

test

 

3.4.2 あてはめとあてはまりの良さ

ポアソン回帰とは、観測データに対するポアソン分布を使った統計モデルのあてはめであり、この統計モデルの対数尤度 \( \log L \) が最大になるパラメータ \( \hat{\beta}_{1} \) と \( \hat{\beta}_{2} \) の推定値を決めることです。

 Rではglm()関数を使うことで手軽にGLMのあてはめができます。pythonでやる方法をググッてみたら、scikit-learnやstatsmodelsというモジュールが有名なようです。statsmodelsはRライクに使えるとのことで、ここではstatsmodelでやってみます。

調べてみると、Statsmodelsには通常の関数型呼び出しに対応するAPIと、R-styleのモデル式を用いる数式対応API(formula API)の2種類が実装されているとのことでした。

Statsmodels formula API と多項式回帰モデル

以下のようにして切片 \( \hat{\beta}_{1} \) と \( \hat{\beta}_{2} \) の推定値を求めてみました。GLMはformula APIで実行し、ポアソン分布の指定に通常のAPIを使っています。Rのように、formulaの指定を’従属変数 ~ 説明変数’という風に指定できました。

 

 3.4.3 ポアソン回帰モデルによる予測

次に、このポアソン回帰の推定結果を使って、様々な体サイズ \( x \) における平均種子数 \( \lambda \) の予測を行います。

 \[ \lambda = \exp (1.29 + 0.0757 x) \]

これを図示します。

predict

 

3.5 説明変数が因子型の統計モデル

施肥効果 \( f_{i} \) を説明変数としてくみこんだモデルも検討します。施肥処理は因子型、カテゴリ型のデータとして格納されています。

Pythonで推定を行います。因子型の説明変数であっても、簡単にモデル式を指定できます。

 

3.6 説明変数が数量型 + 因子型の統計モデル

個体の体サイズ \( x_{i} \) と施肥効果 \( f_{i} \) の複数の説明変数を同時にくみこんだ統計モデルを作ってみます。

statsmodelでは、モデル式の部分をx + fとするだけで適切に処理してくれます。

 


This Post Has 0 Comments

Leave A Reply