scikit-learnで線形単回帰

1.やること

線形回帰を行う。
説明変数が一つの場合、目的変数$y$と説明変数$x$の間に$y=ax+b$的な関係があるとして、$a$と$b$を求めたい。
説明変数が複数の場合はsklearn公式には以下のようにある。

\hat{y}(w, x) = w_0 + w_1 x_1 + ... + w_p x_p

データとして、2019年シーズンの北海道日本ハムファイターズの各試合データを用いる。データは一行＝一試合で、各試合のデータが記されている。
今回は「失点」を目的変数、「被安打」のみを説明変数として線形回帰モデルを構築する。

モデル式$y=ax+b$において、$b$は被安打0の時の失点、$a$は1本ヒットを打たれたときに増える失点数ということになりそうだ。

2.データ

まず、データ（csvファイル）をpd.DataFrameとして読み込む。

2-1

from pathlib import Path

import pandas as pd

p = Path("D:/data/")
df = pd.read_csv(p.joinpath("file.csv"))

df.head()

試合	失点	被安打
1	3	10
2	4	11
3	1	6
4	3	8
5	3	10

つづいて各列を説明変数・目的変数に代入。
ちなみにここで各変数の型をpd.DataFrame（行列的なやつ）ではなくpd.Series（ベクトル的なやつ）にしてしまうと後で面倒ということがわかった。

参考資料

print(type(df[["失点"]]))
print(type(df.loc[:,["失点"]]))
print(type(df.iloc[:, 0:1]))
# <class 'pandas.core.frame.DataFrame'>

print(type(df.loc[:, "失点"]))
print(type(df.iloc[:, 0]))
# <class 'pandas.core.series.Series'>

つまりこうする。1列のpd.DataFrame。

2-2

y = df.loc[:,["失点"]]
x = df.loc[:,["被安打"]]

3.基礎分析

プロットしてみる。

3-1

import matplotlib.pyplot as plt

plt.scatter(x, y)

x軸（横軸）が被安打、y軸（縦軸）が失点。
被安打が多いほど失点も多いという相関がありそう。

4.学習

まずモデルを設定。線形回帰はsklearn.linear_model.LinearRegression()らしい。

4-1

from sklearn import linear_model

model = linear_model.LinearRegression()

つづいてモデルをデータに合わせて学習させる。sklearnではmodel.fit(x, y)でこれを行う。

4-2

model.fit(x, y)

ここでxとyの型がpd.Series（ベクトル的なやつ）だとエラーがでる。行列計算するときにn行1列にすべきところ1行n列になるためと思われる。

4-2b

y = df.iloc[:, 0]
x = df.loc[:, "被安打"]

model.fit(x, y)

ValueError: Expected 2D array, got 1D array instead:　（以下略）

この場合はpd.Series.valuesでnumpy.ndarray型にし、np.reshape(-1, 1)でn行1列に変形するとうまくいった。

4-2c

y = df.iloc[:, 0]
x = df.loc[:, "被安打"]

model.fit(x.values.reshape(-1, 1), y.values.reshape(-1, 1))

# これもあり
model.fit(pd.DataFrame(x), pd.DataFrame(y))

5.結果

学習後、モデル式$ax+b$の、$a$はmodel.coef_、$b$はmodel.intercept_で取得できる。

5-1

print(model.coef_, model.intercept_)

[[0.8755288]] [-3.17832737]

（単回帰でも行列（切片はベクトル）で返ってくるのか……）
ともかくこれによると

失点 = 被安打×0.88 - 3.18

という式に表すことができるという。
被安打がわかったところでこんな計算してられないが、予測失点数はmodel.predict(x)で求めることができるので数式を立てる必要はない。

モデルの決定係数をmodel.score(x, y)で取得できる。
model.score(x, y)は「xからの予測値」と「y」間のR2を算出するもので、sklearn.metrics.r2_score(y, model.predict(x))を返す。引数を間違えそうになった。

5-2

print(model.score(x, y))

0.7048504590959697

ちなみに0.70≒0.84^2。かなり高い。

プロットしてみる。

5-3

lx = [[-100], [100]]

fig = plt.figure(figsize=(4, 4))
ax = fig.add_axes((0, 0, 1, 1), xlim=(-0.5, 23), ylim=(-0.5, 23))
ax.scatter(x, y, color="C0", zorder=10)
ax.plot(lx, model.predict(lx), color="C1", zorder=9)
[ax.spines[side].set_visible(False) for side in ["right", "top"]]
ax.grid(True, axis="both", color="gainsboro", alpha=0.8, zorder=7)

はい。

でもさ

失点はマイナスにならないという文句が聞こえてくる。
そういうときはモデル設定時にfit_intercept=Falseを与えると切片0モデルになる。

5b

model = linear_model.LinearRegression(fit_intercept=False)
model.fit(x, y)

print(model.coef_, model.intercept_, model.score(x, y))

fig = plt.figure(figsize=(4, 4), facecolor="white")
ax = fig.add_axes((0, 0, 1, 1), xlim=(-0.5, 23), ylim=(-0.5, 23))
ax.scatter(x, y, color="C0", zorder=10)
ax.plot(lx, model.predict(lx), color="C1", zorder=9)
[ax.spines[side].set_visible(False) for side in ["right", "top"]]
ax.grid(True, axis="both", color="gainsboro", alpha=0.8, zorder=7)

[[0.54919786]] 0.0 0.5906048285153787

結局

現実的な解釈。

5c

lx = [[(0-model.intercept_)/model.coef_], [100]]

fig = plt.figure(figsize=(4, 4), facecolor="white")
ax = fig.add_axes((0, 0, 1, 1), xlim=(-0.5, 23), ylim=(-0.5, 23))
ax.scatter(x, y, color="C0", zorder=10)
ax.plot([[0]]+lx, [[0]]+model.predict(lx).tolist(), color="C1", zorder=9)
[ax.spines[side].set_visible(False) for side in ["right", "top"]]
ax.grid(True, axis="both", color="gainsboro", alpha=0.8, zorder=7)

ちな

今回はGoogle Colaboratoryを使った。