Python入門トップページ

---

線形回帰分析をしてみよう-3 (sklearn)

問題

次の表はあるスーパーマーケットの売り場面積と売上高の実績値である．ここでは，線形回帰を sklearn.linear_model.LinearRegression を使って行い，新規開店予定の店舗 (P, Q) の売上高を予測しよう．

Python で線形回帰

まずは必要なモジュールをインポートする．

モジュールのインポートimport numpy as np
import pandas as pd
from sklearn import linear_model
import matplotlib.pyplot as plt

# 高解像度ディスプレイ用
from IPython.display import set_matplotlib_formats
# from matplotlib_inline.backend_inline import set_matplotlib_formats # バージョンによってはこちらを有効に
set_matplotlib_formats('retina')

データは GitHub のリポジトリで公開しているので，これを Pandas のデータフレームに直接読み込む．

問題のデータを読み込むurl = "https://github.com/rinsaka/sample-data-sets/blob/master/store-space-sales.csv?raw=true"
df = pd.read_csv(url)
df

Pandas のデータフレーム df から「space」列と「sales」列を取り出してそれぞれ NumPy 配列に変換する．このとき，x に該当する「space」は1次元配列ではなく，2次元配列（つまり行列形式）になるように変換する必要があることに注意する．

Pandas データフレームを NumPy 配列に変換x = df.loc[:, ['space']].values
y = df.loc[:, 'sales'].values
print(x)
print(y)

[[26]
 [28]
 [25]
 [26]
 [28]
 [33]
 [32]
 [29]
 [30]
 [29]
 [38]
 [40]
 [32]
 [26]
 [34]]
[65 65 62 59 69 73 79 71 73 71 86 88 71 63 80]

NumPy 配列の形式を確認する．

x の形式を確認するx.shape

(15, 1)

y の形式を確認するy.shape

(15,)

回帰分析の前に散布図を描いてみよう．

散布図fig, ax = plt.subplots(1, 1, figsize=(6, 4))
ax.scatter(x, y)
ax.set_xlabel('space')
ax.set_ylabel('sales')
ax.set_xlim(20,45)
ax.set_ylim(40,100)
plt.show()

データの準備ができたので，モデルを準備する．

モデルの準備model = linear_model.LinearRegression()

回帰分析は2次元のNumPy配列に変換した x, y を引数に与えるとよい．

回帰分析の実行results = model.fit(x, y)

結果を表示してみよう．

結果の表示results

LinearRegression()

係数を表示してみよう．

係数の表示model.coef_

array([1.83357349])

上の係数は array 形式になっているので，次のような方法で値を取得して表示してみよう．

係数の取得と表示a = model.coef_[0]
print(a)

1.8335734870317013

次は切片を表示してみよう．

切片の表示model.intercept_

15.926032660902955

切片も値を取得して表示してみよう．

切片の取得と表示b = model.intercept_
print(b)

15.926032660902955

決定係数 (R-squared) は score 関数で取得できる（もちろん，ここでの R-squared（係数，傾きなど）の結果と一致するはずである）．

決定係数の表示model.score(x, y)

0.9102297512020177

グラフを表示しよう．

グラフを描画fig, ax = plt.subplots(1, 1, figsize=(6, 4))
ax.scatter(x, y)
ax.plot(x, [a * t + b for t in x])
ax.set_xlabel('space')
ax.set_ylabel('sales')
ax.set_xlim(20,45)
ax.set_ylim(40,100)
ax.set_title('OSL results')
# plt.savefig('ols-2-2021.png', dpi=300, facecolor='white')
plt.show()

新規店舗の売上を予測する．

予測# 予測する
print('P (27) ==> ', a * 27 + b)
print('Q (38) ==> ', a * 38 + b)

P (27) ==>  65.4325168107589
Q (38) ==>  85.60182516810758

なお，sklearn では predict を使って予測ができる．

予測(1)model.predict([[27]])[0]

65.4325168107589

予測(2)model.predict([[38]])[0]

85.60182516810761

さらに，x が 27 と 38 という2つの新規店舗の売上高がそれぞれ 61 と 79 であることがわかっている状況で，これらのデータをテストデータとして，テストデータに対する当てはまりの良さを検証する．まず，NumPy の行列形式でテストデータを準備する．

テストデータの準備x_test = np.array([[27], [38]])
y_test = np.array([61, 79])
print(x_test)
print(y_test)

[[27]
 [38]]
[61 79]

テストデータを NumPy 配列で準備することで，予測もまとめて行うことができる．

y_pred = model.predict(x_test)
print(y_pred)

[65.43251681 85.60182517]

当てはまりの良さは決定係数で検証できる．これは model.score() 関数にテストデータを与えるだけで良い．

スコア：決定係数model.score(x_test, y_test)

0.6096833282227203

なお，この決定係数の計算方法はここでは次のように記載されている． \begin{eqnarray*} R^2 = 1 - \frac{u}{v} \end{eqnarray*} ここで \(u\) は残差の二乗和，\(v\) は標本値と標本平均値の偏差の二乗和である．これらを順に計算して決定係数が正しいことも確認しておこう．

まず，NumPy 配列形式のデータでは要素ごとの演算を簡単に行えるので，残差は次のように求めることができる．

残差y_test - y_pred

array([-4.43251681, -6.60182517])

残差を二乗する．

残差の二乗(y_test - y_pred)**2

array([19.64720528, 43.58409555])

これらの総和を求めるためには NumPy の sum() 関数が利用できる．これで，残差の二乗和 \(u\) を計算できる．

残差二乗和u = ((y_test - y_pred)**2).sum()
print(u)

63.2313008279193

次に，標本平均値は次のように y_test.mean() で求めることができる．

print(y_test)
print(y_test.mean())

[61 79]
70.0

平均値との偏差を求める．

print(y_test)
print(y_test.mean())
print(y_test - y_test.mean())

[61 79]
70.0
[-9.  9.]

偏差の二乗を計算する．

print((y_test - y_test.mean())**2)

[81. 81.]

偏差の二乗和 \(u\) は次のように計算できる．

偏差の二乗和v = ((y_test - y_test.mean())**2).sum()
print(v)

162.0

上の手順をまとめると，決定係数は次のように計算できる．

決定係数u = ((y_test - y_pred)**2).sum()
v = ((y_test - y_test.mean())**2).sum()
r2 = 1 - u/v
print(r2)

0.6096833282227203

したがって，score()で計算した値と一致することが確認できた．

スコア（決定係数）の計算：再掲model.score(x_test, y_test)

0.6096833282227203

---