読みたかった論文を取り寄せたら俺の論文が引用された上でまあまあ正当に批判されていて笑ってしまった。他にも似た形でいろんなことやってる論文がある中でなんでよりにもよって俺のこれを引用するんだよっていうような形にはなってたけど。とはいえみんな各々の哲学で研究して、そしてそのぶつかり合いが技術レベルの押し上げにつながるわけだからまあいいんだけど。俺はデータサイエンスとか語るようになったアレもあって、手法開発ではなんだかんだ実社会における使いやすさを重視してるのでこれでいいんですよ。じゃなかったらClumpinessなんかやんないでしょ。
急にMATLABでRNNがやりたくなってきた（MATLABでRNNシリーズ）
LSTMはじめRNN系の深層学習で何かしらの予測を考えると、そこにはこれまでの分類/回帰という話以外にネットワーク構造におけるSequence-to-Sequence型（一連の流れから一連のデータを吐く） / Sequence-to-One型（一連の流れからその結果としての1つの値を吐く）の違いが出てくる。
前者は時系列データから時系列データを吐くといっても一応いいのかな。たとえば1ヶ月の購買行動をdailyでrecurrentに流し込みながら、その時々の購買行動から次回の購買時期がいつになるのかをその都度予測させるとか。後者だと1ヶ月の購買行動を一通り食わせた上で来月の売上金額を予測するとか。
RNNとかLSTM自体の細かい説明は腐るほどあるので割愛して（ちゃんとやりたければとりあえず”伝説”こと岡谷(2015)『機械学習プロフェッショナルシリーズ 深層学習』から入ればよい）、なおかつTensorflow使った例は腐るほどあるのでMATLABでやる。うちの大学でMATLABと全toolboxが使い放題だということが判明して久々に使いたくなったというのもあるんだけど、個人的にはMATLABでDeep Learningやるのはkerasでlayerをstackしていくのと感覚としてはあまり変わらないと思っている。
ただ、ネットでよくやってるサンプルデータ（Sequence-to-Sequenceで入力の総和を取るだけ）のシミュレーションだとわざわざLSTMなんか使う必要はなくて、いわゆる「浅いネットワーク」の中でも自己回帰系の手法（NAR系のニューラルネットワーク）で十分に対応できる。簡単すぎる問題に難しすぎる手法を当てたところで1匹の蟻に大砲を撃つようなものだし、現象の説明はシンプルであるに越したことはない。
実データ（購買ログとか）への応用はまた近々どっかのpreprintにでも書くことにして今日はサンプルデータだけにする。
サンプルデータの生成
とりあえずtimestepを1〜第T期まで、x, yそれぞれ1変数で、xは毎期乱数で0〜10までの整数として発生、あるt期のy_tは1〜t期までのxの総和に毎期一定のbiasを足していく形として、
$$y_t = \begin{cases} y_{t-1} + x_t + b& (t>1) \\ x_t + b & (t=1) \end{cases} ~~ t \in [1, T]$$
まあつまりは
$$y_t = \sum_{i=1}^{t} x_i + tb$$
の極めてシンプルな形でデータを生成する。
% ---01 Generate Dataset
n_x = 1; n_y = 1;
T = 1000;
X = round(rand(n_x, T)*10);
Y = zeros(n_y, T);
y = 0; bias = 5;
for i = 1:T
y = y + X(i) + bias;
Y(i) = y;
end
X = num2cell(X); Y = num2cell(Y);今回は外部入力を伴う非線形自己回帰ニューラルネットワーク（Nonlinear AutoRegressive models with eXogenous Network, NARXNET）でやる。MATLABだとNonlinear Autoregressive Neural Network with External InputでNARXNETって書いてあるけどいくつか論文みるとメジャーな表現ではないっぽい。
NARXNETはt期以前のOutputが帰還路から隠れ層にもう一度入力される一般的なinput-outputモデル（NARNETとか、線形だとARMAモデルとか）に、毎期説明変数が投入できる構造なので、言ったらこのデータの生成と表現していることはほぼ同じ。非線形のactivationはいらないんだけど、とはいえそれはデータ構造がわかっているからこそ言えることでしかないというのもある。
% ---02 Define Networks
inputDelays = 1:1;
feedbackDelays = 1:1;
hiddenLayerSize = 1;
net = narxnet(...
inputDelays, feedbackDelays, hiddenLayerSize, 'open'...
);
net.trainFcn = 'trainbr';
net.performFcn = 'mse';
[x,xi,ai,t] = preparets(net,X,{},Y);
% Setup Division of Data for Training, Validation, Testing
net.divideParam.trainRatio = 0.5;
net.divideParam.valRatio = 0.25;
net.divideParam.testRatio = 0.25;無駄に構造を複雑にするのは好みではないのと、十分に近似できるはずなので隠れ層内のユニット数（hiddenLayerSize）は1つにしておく。
学習はベイズ正則化で。ただ非線型関数の指定がどこみても見つからない。ここの非線型関数は「MLP（多層パーセプトロン）で近似できさえすればNARX Recurrent Neural Networkになる」[1]という微妙すぎる表現になっている。できれば最初はlinearでシミュレーションしてみたかったんだけど。
% ---03 Train the Network
[net,tr] = train(net,x,t,xi,ai);
yhat = net(x,xi,ai);
e = gsubtract(t,yhat);
performance = perform(net,t,yhat);
se = cell2mat(e).^2;
mse = sum(se)/100;
rmse = sqrt(mse)RMSEで9.0687。散布図を見ると、
分析者としてのこちらはデータの生成構造が完全にわかっているからRMSEが9でもちょっと不満になる（だって線形回帰で説明変数にラグとって投入すれば絶対にもっと綺麗になる）けど、散布図を見ればまあ予測としては決して悪くない。この形で1000個もサンプルを生成させれば最後の方は5桁台になるわけだけど、そこでも誤差は2〜5ぐらいに収まっている。たぶん隠れ層のユニット増やせばもっと綺麗になる。
新規データでテストもする。
% ---04 Cross-validation
Xt = round(rand(n_x, T)*10);
Yt = zeros(n_y, T);
y = 0; bias = 5;
for i = 1:T
y = y + Xt(i) + bias;
Yt(i) = y;
end
Xt = num2cell(Xt); Yt = num2cell(Yt);
[xt,xti,ati,tt] = preparets(net,Xt,{},Yt);
ythat = net(xt,xti,ati);
et = gsubtract(tt,ythat);
performance = perform(net,tt,ythat);
se_t = cell2mat(et).^2;
mse_t = sum(se_t)/100;
rmse_t = sqrt(mse_t)RMSE（テスト）で9.3049、悪くないですね。散布図はほとんど一緒なので割愛。
非線形でもやってみる
やってみるっつってもさっきのデータ生成に非線形関数を挟むだけなんだけど。
% Nonlinear Simulation ------------------------------------
T = 100;
X = round(rand(n_x, T)*10);
Y = zeros(n_y, T);
y = 0; bias = 5;
for i = 1:T
y = y + X(i) + bias;
Y(i) = log(2*y);
end
X = num2cell(X); Y = num2cell(Y);
[x,xi,ai,t] = preparets(net,X,{},Y);
net.trainFcn = 'trainbr';
net.performFcn = 'mse';
% Setup Division of Data for Training, Validation, Testing
net.divideParam.trainRatio = 0.5;
net.divideParam.valRatio = 0.25;
net.divideParam.testRatio = 0.25;
% Train the Network
[net,tr] = train(net,x,t,xi,ai);
% Test the Network
yhat = net(x,xi,ai);
e = gsubtract(t,yhat);
performance = perform(net,t,yhat);
se = cell2mat(e).^2;
mse = sum(se)/100;
rmse = sqrt(mse)
scatter(cell2mat(t),cell2mat(yhat))RMSE=0.0523
まあ標準化してないしデータの形が違うのでこれまでのものと比較はできないのだけれど、意外と綺麗に近似できてるなとは思う。
ちょっと頭の方がうまくできてないし、時系列のエラーを見ると最初の出力が過大だからって補正かけてしばらく過小になってるしと全体的にぐだぐだなのだけれど、とはいえこういう累積していくタイプのモデルというのは t-1 期の予測値に引っ張られるのが帰納的に全値に影響を与えていくから大変なんですよね。それこそ入力を時変にしないと（つまり出力値をrecurrentlyに入力させないと）、
こうやって期首では過大、期末では過小、なのに全体としては誤差が最小化された形になってしまう。機械学習やってる人だとみんな少なからず感じたことがあるだろうと思うけど、機械学習は設定した評価関数を最小化するのに我々が想定していなかったような穴をすり抜け、頑として真面目に働かない時がある（まあ上のは十分に予想できる結果だけどね）。
こんなふうに線形/非線形ともにシミュレーションでは部分的に過大/過小出力になっているのを無理やり調整しているのが見えるんだけど、個人的にはDeep Learningとか含めこういうのやっててなんとなく嫌な気持ちを抱くのは、まさにこういう現象全体のうちあくまで観測できているに過ぎない表層的な部分に合わせて世界を歪めている感じがあるからなんですよね。まあそれでもやるんだけど。
まあ現象の過程を逐一食わせていくとはこういうことなのかな。ここからもう少し複雑な話につなげていく。