ニューラルネットワークでは分類問題を扱うことが多いですが、出力層の活性化関数を変えることで回帰問題も扱えます。この記事では出力層の活性化関数をReLU関数をにして、打ち切りデータの回帰分析を行ってみます。
きっかけはCourseraのDeepLearningのコースを振り返っていたら「線形活性化関数を使うかもしれないのは普通は出力層だけ。ReLU関数にすると出力が正の場合の回帰もできるよ」とAndrew Ng先生が言っていたことより。「あれ、これ打ち切り回帰できるよな」と思って実装してみました。
リポジトリ：https://github.com/koshian2/TruncatedRegressionNN
出力層の活性化関数
出力層の活性化関数はニューラルネットワークの活性化関数の中でも特別な意味を持ちます。隠れ層では適当にReLUやtanhを使っておけば差し障りないことが多いですが、出力層では問題設定が変わります。
例えば、分類問題では出力層の活性化関数をシグモイド関数（2値分類）、ソフトマックス関数（多クラス分類）とします。出力層の活性化関数を変えることで回帰問題にスライドさせることが可能です。出力層を線形活性化関数（活性化関数がない）とすると予測値が実数全体になる一般的な回帰問題になります。同様に出力層をReLU関数にすると予測値が0以上になる打ち切りデータの回帰分析ができます。
打ち切りデータとは
あるところまでは0で、あるところから値が上がり始めるデータです。例えば不動産価格、価格は（普通は）マイナスにはなりません。または降水量や降雪量、こちらは絶対にマイナスにはなりません。
以下は気象庁から取得した、1953年～2017年の旭川の気温と降雪量の関係です。横軸が日の平均気温、縦軸が月間降雪量です。全て月次データです。
このように気温が10度以上では全く降雪がなく、気温が0度～5度を下回ると急に降雪量が増えるという打ち切りデータの典型例です。
今回は特にやりませんが、統計ではこのような打ち切りデータに対しては通常Tobitモデルを使います。ただし、PythonのライブラリではTobitモデルができない（githubから探せばあります）ので、普通はRを使うべきでしょう。Tobitモデルも結構ややこしいので、今回は全く別のニューラルネットワークからのアプローチを行います。ニューラルネットワークとして実装するととてもシンプルです。
ニューラルネットワークにおけるReLU活性化関数は次の通りです。
$$\rm{ReLU}(z) = \max(z,0)$$
なんとなく打ち切りデータに使えそうな気がしますよね。これを実装していきます。
隠れ層なしの場合
asahikawa.csvには1列目に気温が、2列目に降雪量があります。
| 日平均気温 | 降雪深 |
|---|---|
| -10.9 | 98 |
| -9.6 | 102 |
| -2.3 | 64 |
| 3.7 | 7 |
| 10.5 | 0 |
文字列があるとインポートがめんどくさいので日時は省いてしまいましたが、1953年1月、2月…の順に2017年12月までデータが並んでいます。気温から降雪量を推定します。以下のニューラルネットワークを作ります。
ニューラルネットワークとしてはこれ以上ない簡単なものです。
import numpy as np
import matplotlib.pyplot as plt
from keras.models import Sequential
from keras.layers import Dense
from keras.optimizers import Adam
# 1953年～2017年の月次データ（旭川）
data = np.loadtxt("asahikawa.csv", delimiter=",", skiprows=1)
X, y = data[:, 0], data[:, 1]
X = X[:, np.newaxis] #行列にする
model = Sequential()
model.add(Dense(1, activation="relu", input_shape=X.shape[1:]))
model.compile(loss="mean_squared_error", optimizer=Adam(lr=1), metrics=["mse"])
n_epoch = 200
history = model.fit(X, y, epochs=n_epoch, batch_size=128).history
分類問題と異なるところは、損失関数にmean_squared_error（平均二乗誤差）を使っているところです。
pred_X = np.arange(-20, 30, 0.1)[:, np.newaxis]
pred_y = model.predict(pred_X)
print(model.get_weights())
plt.plot(X, y, ".")
plt.plot(pred_X, pred_y)
plt.xlabel("Temperature")
plt.ylabel("Snowfall")
plt.plot()
plt.show()
超浅いモデルなので一瞬で答えが出てきます。いい感じですね。
ニューラルネットワークは非凸最適化なので、収束した解＝大域最適解となる保証がありませんが、今回やった隠れ層のないモデルではほぼ同じような結果となりました。次にやる隠れ層のあるモデルだと答えが毎回変わるので注意が必要です。
係数は以下の通りです。
[array([[-8.899647]], dtype=float32), array([81.43994], dtype=float32)]
これは、降雪量をy、気温をxとして、
$$\hat{y} = \max(-8.899647x+81.43994, 0) $$
と推定するモデルです。maxの中をイコール0として解くと、平均気温が9.15度を下回ると雪が降り始めるということがわかりますね。
では平均気温が-10度のときの降雪量の推定値はどのぐらいでしょうか？x=-10を代入すると、-8.899647*(-10)+81.43994=170.4cmとなります。
隠れ層ありの場合
実はデータをよく見てると、雪が降り始めるのは気温9度ではなく5度ぐらいで、積雪と気温の関係も二段階になっているのでは？という疑問がわきます。そんなときは隠れ層を追加します。
ちなみに回帰分析のときにやりがちな、例えば気温を-5度した別の変数を追加するという作業は必要ありません。隠れ層を追加すれば勝手に学習してくれます。コードも、
model = Sequential()
model.add(Dense(2, activation="relu", input_shape=X.shape[1:]))
model.add(Dense(1, activation="relu"))
model.compile(loss="mean_squared_error", optimizer=Adam(lr=1), metrics=["mse"])
1行追加しただけです。ただし、結果はニューラルネットワークが非凸最適化ゆえに局所解が複数あります。初期値の乱数次第では試行ごとに全く別の結果が出てきて、特に隠れ層を入れた場合に顕著に現れます。たまたま上手くいった例では、
[array([[-1.156841, 9.951706]], dtype=float32), array([ 10.494091, -12.452746],
dtype=float32), array([[ 2.446894 ],
[-2.1623268]], dtype=float32), array([81.09545], dtype=float32)]
となりました。5度以下になると急速に積雪が増え、氷点下では気温が下がっても積雪量の増加はゆるやかになっていく、というのは感覚的にしっくりきます。あまり上手くいかない例だと、
[array([[-3.187832, -5.999301]], dtype=float32), array([ 20.281013, -43.23949 ],
dtype=float32), array([[ 4.608462 ],
[-1.7430465]], dtype=float32), array([0.83780617], dtype=float32)]
失敗例として、隠れ層がない場合とほぼ同じ結果に収束してしまったケースや、y=0という直線に収束したケースもありました。何回かやってみてその中から良さそうなものを選ぶというのをおすすめします。もし結果に再現性をもたせたければ初期値をランダムではなく固定値とすることも検討します。
上手くいったほうの例の係数をもう少し詳しく見ます。ニューラルネットワークのForward Propagationの定義から、次の式で表されます。ここでgはReLU活性化関数です。
\begin{align}
Z_1&= \left[\begin{array}{r}
-1.156841 \\ 9.951706
\end{array}\right] X +
\left[\begin{array}{r}
10.494091 \\ -12.452746
\end{array}\right] \\
A_1&=g(Z_1)\\
Z_2 &= \left[\begin{array}{r} 2.446894 && -2.1623268
\end{array}\right] A_1 + 81.09545 \\
\hat{y} &= g(Z_2)
\end{align}
(1)気温が-10度の場合
\begin{align}
Z_1&= \left[\begin{array}{r}
-1.156841 \\ 9.951706
\end{array}\right] (-10) +
\left[\begin{array}{r}
10.494091 \\ -12.452746
\end{array}\right] =
\left[\begin{array}{r}
22.0625 \\ -111.97
\end{array}\right] \\
A_1&=g(Z_1)=\left[\begin{array}{r}
22.0625 \\ 0
\end{array}\right]\\
Z_2 &= \left[\begin{array}{r} 2.446894 && -2.1623268
\end{array}\right]\left[\begin{array}{r}
22.0625 \\ 0
\end{array}\right] + 81.09545 = 135.0801 \\
\hat{y} &= g(Z_2) = 135.0801
\end{align}
気温が-10度のとき、積雪の推定値は135.1cm。
(2)気温が2度の場合
\begin{align}
Z_1&= \left[\begin{array}{r}
-1.156841 \\ 9.951706
\end{array}\right] 2 +
\left[\begin{array}{r}
10.494091 \\ -12.452746
\end{array}\right] =
\left[\begin{array}{r}
8.180409 \\ 7.450666
\end{array}\right] \\
A_1&=g(Z_1)=\left[\begin{array}{r}
8.180409 \\ 7.450666
\end{array}\right]\\
Z_2 &= \left[\begin{array}{r} 2.446894 && -2.1623268
\end{array}\right]\left[\begin{array}{r}
8.180409 \\ 7.450666
\end{array}\right] + 81.09545 = 85.00127 \\
\hat{y} &= g(Z_2) = 85.00127
\end{align}
気温が2度のとき、積雪の推定値は85.0cm。
(3)気温が10度の場合
\begin{align}
Z_1&= \left[\begin{array}{r}
-1.156841 \\ 9.951706
\end{array}\right] 10 +
\left[\begin{array}{r}
10.494091 \\ -12.452746
\end{array}\right] =
\left[\begin{array}{r}
-1.07432 \\ 87.06431
\end{array}\right] \\
A_1&=g(Z_1)=\left[\begin{array}{r}
0 \\ 87.06431
\end{array}\right]\\
Z_2 &= \left[\begin{array}{r} 2.446894 && -2.1623268
\end{array}\right]\left[\begin{array}{r}
0 \\ 87.06431
\end{array}\right] + 81.09545 = -107.166 \\
\hat{y} &= g(Z_2) = 0
\end{align}
気温が10度のとき、積雪の推定値は0cm。
補足
積雪の例では、特に必要ありませんでしたが、打ち切りの値が0ではなく例えば10で打ち切る場合、または打ち切りがmaxではなくminの場合は、最後の1ユニットのReLUのあとに、1個線形活性化関数の層をつくってこれを出力層とすればよいでしょう（要はReLUの線形変換）。
ReLUという関数の性質から、ニューラルネットワークを用いて打ち切りデータの回帰分析を簡単に実装することができました。