13. Multiple Regression

Multiple Regression

　今日は単純な線型モデルについて説明します。正直解説が必要にないというモデルだと思いますが、記事数が稼げないのでやりたいと思います。

重回帰

　重回帰です。下記の係数aと切片bを学習します。D次元ベクトル$X$と$a$、$b$をもちいて$Y$を求めます。

 Y = a_0 X_0 + ... + a_{D-1} X_{D-1} + b_0 = \sum_{i=0}^D a_i X_i \\ 
\text{where} \ a_D=b_0, X_D = 1

　いろいろなところで計算は見つけることができるので、途中計算はすべてスキップしますが、目的関数は

L = \frac{1}{2} \sum_{i=1}^N |y_i - \tilde{y}_i|^2 

で、これを用いて、

\vec{a} = \left( X^T X \right)^{-1} X^T Y

を導出することができます。愚直に逆行列を計算してもよいのですが、勉強のために修正コレスキー分解をもちいて下三角行列と対角行列に分解した後に逆行列を求めています。特に性能を求めていないので、OpenBlasなどの実装に置き換えたいと思っています。

Ridge正則化

　重回帰にL2正則化項を追加します。これも最終的な部分のみ記載すると、目的関数と係数はそれぞれ下記のようにもとまります。

L = \frac{1}{2} \sum_{i=1}^N |y_i - \tilde{y}_i|^2 + \lambda |a|^2 \\
\vec{a} = \left( X^T X + \lambda I \right)^{-1} X^T Y

$\lambda$が正則化の強さで、大きくすればするほど係数$\vec{a}$の要素の内0に近づく傾向に強くなります。ある程度強くておくことで汎化性能を向上させることができます。$\lambda I$が加わっているのみなので、通常の重回帰とほぼ処理上の変更はありません。

Lasso正則化

　重回帰にL1正則化項を追加します。上記二つは解析的に求めることができましたが、Lassoはそれができません。目的関数は以下で与えられます。通常の重回帰、Ridgeでは第１項の係数の分母にNがありませんでしたが、Lassoでは付け加わっています。

L = \frac{1}{2N} \sum_{i=1}^N |y_i - \tilde{y}_i|^2 + \lambda |a| \\

FortLearnerではCordinate Descentで実装しています。さらに良いアルゴリズムとしてISTAとFISTAがあるようです。ここでは実装できているCordinate Descentのみを見ていくことにします。
　Cordinate Descentでは特徴量ひとつづつの更新を行い、一通り更新が済んだら再度一つ目の特徴量に戻ってという風に更新が行われます。上記の目的関数は３つの場合に分けることができます(参考)。

\frac{\partial L}{\partial a_i} = 
a_i \sum_{n=1}^N x_{ni}^2 
+ \sum_{n=1}^N x_{ni} \left( a_0 + \sum_{k \neq i}^D a_k x_{nk} - y_n \right) + \lambda N \cdot \text{sign} (a_i)

これが０になるような条件で更新を行うので、$a_i$について解いてみると、

a_i = 
\frac{
        \sum_{n=1}^N x_{ni} \left( y_n - a_0 - \sum_{k \neq i}^D a_k x_{nk} \right) - \lambda N \cdot \text{sign}(a_i)
    }{
        \sum_{n=1}^{N} x^2_{ni}
    }

しかしこれでは、問題があります。両辺に$a_i$が存在してしまっています。分母は必ず正であるため、

\sum_{n=1}^N x_{ni} \left( y_n - a_0 - \sum_{k \neq i}^D a_k x_{nk} \right) - \lambda N > 0 \\
\sum_{n=1}^N x_{ni} \left( y_n - a_0 - \sum_{k \neq i}^D a_k x_{nk} \right) + \lambda N < 0

の場合にはそれぞれの値を$a_i$を代入します。最後に分子が、

- \lambda N  < \sum_{n=1}^N x_{ni} \left( y_n - a_0 - \sum_{k \neq i}^D a_k x_{nk} \right) < \lambda N

は、正負どちらでもない、つまり0にしてしまいます。

　FortLearnerでは、高速化しようとして（本当に速くなっているのか確かめてないのですが）、Loop Unrollingをしています。１年以上前のコードをノートも残さずに書いてしまっていたので、自分でも何が何やらわからなくなってしまいました。今後は元のコードは別できちんと残しつつ、拡張していこうと思います。

Logistic Regression

　最後にロジスティック回帰です。重回帰と同様に、係数と切片を学習しますが、出力時にシグモイド関数を通します。

\sigma(x) = \frac{1}{1+\exp(-x)} \\
h_\theta(\vec{x}_i) 
    = \sigma\left(\vec{\theta} \vec{x} + b\right) 

また、この時目的関数は、

J = -\frac{1}{N} \sum^N_{i=1} \left\{ y_i \ln\left(h_\theta(\vec{x}_i) \right) 
      + (1-y_i) \ln \left(1-h_\theta(\vec{x}_i)\right)\right\}

で定義され、GradientとHessianは以下のように計算されます。

\frac{\partial}{\partial \theta_j} J = \frac{1}{N} \sum^N_{i=1} \left( h_\theta(\vec{x}_i) - y_i \right) x_{ij} \\
\frac{\partial^2}{\partial \theta_j \partial \theta_k} J = 
\frac{1}{N} \sum^N_{i=1} h_\theta(\vec{x}_i) \left( 1-h_\theta(\vec{x}_i) \right)  x_{ij} x_{ik}

　単純な勾配降下法で計算してもよいのですが、Scikit-learnのデフォルトではLimited-Memory BFGSが用いられています。FortLearnerでは、BFGS法で実装してみました。ここら辺は本当に試してみただけであまり、最適化されていません。

　実装していて時間がかかったのは、BFGSの読み替え（記載されているアルゴリズムの文字が自分のプログラムのどの変数に対応するのか）に非常に時間がかかりました。いくつかのURLを参考にしたのですが、当然一貫性はないので、毎回頭がこんがらがりました。さらに、Scikit-learnに比べ遅い上に精度が悪いという結果になってしまいました。

　Sklearnではscipy内で実装されたものを利用しています。そしておそらくこのFortranコードを使っていると思うのですが、４０００行弱あってとてもではないですが、読み解くことができませんでした。おそらく単純なPseudoコードでは表せない様々な経験的な改善点があるのではないかと思います。

以上