多変量解析

線形回帰・Lasso・Ridgeを実装してみよう!

こんにちは!

前に、正則化スパース推定を紹介しました。

でも、式や概念の説明だといまいちどんなことが起きているのかわからない・・・

そこで、今回はR言語を用いて実際に解析してどんな違いがあるのかを確認してみます!

用いる手法

今回比較するのは、通常の線形回帰・Lasso回帰・Ridge回帰です。
簡単に、三つの手法についておさらいします。

通常の線形回帰は残差\({\bf y}-{\bf X}{\bf \beta}\)が最小になるように\(\beta\)を推定します。
つまり、次の式を最小化するように\(\beta\)を推定します。
\begin{eqnarray*}
\|{\bf y}-{\bf X}{\bf \beta}\|_2^2\\
\end{eqnarray*}

Ridge回帰は正則化手法であり、多重共線性などが起きていて通常の線形回帰ではパラメータが推定できないような場面でも罰則を与えることで推定できるようにしたものです。
次の式を最小化するように\(\beta\)を推定します。
\begin{eqnarray*}
\|{\bf y}-{\bf X}{\bf \beta}\|_2^2+\lambda\|{\bf \beta}\|_2^2\\
\end{eqnarray*}
第2項は正則化項と呼ばれ、なるべくパラメータ\(\beta\)の絶対値を小さくするような働きをします。
ここで、\(\lambda\)は正則化パラメータであり、解析者が設定するかクロスバリデーションなどで決定する必要があります。

Lasso回帰も正則化手法です。Ridge回帰とは異なりパラメータのいくつかを0に推定できるという特徴があります。
これをスパース推定と呼びます。次の式を最小化するように\(\beta\)を推定します。
\begin{eqnarray*}
\|{\bf y}-{\bf X}{\bf \beta}\|_2^2+\lambda\|{\bf \beta}\|_1\\
\end{eqnarray*}

今回は、線形回帰はglm()関数、LassoとRidgeはglmnetパッケージのglmnet()関数を用います。

人工データ生成

データは、実データではなく自分で発生させる人工データを用います。
今回のデータの生成モデルを次に示します。

\begin{eqnarray*}
y_i&=&1+{\bf x}_i^T{\bf \beta}+\epsilon_i\\
{\bf \beta}&=&(1,0,-1,0,0)^T\\
{\bf x_i}&\sim& N({\bf \mu},{\bf \Sigma})\\
{\bf \mu}&=&{\bf 0}_5\\
{\bf \Sigma}&=&{\bf I}_5\\
\epsilon_i&\sim& N(0,1^2)\\
i&=&1,2,\cdots,n\\
\end{eqnarray*}

説明変数は、平均0・独立で発生させています。
また、\(\beta\)の設定から変数2・変数4・変数5は目的変数に寄与していません。

データ解析

先ほど説明したデータで解析をしてみましょう!

サンプルサイズ\(n\)は1100とします。
100個のデータを学習データ、残りの1000個のデータをテストデータとして実験します。

それぞれ3手法で推定された回帰係数と予測誤差(MSPE;Mean Squared Prediction Error)を計算してみましょう。

通常の線形回帰で推定された回帰係数
(Intercept) 1.110727505
x1 0.862992782
x2 -0.024437489
x3 -0.757605621
x4 -0.098294387
x5 -0.008577822

Lasso回帰で推定された回帰係数
(Intercept) 1.0713394
x1 0.7781268
x2 .
x3 -0.6778876
x4 .
x5 .

Ridge回帰で推定された回帰係数
(Intercept) 1.08903319
x1 0.81095492
x2 -0.02288852
x3 -0.71278554
x4 -0.09335028
x5 -0.01748969

(Intercept)は定数項を表しています。
どうでしょうか、Lasso回帰は変数1と変数3以外は完全に0になっていることがわかります。
つまり、自動的に変数選択ができました!これは嬉しいですね!

また、推定された回帰係数の2乗和を計算すると次のようになりました。
通常の線形回帰
> sum(as.matrix(result_glm$coefficients)^2)
2.562771

Lasso
> sum(coefficients(result_lasso,s=”lambda.min”)^2)
2.212781

Ridge
> sum(coefficients(result_ridge,s=”lambda.min”)^2)
2.361248

やはり、通常の線形回帰と比較して小さくなっています。

最後に、MSPEを比較します。
通常の線形回帰
> mean((data_test$y-yhat_glm)^2)
1.119285

Lasso
> mean((data_test$y-yhat_lasso)^2)
1.171355

Ridge
> mean((data_test$y-yhat_ridge)^2)
1.153019

最も良いのは通常の線形回帰になりました。
これはサンプルをたくさんあって、説明変数に相関もなく多重共線性が起きていないから当たり前ですね!

もっと、変数を増やしてみたりするとLassoやRidgeがMSPEで勝つこともあるでしょう。

ちなみに、このプログラムでは一回の計算なので、結果の信頼性には欠けます。
手法の精度比較を行うのであれば、シミュレーションを何回も繰り返して、その平均やばらつきをみて比較することが望ましいです。このことに関してはまた違う記事で紹介します!