Pythonで回帰分析をちょっとやってみましたのでメモです。
↓このまんまです。
http://www.atmarkit.co.jp/ait/articles/1311/05/news003.html
import pandas as pd
import numpy as np
df = pd.read_csv('children_data2005_08_130819.csv', header=-1, skiprows=2, encoding='Shift_JIS')
#データが欠けている行は削除する。
df = df[np.isfinite(df[9])]
df = df[np.isfinite(df[10])]
#身長と体重のデータがある10列目と11列目のみを抽出する
height = df[9]
weight = df[10]
#プロットしてみる
plot(height, weight, 'bo')
Out[8]: []
#単回帰分析を用いてこのモデルの妥当性を検討する。
#pandasにはols関数が用意されていて
#ols関数は
# 目的変数 y
# 説明変数 x
# 直線が原点(0,0)を通るか切片を持つかどうかのパラメータ(intercept)
#を指定すると、最小二乗法(Ordinary Least Squares)を用いてモデルの計算を行う。
model = pd.ols(y=weight, x=height, intercept=True)
plot(model.x['x'], model.y_fitted, 'g-')
Out[17]: []
#modelの中身のポイントを確認してみる
model
Out[20]:
-------------------------Summary of Regression Analysis-------------------------
Formula: Y ~ +
Number of Observations: 1602
Number of Degrees of Freedom: 2
R-squared: 0.8445
Adj R-squared: 0.8444
Rmse: 2.1957
F-stat (1, 1600): 8691.3653, p-value: 0.0000
Degrees of Freedom: model 1, resid 1600
-----------------------Summary of Estimated Coefficients------------------------
Variable Coef Std Err t-stat p-value CI 2.5% CI 97.5%
--------------------------------------------------------------------------------
x 0.0336 0.0004 93.23 0.0000 0.0329 0.0343
intercept -17.2343 0.3742 -46.06 0.0000 -17.9677 -16.5009
---------------------------------End of Summary---------------------------------
このモデルの数式
[体重]=[身長]×0.0336-17.2343
このモデルが統計的にどのくらい正しいのか検証する。
□決定係数(R-squared: 0.8445)
説明変数が目的変数のどれくらいを説明できるかを表す値で、0~1の間の値を取り、寄与率と呼ばれることもあります。
今回のモデルでは0.8445となっており、このモデルで84%以上説明できているということができます。
□F検定(F-stat (1, 1600): 8691.3653, p-value: 0.0000)
F検定(F-test)と呼ばれるモデル全体の妥当性を検討する際に使う値です。F値の有意確率(p-value)を判断基準として用います。
今回のモデルの場合は、F値(8691.3653)が十分に大きく、有意確率(0.0000)も十分に小さく0.01以下であるため99%以上の確率で妥当であると言えます。
有意確率の比較は、99%の場合には0.01ですが、95%の確率で検定する場合には、0.05と比較して妥当性を検証します。
□t検定
モデル全体の妥当性をF検定で判断した後は、それぞれのパラメータ(係数aと定数b)の妥当性を検証します。
係数a(x):t-stat=93.23、p-value=0.0000
定数b(intercept):t-stat=-46.06、p-value=0.0000
このモデルの場合、xの有意確率もinterceptの有意確率も0.01よりも十分に小さい値なので、両方の値は妥当であると判断できます。
■PVとセッション数の関係
せっかくなので上のやり方を真似してあるサイトのPVとセッションの関係を
回帰分析してみます。
import pandas as pd
import numpy as np
%pylab
df = pd.read_csv('data.csv', names=['PV','SESSION'], header=-1, encoding='Shift_JIS')
df = df.replace(',','', regex=True)
## 数字に「,」がはいっているので取り除いています。
pv = df['PV'].astype(float)
se = df['SESSION'].astype(float)
plot(pv, se, 'bo')
model = pd.ols(y=se, x=pv, intercept=True)
plot(model.x['x'], model.y_fitted, 'g-')
model
Out[32]:
-------------------------Summary of Regression Analysis-------------------------
Formula: Y ~ +
Number of Observations: 61
Number of Degrees of Freedom: 2
R-squared: 0.2639
Adj R-squared: 0.2515
Rmse: 16730.9688
F-stat (1, 59): 21.1571, p-value: 0.0000
Degrees of Freedom: model 1, resid 59
-----------------------Summary of Estimated Coefficients------------------------
Variable Coef Std Err t-stat p-value CI 2.5% CI 97.5%
--------------------------------------------------------------------------------
x 0.0319 0.0069 4.60 0.0000 0.0183 0.0455
intercept 192177.6432 10286.0067 18.68 0.0000 172017.0700 212338.2163
---------------------------------End of Summary---------------------------------
このモデルの数式
[SESSION]=[PV]× 0.0319+192177.6432
このモデルが統計的にどのくらい正しいのか検証する。
□決定係数(R-squared: 0.2639)
今回のモデルでは0.2639となっており、このモデルで26%以上説明できているということでしょうか。
□F検定(F-stat (1, 59): 21.1571, p-value: 0.0000)
今回のモデルは、F値( 21.1571)が小さいですね。
なので、このモデルはだどうではなさそうですね。もっとサンプルがあれば変わってくるのでしょうか。
ということで、あんまりPVとセッション数の関係は妥当性のないモデルのようです。
まあでもいろいろと勉強になりました。
こいつはいろんな場面で使えそうです。
数字のコンマ(,)を取り除く方法
df = df.replace(',','', regex=True)
以上