【機械学習】単回帰分析をわかりやすく解説|python
単回帰分析
教師あり学習
今回は、教師あり学習の基礎中の基礎である「単回帰分析」を実装します。
教師あり学習とは、説明変数（インプット）から目的変数（アウトプット）を予測するモデルを求める手法です。
訓練データには目的変数や説明変数があり、あらかじめ作ったモデルに訓練データの説明変数を入力し、そのモデルからの出力が訓練データの目的変数に近づくようにモデルのパラメータを調整することで学習していきます。
さて、単回帰分析の説明の前に、構造を理解しておきましょう。
まず、教師あり学習は、目的変数（予測したい変数）のデータ形式によって主に2つのカテゴリに分類されます
- 回帰（regression）: 目的変数が連続的な数値を取る場合（例：株価、身長、気温）
- 分類（classification）: 目的変数が離散的なカテゴリを取る場合（例：「男性・女性」、「幼児・小学生・学生・大人」）
教師あり学習には様々なアルゴリズムがあり、それぞれ回帰問題や分類問題、あるいはその両方に適用できます。
代表的なものには以下があります。
- 重回帰分析（multiple linear regression）
- ロジスティック回帰（logistic regression）
- k近傍法（k-Nearest Neighbors）
- 決定木（Decision Tree）
単回帰分析は、中でも最も基本的な回帰分析の手法で、1つの説明変数（独立変数）を用いて1つの目的変数（従属変数）を予測するモデルです。見たことがありますかね？
$${y = \beta_0 + \beta_1x + \epsilon}$$
- ${y}$ は目的変数
- ${x}$ は説明変数
- ${\beta_0}$ は切片（定数項）
- ${\beta_1}$ は傾き（回帰係数）
- ${\epsilon}$ は誤差項
${\beta_0}$ と ${\beta_1}$ はモデルのパラメーターと呼ばれ、データから推定されます。
上記の通り、重みをもとに目的変数に対する説明変数の寄与を考えられる+線形モデルのため、根拠が明確で扱いやすいなどの利点があります。
説明変数を2個以上にすると重回帰分析になり、目的変数に対して非線形活性化関数(例えばシグモイド関数)を加えると単純パーセプトロン、何層も非活性化関数を重ねたものが多層パーセプトロン(深層学習)になります。
補足｜相関係数と回帰係数は何が違うのか？
ここで少し寄り道です。
特に単回帰分析においては、1つの説明変数と目的変数の関係を図るという点で相関係数との違いをうまく説明できないのではないのでしょうか？
結論違うのですが、定数項と誤差項がない線形モデルに限り適用できます。多変量回帰では単純には適用できません。
実は、誤差項と定数項を除外し、回帰係数を求めた後、標準偏差でのスケーリングを行うと相関係数と一致します。
ここで、それぞれの定義を今一度見てみましょう。
相関係数（Pearsonの相関係数）
$${r = \frac{\sum_{i=1}^{n}(x_i – \bar{x})(y_i – \bar{y})}{\sqrt{\sum_{i=1}^{n}(x_i – \bar{x})^2}\sqrt{\sum_{i=1}^{n}(y_i – \bar{y})^2}}}$$
単回帰の回帰係数（定数項と誤差項がない場合）
$${\beta = \frac{\sum_{i=1}^{n}x_i y_i}{\sum_{i=1}^{n}x_i^2}}$$
目的変数と説明変数を、それぞれの標準偏差で割ることで、相関係数と回帰係数が一致します。
標準化すると、$$\beta = r$$となります。
では、なぜ一致するのでしょうか？
標準化によって、変数の尺度の影響を取り除くことができるからですね。
標準化すると、${x_i = \frac{x_i – \bar{x}}{s_x}}$、${y_i = \frac{y_i – \bar{y}}{s_y}}$となります。
標準化した回帰係数を計算すると、相関係数の式と同一になりました。
こんな感じで切り口を変えて勉強してみると意外な発見があると思います。
相関係数や標準偏差については、こちらのコンテンツがおすすめです。
ベクトルでの表現と損失関数について
では、線形代数を使って単回帰分析の数学的背景をご説明します。
$$f(x)=w_{1}x+w_{0}$$
単回帰では上のような直線を想定します。
\(x\)については特徴ベクトル、\(w\)については重みベクトルなので、以下のように整理できます。
$$\boldsymbol{x}=[x,1]^T$$
$$\boldsymbol{w}=[w_{1},w_{0}]^T$$
よって上のようなベクトルの内積の形で直線を表すことができます。
重回帰分析なら高次元での平面または超平面ですね。
ちなみに\(R^n\) において次元がちょうど 1 だけ小さい \(n-1\) 次元 の空間のことを超平面 (hyperplane) といいます。
$$\begin{equation} \mathbf{y} = \mathbf{X}\mathbf{w} + \boldsymbol{\epsilon} \end{equation}$$
さて、線型回帰分析では予測値\(\hat{y}\)と教師データ\(y\)の二乗誤差を最小化するような重みベクトル\(\boldsymbol{w}\)を探したいです。この損失関数を残差平方和と呼びます。
$$\sum_{i}^N||y_{i}-\hat{y}_{i}||^2$$
手法によって損失関数の種類は違ってくるし、目的によっても変わってきます。
RSSは全ての予測に対する誤差の合計を示すため、データの数が多いほど大きな値になる傾向があります。
平均二乗誤差と残差平方和は違うのか
平均二乗誤差（MSE）: MSEはRSSをデータの数で割ったもので、予測誤差の平均を示します。
MSEはデータの数によらず誤差の度合いを示すため、データセットの大きさが異なる場合でもモデルの性能を比較しやすいです。
定義式は以下の通りです。
$$\begin{equation} MSE = \frac{1}{n}\sum_{i=1}^{n} (y_i – \hat{y}_i)^2 \end{equation}$$
主な違いはRSSが誤差の総量を、MSEが誤差の平均を表している点です。
詳しくは以下のコンテンツをご覧ください
【MSEを最小化】ガウス・マルコフの定理と最良線形不偏推定量について
この重みベクトルは解析的に求めることができます。
まず損失関数であるRSS(残差平方和)を表します。
$$\begin{equation} \text{RSS}(\mathbf{w}) = (\mathbf{y} – \mathbf{X}\mathbf{w})^T(\mathbf{y} – \mathbf{X}\mathbf{w}) \end{equation}$$
次に重みベクトル\(w\)について微分してイコール\(0\)にします。
$$\begin{equation} \frac{\partial \text{RSS}(\mathbf{w})}{\partial \mathbf{w}} = -2\mathbf{X}^T(\mathbf{y} – \mathbf{X}\mathbf{w}) = 0 \end{equation}$$
$$\begin{equation} \mathbf{X}^T\mathbf{y} – \mathbf{X}^T\mathbf{X}\mathbf{w} = 0 \end{equation}$$
$$\begin{equation} \mathbf{w} = (\mathbf{X}^T\mathbf{X})^{-1}\mathbf{X}^T\mathbf{y} \end{equation}$$
このように解析的に解を求めることができました。ただし、\(\mathbf{X}^T\mathbf{X}\)が正則な場合に限ります。
では、pythonで実装していきます。
CODE
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
# irisデータセットのロード
from sklearn.datasets import load_iris
iris_dataset = load_iris()
df = pd.DataFrame(iris_dataset.data, columns=iris_dataset.feature_names)
df.head()scikitlearnのirisデータセットを使います。
このままでは扱いづらいので上のようにデータフレーム型に変換しています。
df.head()で数行だけデータフレームを見てみましょう。
| sepal length (cm) | sepal width (cm) | petal length (cm) | petal width (cm) | |
| 0 | 5.1 | 3.5 | 1.4 | 0.2 |
| 1 | 4.9 | 3.0 | 1.4 | 0.2 |
| 2 | 4.7 | 3.2 | 1.3 | 0.2 |
| 3 | 4.6 | 3.1 | 1.5 | 0.2 |
| 4 | 5.0 | 3.6 | 1.4 | 0.2 |
一旦欠損値があるかどうかを確認します。
df.isnull().values.sum()こちら結果は0でした、欠損値はないようです。
次は相関行列を見て、各カラム同士の関係を見てみましょう。
corr_mat = df.corr(method='pearson')
corr_mat| sepal length (cm) | sepal width (cm) | petal length (cm) | petal width (cm) | |
| sepal length (cm) | 1.000000 | -0.117570 | 0.871754 | 0.817941 |
| sepal width (cm) | -0.117570 | 1.000000 | -0.428440 | -0.366126 |
| petal length (cm) | 0.871754 | -0.428440 | 1.000000 | 0.962865 |
| petal width (cm) | 0.817941 | -0.366126 | 0.962865 | 1.000000 |
df.corrの引数にはどの方法で相関係数を求めるかを入力することができます。
今回は一般的なピアソンの相関係数にしましたが、スピアマン(df.corr(“spearman”))の順位相関係数やケンドールの順位相関係数(df.corr(“kendall”))にすることもできます。
順位相関係数について詳しく知りたい方は、【外れ値に対処】順位相関係数と相関係数の違いについて | pythonをご覧ください。
さて、相関行列を見ると、sepal lengthとpetal lengthには強い相関があるようにみられます。
今回は、sepal lengthを説明変数にし、petal lengthを目的変数とした単回帰分析を行ってみましょう。
from sklearn import linear_model
X = df[["sepal length (cm)"]]
Y = df[["petal length (cm)"]]
model_lr = linear_model.LinearRegression()
model_lr.fit(X, Y)
#決定係数を算出する
my_result = model_lr.score(X, Y)
my_result決定係数は、0.759945‥となりました。まずまずのようです。
決定係数については、決定係数とは？説明変数の確らしさを図る指標の一つ。をご覧ください。
最後に、係数と定数項を確認しましょう。
w = model.coef_
b = model.intercept_
print(w)
print(b)\(w = 1.85843298, b = -7.10144337\)となりました。
これは、\(petal length = -7.10144337 + 1.85843298*sepal length +u_{i}\)という線形の式ができるということを意味します。
