パイソンメルマガ「Python で、情報、数学再学習」 第２１号 回帰分析
２０２３年１１月２０日 ２１号
このメルマガ(パイマガ)では、Python を使って、数学の再学習、高校「情報」を勉強していきます。
今回は、回帰分析、重回帰分析を勉強します。
今回の内容は、
「NumPyで行列 回帰分析 その8 最小2乗法｜Pythonで数学を学ぼう！ 第29回」
を参考させてもらいました。
この記事のプログラムを実行するために、以下の初期プログラムを実行してください。
#今回のプログラムを動作させるための初期設定プログラム
#初期設定
!pip install -q japanize_matplotlib
import numpy as np
import matplotlib.pyplot as plt
import japanize_matplotlib
from mpl_toolkits.mplot3d import Axes3D #3Dplot
# データ用意
X1 = [165, 170, 172, 175, 170, 172, 183, 187, 180, 185] # 身長
X2 = [23, 24, 25, 27, 25, 24, 28, 29, 28, 29] # 足の長さ
Y = [52, 60, 63, 65, 78, 70, 72, 85, 90, 94] # 体重
回帰分析 Regression Analysis
身長から体重などを予想する式を求めることを、回帰分析といいます。その式を回帰直線といいます。
回帰分析は、人工知能の基本理論の一つです。今回は、その基本を説明します。
下のプログラムを実行して、赤く表示されている線が回帰直線です。
# 回帰直線 y = ax + b
x1 = np.array(X1)
y = np.array(Y)
# 回帰直線を計算するため準備
ones = np.ones(len(x1))
X = np.array([ones, x1]).T
# 回帰直線の係数を求める
theta = np.linalg.inv(X.T @ X) @ X.T @ y
print(f'回帰直線 y = {theta[1]:.2f}x + ({theta[0]:.2f})')
# グラフ描画
margin = 1
x1_min, x1_max = x1.min()-margin, x1.max()+margin
xx = np.arange(x1_min, x1_max, 1)
yy = theta[1] * xx + theta[0]
plt.scatter(x1, y, c='Blue') #散布図
plt.plot(xx, yy, c='Red') #直線
plt.xlabel('x:身長（cm）')
plt.ylabel('y:体重（kg）')
青い点は、１０人分の身長と体重のデータです。赤い直線は、回帰直線といい、身長と体重の関係を表し、身長のから体重を予想するときに、使用します。
青い点からこの直線の式を求めることを回帰分析といい、この直線を回帰直線といいます。
求め方は、直線を仮定し、それぞれの青い点から、直線へのｙ軸での距離(縦の距離)の２乗を合計し、その最小値から計算します。この計算方法は、最小２乗法といいます。
直線の式、y = ax + b の a, b を最小２乗法を使って計算します。上のプログラムでは、10 行目の
theta = np.linalg.inv(X.T @ X) @ X.T @ y
です。X の中に、１０人分の身長と体重のデータが入っており、np.linalg.inv( ) 関数内で、最小二乗法で計算し、直線の a, b にあたるものを、theta[1], theta[0] で返しています。
重回帰分析 Multiple regression analysis
上のグラフでは、身長から体重を求める回帰直線を表示しました。これに足の長さの要素を追加し、身長と足の大きさから、体重を予想することにします。数学的には、２次元から３次元となります。
予想に使う変数が２つ以上になると、重回帰分析といいます。
原理は、いっしょで、次元が増えるだけです。通常、３次元は、x 軸, y 軸, z 軸 で表現されますが、x1 軸, x2 軸, y 軸 とします。４次元、５次元と次元が増えると、x4, x5 軸となります。
ここでは、３次元ですので、計算で求まる式は、回帰直線でなく、回帰平面となります。
以下のプログラムで実行しています。
平面の式の形式は
y = ax1 + bx2 + c
です。
この a, b, c を計算して、theta[1],theta[2], theta[0] で返してます。
# 回帰平面 y = ax1 + bx2 + c
x1 = np.array(X1)
x2 = np.array(X2)
y = np.array(Y)
# 回帰平面を計算するための準備
ones = np.ones(len(x1))
X = np.array([ones, x1, x2]).T #重回帰だと次元が増える
# 回帰平面の係数を求める
theta = np.linalg.inv(X.T @ X) @ X.T @ y
print(f'回帰平面 y = {theta[1]:.2f}*x1 + {theta[2]:.2f}*x2 + {theta[0]:.2f}')
# グラフ描画
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
mesh_size = 1
margin = 0.01
x1_min, x1_max = x1.min()-margin, x1.max()+margin
x2_min, x2_max = x2.min()-margin, x2.max()+margin
x1_range = np.arange(x1_min, x1_max, mesh_size)
x2_range = np.arange(x2_min, x2_max, mesh_size)
xx1, xx2 = np.meshgrid(x1_range, x2_range)
yy = theta[1] * xx1 + theta[2] * xx2 + theta[0]
ax.scatter3D(x1, x2, y, color='Blue') #散布図
ax.plot_surface(xx1, xx2, yy, color='red', alpha=0.5) #平面
ax.set_xlabel("x1:身長（cm）")
ax.set_ylabel("x2:足の長さ（cm）")
ax.set_zlabel("y:体重（kg）")
上で、求めた式を使って、実際に試してみましょう、身長 177 cm で、足が、27 cm の人の体重を予想すると
回帰平面 y = 0.41x1 + 3.60x2 + -93.37
を使って
x1 = 177
x2 = 27
y = 0.41*x1 + 3.60*x2 + -93.37
print(y) #76.39999999999998
となり、約76kgと予想できます。
回帰平面の係数を計算するのに使った関数
theta = np.linalg.inv(X.T @ X) @ X.T @ y
の数学の分野は、線形代数(行列演算)と言われています。こちらは、また、別の号で扱います。
Discussion