hidekatsu-izuno 日々の記録

プログラミング、経済政策など伊津野英克が興味あることについて適当に語ります(旧サイト:A.R.N [日記])

Python StatsModels メモ

R と比較すると微妙にサポートされていない機能があって困ることが多い StatsModels ですが、Python に寄せていきたいので、できるだけ使ってみてます。

ライブラリのロード

import statsmodels.api as sm

# R互換の関数方式を使う場合はこっち
import statsmodels.formula.api as smf

データのロード

import pandas as pd
data = pd.read_csv('data.csv')
# TSVは、data = pd.read_csv('data.tsv', delimiter='\t')
data = data.fillna(0) # 欠損値がある場合は 0 で埋める

モデルの定義

単回帰

model = sm.OLS(data.y, sm.add_constant(data.x))
#  あるいは、 model = smf.ols('y ~ x', data)
results = model.fit()
results.summary()

重回帰

model = sm.OLS(data.y, sm.add_constant(data[ ['x1', 'x2', 'x3'] ]))
#  あるいは、 model = smf.ols('y ~ x1 + x2 + x3', data)
results = model.fit()
results.summary()

ロバスト回帰

model = sm.RLM(data.y, sm.add_constant(data.x))
#  あるいは、 model = smf.rlm('y ~ x', data)
results = model.fit()
results.summary()

分位点回帰

from statsmodels.regression.quantile_regression import QuantReg
model = QuantReg(data.y, sm.add_constant(data.x))
#  あるいは、 model = smf.quantreg('y ~ x', data)
results = model.fit(q = 0.5) # 50%分位点
results.summary()

リッジ回帰、ラッソ回帰、Elastic Net

model = sm.OLS(data.y, sm.add_constant(data.x))
results = model.fit_regularized(L1_wt=0) # L1_wt が 0 の時 リッジ、1の時ラッソ、その間が Elastic Net
results.params

OLSと同じだが、fit の代わりに fit_regularized を使う。ただ、fit_regularized だと、なぜか summary() が None を返すので、params の値を参照する必要がある。

一般化線形モデル (GLM)

model = sm.GLM(data.y, data.x, family=sm.families.Gamma())
#  あるいは、 model = smf.glm('y ~ x', data, family=sm.families.Gamma())
results = model.fit()
results.summary()

一般化線形混合モデル (GLMM)

model = sm.MixedLM(data.y, sm.add_constant(data.x), groups=data.g)
#  あるいは、 model = smf.mixedlm('y ~ x', data, groups=data.g)
results = model.fit()
results.summary()

グラフの表示

import matplotlib.pyplot as plt
plt.scatter(data.x, data.y) # 散布図
plt.plot(range(int(data.x.min()), int(data.x.max())), results.predict(range(int(data.x.min()), int(data.x.max())))) # 折れ線グラフ

# 図表をクリアする
plt.clf()