PyTorch深層学習⑦線形回帰:実装編
前回は線形回帰の理論的な話をしました。今回は、線形回帰(Linear Regression)をPyTorchのモジュールを使って実装します。
実装することでモデル学習の仕組みの理解が深まります。
今回扱うのは、もっとも単純な線形回帰で、入力と出力の関係を直線で予測します。
例えば、身長から体重を予測する場合を考えます。身長が入力で、体重が出力とします。この身長の数値と体重の数値をグラフにプロットします。そして、その点たちをなるべくうまく通る直線を見つけることを目指します。
線形回帰は多変数の場合でも応用ができるので、シンプルながらも強力な予測手法です。
Python環境の設定
Pythonの仮想環境を作ってPyTorchとJupyterなどをインストールします。
mkdir linear_model
cd linear_model
python3 -m venv venv
source venv/bin/activate
# pip をアップグレードしておく
pip install --upgrade pip
# 必要なライブラリをインストール
pip install torch matplotlib jupyterいつものようにJupyterノートブックを立ち上げてPython3のノートブックを作成してください。
データの準備
直線のデータ
訓練に使うデータを準備します。
まず、必要なライブラリをインポートします。
import numpy as np
import torch
import torch.nn as nn
import matplotlib.pyplot as pltここで torch.nn を nn としてインポートしています。nn はニューラルネットワーク(Neural Network)の意味で、モデル作成に役立つモジュールなどがあります。あとで使います。

まずは、直線のデータを作ってみます。直線は傾きと切片だけ決まるので次のようにデータを作ることができます。
# 傾きと切片
a = 0.5
b = 1.0
# xの値
x = np.arange(-10, 20) # -10から20未満までの整数
# yの値
y = a * x + bこれを表示します。
plt.scatter(x, y)
plt.axvline(0) # 縦軸の線
plt.axhline(0) # 横軸の線
plt.show()x と y が線形の関係にあるのがわかります。当然ではありますが。

これだけだと、面白くないのでデータにノイズを追加して予想し難くしましょう。
ノイズありのデータ
次にノイズを追加したデータを作ります。再現性を確保するためにNumPyの乱数シードを設定します。
# 乱数のシードを設定
np.random.seed(123)
# ノイズを作る
noises = np.random.randn(len(x)) * 2
# 直線にノイズを追加
y = a * x + b + noisesx の個数分だけ正規分布からのランダムな値を作り直線のデータに追加しています。
これを表示します。
plt.scatter(x, y)
plt.axvline(0) # 縦軸の線
plt.axhline(0) # 横軸の線
plt.show()先ほどとは違って直線の関係がわかり難くなっています。

このデータを使って傾きと切片のパラメータを推測できるようにモデルを学習させましょう。
訓練用データ
x をプリントしてみます。
x
1次元のNumPy配列です。y も同様です。
y
では、x と y をPyTorchのTensorに変換します。
x_train = torch.FloatTensor(x).unsqueeze(1)
y_train = torch.FloatTensor(y).unsqueeze(1)unsqueeze(1) としているのは2番目の次元を付け足すためです。シェイプを見てみましょう。
x_train.shape
30個のデータが2次元のテンソルになっています。データをプリントしてみます。
x_train
このように、行列(2次元テンソル)になっており、各行が一つの要素を持つベクトル(1次元テンソル)になっています。
これがPyTorchのモデルに渡す入力データになります。つまり、入力データは一行が一つの入力になっています。
今回は入力変数が x だけなのですが、もし入力変数が2つあれば、各行の値は2つの要素を持つベクトルとなります。
また、y_train はモデルの出力に対する正解データとなるので、やはり一行に一つの値になります。
y_train.shape
y_train のデータを見てみます。

以上で訓練用のデータの準備が完了しました。
今回はバッチやデータローダーは使いません。
線形モデル
nn.Linear
PyTorchの nn にはたくさんのモジュールがあります。モジュールとはデータを入力すると何らかの計算をしてくれるものです。

今回は線形モデルを使いたいので、nn.Linear を使います。
# PyTorchの乱数シード固定
torch.manual_seed(123)
# 入力1変数、出力1変数
model = nn.Linear(1, 1)
torch.manual_seed(123) としてPyTorchの乱数シードを固定しています。これは nn.Linear の初期化においてパラメータの値がランダムに設定されるからです。乱数シードを固定すれば、パラメータの値の再現性を保つことができます。
model = nn.Linear(1, 1) では、入力変数が1つで出力変数が1つと指定しています。よって、model は、入力に傾きを掛けて切片を足して計算したものを出力します。つまりは線形関数としての計算を行います。
ただし、先ほど述べたようにパラメータ(傾きと切片)の値はランダムに初期化されています。
次のようにしてパラメータの値を確認できます。
# generatorからリストにする
params = list(model.parameters())
# スカラーとして取り出す
a = params[0].item()
b = params[1].item()
print(a, b)model.parameters() は generator を返すので一旦リストに変換します。そこから1番目と2番目のパラメータをスカラーとして取り出します。

このセルを何度実行してもパラメータが同じ値で初期化されているのがわかります。
フィード・フォワード
モデルに入力データを挿入して出力を得ることをフィード・フォワードと呼びます。
フィード(Feed)は「餌を与える」といった意味があります。プリンターに用紙を挿入するときもフィードといったりします。
さらに、フォワードは「前へ」という意味になります。これは、後で見るように、ニューラルネットワークの計算グラフでの順方向を意味します。
モデルに入力データをフィード・フォワードすることで、現在のモデルのパラメータによる予測値を計算できます。
# フィード・フォワード
y_hat = model(x_train)
y_haty_hat とは$${\hat{y}}$$の意味で、x_train からの入力値に対して y_train の値を予測したものという意味です。

訓練されていないモデル、あるいは訓練初期のモデルでは、この予測値は y_train とは全く異なる値になっています。
なので訓練を通してパラメータの値を調節していきます。そのため、損失の大きさを図る必要があります。それよってパラメータの更新の度合いや方向が変わってくるからです。
損失関数
MSE
正解値である y_train と予測値 y_hat の差は誤差になります。図で表すと以下のようなイメージになります。

この誤差を二乗して平均を計算したものを全体の誤差の代表値として使います。これを平均二乗誤差(Mean Squared Error、MSE)と呼び、次のように計算します。
$$
\text{MSE} = \dfrac{1}{N} \sum\limits_{i=1}^N (y_i - \hat{y}_i)^2
$$
これによって、誤差の正負に関わらず正解値と予測値の差の大きさを表現できます。
PyTorchを使うとMSEは次のように計算できます。
mean_squared_error = torch.mean((y_hat - y_train)**2)
mean_squared_error.item()
また、nn.MSELoss というモジュールがあり、同じ計算をすることができます。今回はこれを損失関数(Loss Function)として使います。
loss_function = nn.MSELoss()
loss = loss_function(y_hat, y_train)
loss.item()
同じ値が返されました。
勾配
パラメータを調節すると、損失値の増減に影響します。この度合いと方向を知るために各パラメータに対して損失関数の偏微分を計算します。また、すべてのパラメータに対する損失関数の偏微分をまとめたものを勾配と呼びます。
$$
\nabla \text{MSE} = \left[ \dfrac{\partial \text{MSE}}{\partial a}, \dfrac{\partial \text{MSE}}{\partial b} \right]
$$
この辺は前回詳しく解説しました。
勾配がわかるとパラメータを更新する方向がわかります。
例えば、下図では赤の線が正解で、緑の線が予測だとします。ここでパラメータ a(傾き)に対するMSEの偏微分が正だとすると、パラメータ a を増加させると損失がより大きくなることになります。

よって、パラメータ a を減少させて、青の線ように傾きを調節すれば損失値が減少することになります。
また、下図では緑の線が予測であり、パラメータ b(切片)に対するMSEの偏微分が正だとすると、パラメータ b を増加させると損失がより大きくなることになります。

よって、パラメータ b を減少させて、青の線のように切片を調節すれば損失値が減少することになります。
このようにパラメータの値を少しずつ最適化していくのですが、そのためにオプティマイザを使います。
オプティマイザ
SGD
今回は、オプティマイザ(Optimizer)として、SGD(Stochastic Gradient Descent、確率的勾配降下法)を使います。これも前回解説しましたが、次のような更新式を使います。
$$
\boldsymbol{w} \leftarrow \boldsymbol{w} - \alpha \nabla \text{MSE}
$$
ここで$${\boldsymbol{w} = (a, b)}$$で、$${\alpha}$$は学習率です。
つまりは、勾配と逆方向にパラメータを更新するようになっています。また、学習率によって更新の度合いを調整しています。
torch.optim.SGD
PyTorchのオプティマイザは torch.optim の中にあります。次のようにモデルのパラメータを渡し、学習率を指定してSGDを作成します。
# 学習率
learning_rate = 0.01
# オプティマイザの作成
optimizer = torch.optim.SGD(model.parameters(), lr = learning_rate)
実際にどのようにオプティマイザを使ってパラメータの更新をするのかは少し後に解説します。ここでは、オプティマイザが勾配と学習率にしたがってモデルのパラメータの値を更新するとだけ理解してください。
これまで解説したものを図にまとめると以下になります。

これを繰り返すことでパラメータが損失値を小さくするように更新されていきます。
なお、損失関数から勾配を計算するプロセスを誤差逆伝播(Back Propagation)と呼びます。これを実現するためにPyTorchは計算グラフを作ります。
計算グラフ
順方向接続
計算グラフ(Computational Graph)とは、ニューラルネットワークのモデル内の計算を下図のようなグラフの形式で表現したものです。

今回は、入力変数が一つのモデルなので、パラメータとしては「傾き」と「切片」の二つしかありません。より一般に、「重み」(weight)と「バイアス」(bias)と呼びます。なお、「重み」を「パラメータ」の意味で使うこともあります。この場合は、バイアスも重みの一部になります。文脈で判断してください。
計算グラフを順番に辿っていくと、x_train からの入力値にパラメータの重みを掛けて(multiply)してから、バイアスを加算(add)します。これが、予測値 y_hat となります。この予測値 y_hat と正解値 y_train を使ってMSEによる損失値 loss を計算します。
先ほども述べたように、入力値を与えた順方向に計算することをフィード・フォワードと呼びます。
逆方向接続
損失値が計算された後に、各パラメータに対して損失関数の偏微分を行います。つまりは、パラメータに対する損失関数の勾配を計算します。
これは微分の連鎖律を利用しており、損失値から始まって計算グラフを逆方向になぞっていきます。

y_hatをMSEの変数の一つと捉えると、y_hatに対するMSEの偏微分は、$${ \frac{\partial \text{MSE}}{\partial \hat{y}}}$$となります。さらに、パラメータである傾き a と切片 b を y_hat の変数と捉えると、それぞれの偏微分$${\frac{\partial \hat{y}}{\partial a}}$$と$${\frac{\partial \hat{y}}{\partial b}}$$を計算することがでいます。
これらの偏微分を連鎖させることで各パラメータに対するMSEの偏微分を次のように計算することができます。
$$
\frac{\partial \text{MSE}}{\partial a} = \frac{\partial \text{MSE}}{\partial \hat{y}} \frac{\partial \hat{y}}{\partial a}
$$
$$
\frac{\partial \text{MSE}}{\partial b} = \frac{\partial \text{MSE}}{\partial \hat{y}} \frac{\partial \hat{y}}{\partial b}
$$
このように、計算グラフを逆方向へ向かって偏微分を連鎖させていくことで、損失関数から勾配を求める手法を誤差逆伝播法(Back Propagation)と呼びます。
誤差逆伝播法
PyTorchで誤差逆伝播を行うには、次のようなコードを使います。
# 損失値の計算
loss = loss_function(y_hat, y_train)
# パラメータの勾配をリセットする
optimizer.zero_grad()
# 誤差逆伝播を行う
loss.backward()
# パラメータの更新をする
optimizer.step()一つずつ解説します。まず、損失値の計算をします。
# 損失値の計算
loss = loss_function(y_hat, y_train)
loss
ここで損失値 loss に grad_fnという属性がついています。これは、そのテンソルがどの操作(関数)によって作成されたかを示すためです。逆伝播のための計算グラフを構築する際に使用されます。
次に、パラメータに関連する勾配情報をリセットします。
# パラメータの勾配をリセットする
optimizer.zero_grad()具体的には、各パラメータの grad 属性に格納されている値がゼロになります。実は、この grad 属性に各パラメータに対する損失関数の偏微分の値を設定するのが次のコマンドになります。
# 誤差逆伝播を行う
loss.backward()loss.backward()は計算グラフを逆に辿って、.grad_fn を参照しながら微分の連鎖律を使って、各パラメータに対する損失関数の偏微分の値を収納します。
このように自動的に微分の計算を行うことを自動微分(automatic differentiation)と呼びます。逆伝播では計算グラフを使って自動微分を行うことで損失関数の勾配を計算しています。勾配の値は、各パラメータの偏微分の値としてパラメータの .grad 属性に設定されます。
この勾配の値を利用してパラメータを更新するのが次のコマンドになります。
# パラメータの更新をする
optimizer.step()ステップという名前は、「次のステップへの移行」といったところでしょうか。これが呼ばれた時に、パラメータが更新されます。その際に勾配と学習率を使うのがSGDオプティマイザになります。
誤差逆伝播を含めて、これまで解説したものを図にまとめると以下になります。

モデルの学習
訓練ループ
訓練ループ(training loop)では、フィード・フォワードと誤差逆伝播を繰り返しながらパラメータの更新を行います。
一回のループをエポック(epoch)と呼びます。今回は300エポックの訓練を行います。
epochs = 300
for epoch in range(epochs):
# フィード・フォワード
y_hat = model(x_train)
# 損失値の計算
loss = loss_function(y_hat, y_train)
# パラメータの勾配をリセットする
optimizer.zero_grad()
# 誤差逆伝播を行う
loss.backward()
# パラメータの更新をする
optimizer.step()
if epoch % 30 == 0:
print(f'epoch: {epoch:3d} loss: {loss.item():3.2f}')実行結果
実行結果は次のようになります。

180エポックあたりから損失値の値は変わりません。データにノイズがあるので損失値がゼロにならないのは想定内です。
次にパラメータの値を見てみます。
a = model.weight.item()
b = model.bias.item()
print(f'傾き: {a:.2f}')
print(f'切片: {b:.2f}')
これは期待される値(0.5と1.0)にかなり近くなっています。これもノイズのせいで完璧に元の直線を再現できていませんが、訓練がうまくいっているのは確認できます。
この傾きと切片を使ってグラフにしてみましょう。
plt.scatter(x, y)
plt.axvline(0)
plt.axhline(0)
plt.axline((0, b), slope=a, color='red')
plt.show()
いかがでしょうか。
まとめ
今回は線形回帰をPyTorchで実装しました。入力と出力ともに一変数の単純なものですが、色々と応用が効きます。
例えば、nn.Linear(10, 2) とすれば、10個の入力変数から2個の出力変数へと線形変換することができます。
また、データの量が多いときは、バッチやデータローダーを組みわせることで限られたメモリ上でも訓練を行うことができます。
さらに、線形関数による出力に対して非線形の活性化関数を適用するば、分類などに使うこともできます。
次回は分類を行うロジスティック回帰について解説します。
(続く)
