目次
線形代数はタイヘン
前回・前々回、ベクトルを取り上げました。ベクトルの次ということで行列をはじめていきます。ベクトルと行列は線形代数学という大きな数学へと結実します。
1冊にまとめられた線形代数のテキストを自力で終わらせるには、相当な根性と実力が必要になります。
前回取り上げたコサイン類似度とは、データをベクトル化することでn次元空間における2つのn次元ベクトルのなす角（実際にはそのコサイン）を計算することで比較するというアイディアです。
ベクトルは3次元が4次元に次元が大きくなってもさほど計算が面倒なことはありません。行列の場合、次元が1つ大きくなるだけで行列式・逆行列の計算はタイヘンになります。スラスラと手計算できるのは2次元どまりです。3次元、4次元の行列の計算ですら容易に手計算できません。
手計算ができないまま行列式・逆行列・連立分程式の仕組みをものにすることはできません。計算量の膨大さこそが行列にとって最大の壁です。
そこでPythonの登場です。
ベクトルと同様にNumPyライブラリを用いて、行列は容易に定義でき、NumPyに用意された関数を用いることで行列式・逆行列・連立分程式といった行列の様々な計算を1行のコードで済ませることができます。
これまでの連載同様にPythonで線形代数を学ぶことはもちろん可能です。しかし、それでも線形代数のテキストの章立ての長さは長いまま変わることはありません。線形代数とは、連立方程式の解法をいかにシステマティックに行うかという問題が原点にあります。
はたして行列という新しい概念が誕生し、そこに見事な美しい法則が発見されていくことになったのです。行列の美しい法則の結晶が線形代数・線形代数学・線型代数・線型代数学（用語の組み合わせは4通り）です。
「NumPyで行列」として行列の計算を紹介していきます。行列は面白そうと思ってもらえるプログラムを用意しました。それが回帰分析です。
多数のデータをもっともうまく近似する直線（曲線）を求めることで、予測シミュレーションを可能にする回帰分析はデータサイエンスの基本です。
Pythonによる回帰分析
さっそくPythonによる回帰分析のコードを見ていただきましょう。
用意したデータは10人分の身長（cm）と体重（kg）です。まずこれらをnp.array()関数を用いて説明変数ベクトル（配列）xと目的変数（体重）ベクトルyを定義します。
次に説明変数（身長）から説明変数行列Xを定義します。はたして行列Xに対して転置行列、逆行列の計算をして最後に目的変数（体重）ベクトルyをかけることで回帰式の係数（回帰係数）が得られます。この回帰係数の計算コードはたった1行です。
最後に計算結果を出力します。
（身長と体重）回帰分析プログラム 「regression1.py」
>>> import numpy as np
>>>
>>> # データ定義
>>> x = np.array([170, 166, 157, 164, 157, 160, 163, 160, 153, 165]) # 説明変数 身長
>>> y = np.array([57, 53, 47, 52, 47, 53, 48, 48, 43, 57]) # 目的変数 体重
>>>
>>> # 説明変数行列Xの定義
>>> ones = np.ones(len(x))
>>> X = np.array([ones, x]).T
>>>
>>> # 回帰係数の計算
>>> theta = np.linalg.inv(X.T @ X) @ X.T @ y
>>>
>>> print(f'説明変数行列 X =\n {X}')
>>> print(f'回帰係数 θ_0 = {theta[0]}')
>>> print(f'回帰係数 θ_1 = {theta[1]}')
>>> print(f'回帰式 y = {theta[0]} + {theta[1]} * x')
このプログラムファイルは次からダウンロードできます。
https://drive.google.com/file/d/1U1m17D--NUpoZ7wTijzbXlHCPfLb6Ibv/view?usp=sharing
実行結果
次がregression1.pyの実行結果です。最後に回帰式がxの1次式として出力されています。
これより身長が175cmであれば体重はいくつになるか計算できます。回帰式のxに175を代入して計算してみます。結果は約61.3kgです。
>>> -78.77006507593245 + 0.8004338394793968 * 175
>>> 61.30585683296201
この例のように説明変数が1つの場合には回帰式はxの1次式となるので、グラフは直線です。これを回帰直線といいます。
10人分のデータの散布図と回帰直線を描くプログラムは、回帰分析プログラム「regression1.py」にmatplotlibライブラリを用いたグラフ描画部分を付け加えたのが回帰分析プログラム「regression2.py」です。
（身長と体重）回帰分析プログラム「regression2.py」
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>>
>>>
>>> # データ定義
>>> x = np.array([170, 166, 157, 164, 157, 160, 163, 160, 153, 165]) # 説明変数 身長
>>> y = np.array([57, 53, 47, 52, 47, 53, 48, 48, 43, 57]) # 目的変数 体重
>>>
>>>
>>> # 説明変数行列Xの定義
>>> ones = np.ones(len(x))
>>> X = np.array([ones, x]).T
>>>
>>> # 回帰係数の計算
>>> theta = np.linalg.inv(X.T @ X) @ X.T @ y
>>>
>>> print(f'説明変数行列 X =\n {X}')
>>> print(f'回帰係数 θ_0 = {theta[0]}')
>>> print(f'回帰係数 θ_1 = {theta[1]}')
>>> print(f'回帰式 y = {theta[0]} + {theta[1]} * x')
>>>
>>> # データ散布図と回帰直線 グラフ描画
>>> margin = 1
>>> x_min, x_max = x.min()-margin, x.max()+margin # xの範囲
>>> xx = np.arange(x_min, x_max, 1) # 回帰直線グラフ用変数:xx
>>>
>>> # yy = θ_0 + θ_1 * xx # 回帰直線 グラフ用変数:xx,yy
>>> yy = + theta[0] + theta[1] * xx # 回帰直線 グラフ用変数:xx,yy
>>>
>>> plt.scatter(x, y, c='Blue') # データを青丸でプロット
>>> plt.plot(xx, yy, c='Red') # 回帰直線を赤線でプロット
>>> plt.xlabel('x:height(cm)')
>>> plt.ylabel('y:weight(kg)')
>>> plt.show()
このプログラムファイルは次からダウンロードできます。
https://drive.google.com/file/d/1ZE5B6X7Sp5ugeBgvb7DZgyEsjsZIE_mg/view?usp=sharing
実行結果
numpy.polyfit関数を使った回帰分析
NumPyライブラリには回帰分析のためのpolyfit関数が用意されています。これを用いればコードはもっと簡単になります。
>>> theta_1 , theta_0 = np.polyfit(x, y, 1)
第3引数を1とすることで回帰式を1次式と指示します。
回帰係数はxの係数（theta_1）、切片（theta_0）の順序で出力されます。
（身長と体重）回帰分析プログラム「regression3.py」
>>> import matplotlib.pyplot as plt
>>>
>>> # データ定義
>>> x = np.array([170, 166, 157, 164, 157, 160, 163, 160, 153, 165]) # 説明変数 身長
>>> y = np.array([57, 53, 47, 52, 47, 53, 48, 48, 43, 57]) # 目的変数 体重
>>>
>>> # 回帰係数の計算
>>> theta_1 , theta_0 = np.polyfit(x, y, 1)
>>>
>>> print(f'回帰係数 θ_0 = {theta_0}')
>>> print(f'回帰係数 θ_1 = {theta_1}')
>>> print(f'回帰式 y = {theta_0} + {theta_1} * x')
>>>
>>> # データ散布図と回帰直線 グラフ描画
>>> margin = 1
>>> x_min, x_max = x.min()-margin, x.max()+margin # xの範囲
>>> xx = np.arange(x_min, x_max, 1) # 回帰直線グラフ用変数:xx
>>>
>>> # yy = θ_0 + θ_1 * xx # 回帰直線 グラフ用変数:xx,yy
>>> yy = + theta_0 + theta_1 * xx # 回帰直線 グラフ用変数:xx,yy
>>>
>>> plt.scatter(x, y, c='Blue') # データを青丸でプロット
>>> plt.plot(xx, yy, c='Red') # 回帰直線を赤線でプロット
>>> plt.xlabel('x:height(cm)')
>>> plt.ylabel('y:weight(kg)')
>>> plt.show()
このプログラムファイルは次からダウンロードできます。
https://drive.google.com/file/d/1X5KCOubqdvUyl5NPIBMthXbqdfTFfjwj/view?usp=sharing
実行結果
回帰係数は先の結果とほぼ同じ結果になることが確認できます。
回帰分析の鍵「最小2乗法」
いかがでしたでしょうか。身長、体重のデータを近似する回帰直線を求めるコードはこのように数行でできてしまいます。それもこれも次の1行のおかげです。
>>> # 回帰係数の計算
>>> theta = np.linalg.inv(X.T @ X) @ X.T @ y
なぜこのような行列の計算（行列の積、転置行列、逆行列）によって回帰係数が得られるのか。今後の連載でじっくり見ていきたいと思います。そこで展開されるのは驚きの連続です。
時代は今から200年前のフランスに遡ります。1789年7月14日、バスティーユ襲撃が起きます。所謂、フランス革命です。フランス革命の中で生まれたのがメートルです。
登場人物は全員数学者、ラグランジュ、ラプラス、ドランブルそしてガウス。メートルという統一単位をつくるために数学は総動員されました。はたして彼らの闘いの中から新しい数学──最小2乗法が誕生しました。最小2乗法のおかげで回帰係数の計算が可能になります。
上記の1行の数式は、行列とメートルとフランス革命という壮大な物語から生みだされました。次回以降をご期待ください。
