My Notes

統計学とかR(R言語)とかPython3の覚え書きとか走り書きとか。 座右の銘にしたい: All work and no play makes Jack a dull boy.

Python3(pandas, Scipy, StatsModels)で、要約統計量、相関係数、p値、単回帰分析。statsmodels.formula.api.ols(), statsmodels.formula.api as smf, smf.ols()

Python3コード

#!/usr/bin/env python3


"""(docstring)
"""


# 使用したデータは『マンガでわかる統計学 [回帰分析編]』第2章


import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import scipy as sp
from scipy import stats
import statsmodels.formula.api as smf
# statsmodels.formula.api.ols()
#
# patsyが内部で使われてるようなので、動かない場合、
# patsyがinstallされてるかどうか確認。
# statsmodelsが入っていれば、patsyも入っているはずだが。


def main():
    """(docstring)
    """

    # matplotlib
    # macOSやOS Xで文字化けするなら。
    # (Win機ならフォント名をWin用のに書けばいけるかもしれない(試してない))。
    font = {'family' : 'Osaka'}

    dat = {'最高気温': [29, 28, 34, 31, 25,
                      29, 32, 31, 24, 33,
                      25, 31, 26, 30],
           'アイスティー': [77, 62, 93, 84, 59,
                          64, 80, 75, 58, 91,
                          51, 73, 65, 84]}

    df = pd.DataFrame(dat, columns = ['最高気温', 'アイスティー'])
    print('df')
    print('')
    print(df)
    print('')

    # 要約統計量(pandas)
    print(df.describe())
    print('')

    # x, yとしたいならば。
    # x = df.ix[:, '最高気温']
    # y = df.ix[:, 'アイスティー']

    # matplotlib
    # 散布図。
    plt.scatter(df.ix[:, '最高気温'], df.ix[:, 'アイスティー'])
    plt.grid()
    plt.xlabel('最高気温')
    plt.ylabel('アイスティー')
    plt.show()

    # 相関係数(ピアソン)とp値。
    print('scipy.stats.pearsonr()')
    r, p = sp.stats.pearsonr(df.ix[:, '最高気温'], df.ix[:, 'アイスティー'])
    print('{0} {1}{2}    {3} {4}'.format('相関係数r', r, ',', 'p値', p))
    print('')

    # statsmodels.formula.api.ols()
    # statsmodels.formula.api as smf
    # smf.ols() として使用。
    # (UserWarning: kurtosistest ... とりあえずこの警告は無視してもいいと思う)。
    #
    # 基本的には、R(R言語)と似ている。
    # (似せて作ったんだろうから当たり前だろうけど)。
    mod = smf.ols(formula="df.ix[:, 'アイスティー'] ~ df.ix[:, '最高気温']", data=df)
    res = mod.fit()
    print(res)
    print('')

    print(res.summary())
    print('')

    print(res.params)
    print('')

    print(res.predict())
    print('')

    # 散布図と回帰直線
    plt.scatter(df.ix[:, '最高気温'], df.ix[:, 'アイスティー'])
    # res.predict()には引数不要。
    plt.plot(df.ix[:, '最高気温'], res.predict(), color='red')
    # 回帰直線は自分で計算式を書いても描ける。
    # b, a = res.params # 切片と回帰係数
    # print(b, a)
    # 第1引数にx(ここでは、xは具体的にはすべてdf.ix[:, '最高気温']),
    # 第2引数に y = a*x+b の計算式部分。
    # plt.plot(df.ix[:, '最高気温'], a*df.ix[:, '最高気温']+b, color='red')
    # y = a*x+b
    # で代入して、yを第2引数にしてもいい。
    # plt.plot(df.ix[:, '最高気温'], y, color='red')
    plt.grid()
    plt.xlabel('最高気温')
    plt.ylabel('アイスティーの注文数')
    plt.show()


if __name__ == '__main__':
    main()


# Fitting models using R-style formulas
# statsmodels.formula.api as smf
# のols()を使ったほうが楽だった。

出力

df

    最高気温  アイスティー
0     29      77
1     28      62
2     34      93
3     31      84
4     25      59
5     29      64
6     32      80
7     31      75
8     24      58
9     33      91
10    25      51
11    31      73
12    26      65
13    30      84

            最高気温     アイスティー
count  14.000000  14.000000
mean   29.142857  72.571429
std     3.158801  13.019006
min    24.000000  51.000000
25%    26.500000  62.500000
50%    29.500000  74.000000
75%    31.000000  83.000000
max    34.000000  93.000000

scipy.stats.pearsonr()
相関係数r 0.9069229780508896,    p値 7.66141280445014e-06

<statsmodels.regression.linear_model.RegressionResultsWrapper object at 0x10c08d4e0>

                            OLS Regression Results                            
==============================================================================
Dep. Variable:     df.ix[:, 'アイスティー']   R-squared:                       0.823
Model:                            OLS   Adj. R-squared:                  0.808
Method:                 Least Squares   F-statistic:                     55.61
Date:                Sun, 13 Aug 2017   Prob (F-statistic):           7.66e-06
Time:                        13:34:54   Log-Likelihood:                -43.174
No. Observations:                  14   AIC:                             90.35
Df Residuals:                      12   BIC:                             91.63
Df Model:                           1                                         
Covariance Type:            nonrobust                                         
====================================================================================
                       coef    std err          t      P>|t|      [0.025      0.975]
------------------------------------------------------------------------------------
Intercept          -36.3612     14.687     -2.476      0.029     -68.362      -4.360
df.ix[:, '最高気温']     3.7379      0.501      7.457      0.000       2.646       4.830
==============================================================================
Omnibus:                        5.158   Durbin-Watson:                   1.669
Prob(Omnibus):                  0.076   Jarque-Bera (JB):                1.433
Skew:                          -0.180   Prob(JB):                        0.488
Kurtosis:                       1.474   Cond. No.                         282.
==============================================================================

Warnings:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

Intercept          -36.361233
df.ix[:, '最高気温']     3.737885
dtype: float64

[ 72.03744493  68.29955947  90.72687225  79.51321586  57.08590308
  72.03744493  83.25110132  79.51321586  53.34801762  86.98898678
  57.08590308  79.51321586  60.82378855  75.7753304 ]

散布図と回帰直線のスクリーンショット

f:id:my_notes:20170813133802p:plain

f:id:my_notes:20170813133823p:plain

参考文献

マンガでわかる統計学 回帰分析編

マンガでわかる統計学 回帰分析編