学生による学生のためのデータサイエンス勉強会

【さきどりPython#9】非線形回帰モデルを作ってみよう

https://bdarc.net/wp-content/uploads/2020/05/icon-1.pnghistoroid

こんにちは。清野(@historoid1)です。

kaggleデータで線形回帰モデルをつくる話をしましたが、ちょっと書き加えたいところがあったので、まとめておきます。

今回は非線形回帰がテーマです。

前回の記事(以下に示します)では、kaggleデータで線形回帰モデルを動かしました。

【さきどりPython#9】kaggleのデータで線形回帰モデルをつくろう

そのモデルは単回帰(説明変数が1つ)で、線形回帰(一次関数)でした。

今回は重回帰(変数が2つ以上)モデルや、非線形(二次以上の関数)モデルを作ってみたいと思います。

非線形回帰モデルを作ってみよう

scikit-learnの線形モデルの回帰関数

scikit-learnの線形回帰モデルはデフォルトでは2つの変数をとることができます。つまり重回帰ですね。

この重回帰関数を式で表すと以下のようになります。

$$\hat{y} = \hat{\beta_0} + \hat{\beta_1}x_1 + \hat{\beta_2}x_2 + \varepsilon$$

\(\hat{y}\)は推測値、\(\hat{\beta_0}\)は説明変数\(x_1\)に対する係数(重み)、\(\hat{\beta_2}\)は説明変数\(x_2\)に対する係数(重み)です。

なお\(\beta_0\)はバイアス(切片)、\(\varepsilon\)は誤差といわれます。

前回の記事では、実はむりやり変数を1つにして実装していました。この式でいえば\(x_2\)の項がないということです。

前回のおさらい

使用するデータ

kaggle

kaggleのアボカド価格予測のデータです。

目的変数Total Volumeで、説明変数Total Bagsでした。

データの分布

この2つをプロットすると以下の散布図が得られます。

線形回帰モデルの実装と可視化

そして線形回帰モデルを実装して、線形回帰直線を得ることが出来ました。

# ライブラリをロード
from sklearn.linear_model import LinearRegression

# 線形回帰器を作成
regression = LinearRegression()

# データ
feature = df[['Total Bags']].values
target = df['Total Volume'].values

# 線形回帰器を訓練
model = regression.fit(feature, target)

ここまでが前回までの大まかな流れです。

説明変数を2つにして実装する

教科書ではもともと2変数だったので、それにならって2つに増やしてみたいと思います。

説明変数を2つにして実装する

今回は4770列を説明変数に加えます。4770の数だけ列を加えるのではなくて、そういう名前の列ですよ。

4770列とTotal Volume列の散布図
Total Bags4770列の散布図

2つのグラフは似てるけど別物ですからね

では以上の3変数を可視化しましょう。

from mpl_toolkits.mplot3d import Axes3D
# seaborn使ってみます。
import seaborn as sns

fig = plt.figure(figsize=(10,10))
ax = Axes3D(fig)
ax.scatter3D(df['Total Bags'], df['4770'], t)

ax.set_xlabel('Total Bags')
ax.set_ylabel('4770')
ax.set_zlabel('Total Volume')

plt.show()

うーん……。なんかよくわかりませんが、Vの字に分かれてるんですかね?たぶん。

重みとバイアスの確認

# 線形回帰器を作成
lr = LinearRegression()

# データ
f = df[['Total Bags', '4770']]
t = df['Total Volume']

# 線形回帰器を訓練
mdl = lr.fit(f, t)


# 係数(傾き、重み)の確認
mdl.coef_
# array([2.55865479, 9.42474542])

# バイアスの確認
mdl.intercept_
# 22231.323813074036

重みは配列で返ってきましたが、これは\(\hat{\beta_1}\)と\(\hat{\beta_2}\)に相当します。

つまりTotal BagsTotal Volumeの線形回帰直線の傾きが\(2.55865479\)で、4770Total Volumeでの傾きが\(9.42474542\)というわけです。

3変数による3次元空間を線形回帰するので、平面で分離することになります。

いかがでしょうか。これでわかりましたね。

しかも結構きれいに線形回帰されてる感じですね。

【応用編】説明変数同士の影響を考慮する

交互作用 interaction とは

交互作用は2つの因子が組み合わさることで初めて現れる相乗効果のことです。「肥料の量×土の種類」の場合、肥料の量と土の種類が相互に影響を及ぼし合っていることを表します。また、交互作用による効果のことを「交互作用効果」といいます。

統計WEBより引用

というわけで、交互作用項というものを設定します。

$$\hat{y} = \hat{\beta_0} + \hat{\beta_1}x_1 + \hat{\beta_2}x_2 + \hat{\beta_3}x_1x_2 + \varepsilon$$

\(\hat{\beta_3}x_1x_2\)が交互作用項です。

交互作用項を設定する

交互作用があるかは不明ですが、今回の例でも試してみましょう。

from sklearn.linear_model import LinearRegression
# 新しいライブラリ
from sklearn.preprocessing import PolynomialFeatures

# データ
f = df[['Total Bags', '4770']]
t = df['Total Volume']

# 交互作用項をつくる
polynomial = PolynomialFeatures(degree=2, interaction_only=True)
features_polynomial = polynomial.fit_transform(f)

# 線形回帰器を作成
reg = LinearRegression()

# 線形回帰器を訓練
interact_model = reg.fit(features_polynomial, t)

degree=2では、2つの変数の組合せで交互作用項を作ることを指定しています。

またinteraction_only=Trueでは、多項式特徴量(\(x_1^2\)のようなべき乗の項)を作らないと指定しています。

interact = np.multiply(df['Total Bags'], df['4770'])

上記のように自分で交互作用項を作っても構いません。

つぎに、polynomial.fit_transform(f)で\(\hat{y} = \hat{\beta_0} + \hat{\beta_1}x_1 + \hat{\beta_2}x_2 + \hat{\beta_3}x_1x_2 + \varepsilon\)を作っています。

重みの確認

interact_model.coef_
# array([ 0.00000000e+00,  2.65331176e+00,  1.05499410e+01, -2.36299414e-07])

これで\(\hat{\beta_0}, \hat{\beta_1}, \hat{\beta_2}, \hat{\beta_3}\)が返されました。

この式は4次元空間なので可視化はできません

非線形回帰を行う

非線形回帰関数

さて、いよいよ非線形回帰を行います。

まずは式を確認しましょう。

$$\hat{y} = \hat{\beta_0} + \hat{\beta_1}x_1 + \hat{\beta_2}x_1^2 + \dots + \hat{\beta_d}x_i^d$$

このようになっています。

ちょっとわかりにくいのですが、例えば説明変数が3個、次数が2なら以下のようになります。

$$\hat{y} = \hat{\beta_0} + (\hat{\beta_1}x_1 + \hat{\beta_2}x_1^2) + (\hat{\beta_1}x_2 + \hat{\beta_2}x_2^2) + (\hat{\beta_1}x_3 + \hat{\beta_2}x_3^2) + \varepsilon $$

非線形回帰関数の実装

上記のような回帰式だと複雑なので、変数1つに戻ってみます。

もともとはこんな感じの回帰直線でした。

これを3次式で表してみたいと思います。

$$\hat{y} = \hat{\beta_0} + \hat{\beta_1}x_1 + \hat{\beta_2}x_1^2 + \hat{\beta_3}x_1^3 + \varepsilon$$

from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures

# データ
f = df['Total Bags'].values.reshape(-1, 1)
t = df['Total Volume'].values
dimension = 3

# 多項式をつくる
polynomial = PolynomialFeatures(degree=dimension)
features_polynomial = polynomial.fit_transform(f)

# 線形回帰器を作成
lrp = LinearRegression()

# 線形回帰器を訓練
lrp_model = lrp.fit(features_polynomial, t)

# 回帰式
res = []
for x in f:
    tmp = 0
    for i in range(dimension):
        tmp += lrp_model.coef_[i]*x**i
    res.append(tmp)

# 散布図を描く
fig = plt.figure(figsize=(10, 5))
ax = plt.axes()
plt.plot(f, t, 'o', color='blue')

# 回帰関数をプロット
plt.plot(f, res,'o', color='orange')

なんかあんまり一次式と変わりませんでしたね。

4次式にすると以下のようになります。

誤差が大きすぎて、全体が潰れて見えています。

つまりこのデータでは1次から3次程度が良い次数であることが分かります。


https://bdarc.net/wp-content/uploads/2020/05/icon-1.pnghistoroid

今回はここまでです。

あまり良い例を見せることが出来ませんでしたが、自分で適したデータを見つけて試してみてください。

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です