ホーム » 「正則化」タグがついた投稿

タグアーカイブ: 正則化

Lassoの最適化アルゴリズムの考察

2019年2月25日 / コメントする

こんにちは。

前回の最後に、
– Coordinate descent（座標降下法）
– ISTA (メジャライザー最小化)
– FISTA (高速化)
の3つの最適化方法を紹介しました。今回は上の2つを解説します。そして次回実装へうつりましょう。

まずは簡単なLassoからです。

Coordinate descent

いきなり解説に入る前に何ですが、座標降下法ってメチャクチャシンプルなんですよ。
けど全く有名じゃないんですよね。（多分、理由は記事の最後）まあ、アルゴリズムみましょう

目的関数 $f(x) = g(x) + \sum_{i=1}^{n} h_i(x)$ の最適化を考えます。
ここで $g$ はconvexかつdifferentiabl、 $h_i(x)$ はconvexとします。
（例えばconvexかつundifferentiablな $f(x)$ を考えた場合、反例があります。詳しくはREADMORE）
$k=1, 2, \cdots$ 回目の更新を右肩の添え字で表すと、各変数に対して

$x_1^{(k)} = \texttt{argmin } f(x_1, x_2^{(k-1)}, x_3^{(k-1)}, \cdots, x_n^{(k-1)})$

$x_2^{(k)} = \texttt{argmin } f(x_1^{(k-1)}, x_2, x_3^{(k-1)}, \cdots, x_n^{(k-1)})$

$x_3^{(k)} = \texttt{argmin } f(x_1^{(k-1)}, x_2^{(k-1)}, x_3, \cdots, x_n^{(k-1)})$

のように更新するアルゴリズムです。

Pointは１つ、変数を１つずつ更新するんです。その他の変数は定数とみるんです。なので

この図の感じですね。しかし次の問題を考えてみましょう $x \in \mathbb{R}^2$ に対して

$\texttt{Minimize } x_1^2 + x_2^2 + 5|x_1 + x_2|$

です。例えば初期値を $x=(1,-1)$ で座標降下法をアプライしてみると
$x_1,　x_2$ のどちらから更新しても $x=(1,-1)$ で停滞することがわかります。しかし
これは明らかに最適解は $x=(0,0)$ ですよね。こういったことからスパースな問題にはあまり良くないようです。

次にメジャライザーの方を行きましょう。

Iterative Shrinkage Thresholding

ちなみにこのアルゴリズムは初めて聞きました。EMアルゴリズムでも似たことが行われますね。
こちらもアルゴリズムは単純で

テイラー展開を行い、目的関数を置き換えることで最適化問題を解きやすい形にしよう

というものです。

一次元の場合、 $f$ を $x=a$ 周りでのテイラー展開

$f(x) \approx f(a) + f^{(1)}(a)(x-a) + \frac{1}{2!} f^{(2)}(a)(x-a)^2 + \frac{1}{3!} f^{(3)}(a)(x-a)^3 + \cdots$

多次元の場合、

$f(\vc{x})= f(x_1,x_2, \ldots, x_n)$

$f(\bf{x}) \approx f(\bf{a}) + Df(\bf{a}) (\bf{x}-\bf{a}) + \frac{1}{2} (\bf{x}-\bf{a})^T Hf(\bf{a}) (\bf{x}-\bf{a})$

下の図がイメージになります。上界というのは「より大きいやつ」という意味です。

この図からアルゴリズムがどうやって動いて何で最適化できるのかは直感的にわかったかと思います。
しかし、このメジャライザーが本当に存在するのか？など疑問に思う人がいるでしょう。
では考察しましょう。

まず初期値を $x^{(t)}$ とします。そして目的関数を $\mathcal{L}(x)$ ,メジャライザーを $\mathcal{T}(x)$ とします。
そして、メジャライザーは次のような性質を満たしてほしいです。
– $\mathcal{L}(x^{(t)}) = \mathcal{T}(x^{(t)})$
– $\forall x\neq x^{(t)}$ に対して $\mathcal{L}(x) \leq \mathcal{T}(x)$
すると性質から明らかに $T$ の最小化が $L$ の最小化に繋がることがわかります。
では実際にそんな都合のいい $\mathcal{T}(x)$ を探しましょう。の前にちょっと寄り道します。

どんな形のメジャライザーが都合いい？

メジャライザー $T(x)$ を探すといえどこのままでは全くわかりません。
$T(x) = x^2$ にすればいいのか？いや、それもわかりません。なので

都合のいいメジャライザーの形を予想します。

いきなりですが次の最適解 $x$ を考えましょう。次元は１です。( $\lambda \geq 0$ )

$\texttt{argmin } \frac{1}{2} (y-x)^2 + \lambda |x|$

前回の記事でやったようにこれは場合分けで簡単にできますね
$\forall x \geq 0$ に対して

$\frac{1}{2} (x^2 - 2xy + y^2 + 2\lambda x) = \frac{1}{2} (x - (y-\lambda))^2 + \texttt{const}$

同様にして $\forall x < 0$ に対して

$\frac{1}{2} (x - (y+\lambda))^2 + \texttt{const}$

では最後に $y$ の範囲で場合分けすると前回のSoft thresholding functionが出てきましたね。

$x = S ( y, \lambda )$

詳しく復習すると
$y < - \lambda$ の時

$x= y + \lambda$

$- \lambda \leq y \leq \lambda$ の時

$x = 0$

$y > \lambda$ の時

$x = y - \lambda$

これが大切です。

　つまり平方完成して解が求まりました。

ではLassoの目的を再掲しましょう。

$\texttt{argmin } ||y - Xb||_2^2 + \lambda|b|_1^1$

今の例と似てますよね。

　しかしここままでは平方完成できない

ここで

　平方完成できるやつを考えよう！！！

これがほしいメジャライザーを予測するということです。

では元に戻って都合のいい $\mathcal{T}(x)$ を探しましょう。平方完成できるやるがいいんでしたね。
テイラーの２次近似の形を想像しながら

$\mathcal{T}(\mathbf{x}) = \mathcal{L}(\mathbf{x}^{(t)}) + \frac{d\mathcal{L}(\mathbf{x}^{(t)})}{d\mathbf{x}}(\mathbf{x} - \mathbf{x}^{(t)}) + \frac{\rho}{2}\bigl\vert\mathbf{x}-\mathbf{x}^{(t)}\bigr\vert^2$

という形で考えましょう、( $\frac{\rho}{2}$ は計算の都合上のもの)そうすれば

$= \frac{\rho}{2} \Bigl\vert \bigl( \mathbf{x}^{(t)} - \frac{1}{\rho} \frac{d\mathcal{L}(\mathbf{x}^{(t)})}{d\mathbf{x}}\bigr) - \mathbf{x} \Bigr\vert^2 + \texttt{const}$

とできますね。よって

$\text{argmin } \Bigl(\mathcal{T}(\mathbf{x}) + \lambda |\mathbf{x}|\Bigr)$

は

$x^{(t+1)} = S \Bigl( \frac{1}{\rho} \frac{d\mathcal{L}(\mathbf{x}^{(t)})}{d\mathbf{x}}, \lambda \Bigr)$

となります。

　実はまだ終わっていない

そうなんです。先ほどメジャライザーを考えたときに $\rho$ を使いましたね。
この $\rho$ を決める必要があります。メジャライザーの性質（２）から

$\mathcal{T}(\mathbf{x}) - \mathcal{L}(\bf{x}) = \frac{1}{2} (\bf{x} - \bf{x}^{(t)})^\top (\rho I - X^TX) (\bf{x} - \bf{x}^{(t)})$

が「０以上」を任意の $\bf{x}$ について言える必要があるので

　 $(\rho I - X^TX)$ が半正定値

の必要があります。つまり、 $X^TX$ の最大の固有値が $\rho$ として取りうる最小の値になります。

そこで固有値の大きさを見積れる定理があるのでそれを使いましょう。
– Gershgorin circle theorem
– Gershgorin’s Theorem for Estimating Eigenvalues
この定理からmaxをとることで

$\texttt{eigenvalue } \le \max_j \sum_i\bigl\vert{\mathrm{X^\top X}}_{ij}\bigr\vert$

と固有値の大きさが見積もれます。

長かったですがとりあえずLassoはおしまいです。

でわ。

Lassoの進化、Group Lassoとは

2019年2月21日 / コメントする

こんにちは。

Lassoでスパース、スパースと言ってましたが実際はスパース推定という言葉がよく使われます。これについての説明を軽くしてからGroup lassoの紹介をします。

スパース推定とは

スパースは”疎”を意味します。スパース推定とはどのパラメータが0になるかを推定すること。つまり、データの本質がわかります。直感的には次の図から

もう少し実用的な例は？

右の入力から左の出力を考えるLassoモデルを考えると

のようにアレルゲン反応のスパース推定ができます。

つまり、スパース推定により関与しないパラメータがわかります。

Group Lassoとは

グループラッソとは

「潰れる変数がグループになったLassoモデル」

先ほどのアレルゲンを例にとると例えば、ヒノキなどの個体単位ではなく「花粉」というグループ単位でスパース性が検証できるモデルのことです。Lassoはの正規化項は次のものでした。

$\Omega_{Lasso}(\beta) = \sum_{i=1}^{n} |\beta_i|$

Group Lassoの正規化項は次で定義されます

$\Omega_{Group}(\beta) = \sum_{g=1}^{G} \sqrt[]{p_g} |\beta_{g}|_2$

ここで $g$ はg番目のグループを表すindex。(ただし、 $g = 1,\cdots, G$ 、 $p_g$ はグループgの大きさ)
前述通りGroup Lassoでは特徴をグループ化します。よって、事前に類似の傾向がありそうな特徴の情報を考慮します。

なるほど、と思った方と、ん？、と思った方がいると思います。

グループ単位で本当に潰れるのか？

仮にそうならどうやって変数を扱うのか？

これは図を用いてチェックしましょう。

Group Lassoの解

参考文献の論文に従い、解説します。

まず $J$ 個の変数からなるもっとも一般的な回帰問題を考えます

$Y = \sum_{j=1}^{J} X_j \beta_j + \epsilon$

Yは $n\times 1$ ベクトル、 $\epsilon \sim \mathcal{N}(0, \sigma^2I)$ 、 $X_j$ は $j$ 番目のデータに対応した $n\times p_j$ 行列で $\beta_j$ はサイズ $p_j \times 1: j=1,\cdots,J$ の係数からなるベクトルとする。さらに各 $X_j$ は直交行列であると仮定する。
すなはち $X_j^t X_j = I_{p_j}; j=1,\cdots ,J$ とする。さらに、 $X = (X_1, X_2, \cdots, X_J)$ 、 $\beta = (\beta_1^t, \beta_2^t, \cdots, \beta_J^t)^t$ とすると上式は $Y=X\beta + \epsilon$ かける。

長ったらしく書きましたが要は「グループの大きさが各 $p_i$ のサイズ」です。なので $p_1 = p_2 = \cdots = p_J = 1$ を考えるとこれはLassoそのものとなります。

$\eta \in \mathbb{R}^d, d \geq 1$ と $d \times d$ の正定値対称行列 $K$ に対して次を定める。

$|\eta|_K = (\eta^t K \eta)^{1/2}$

ただし、 $|\eta| = |\eta|_{I_d}$ とする。正定値行列 $K_1,\cdots,K_J$ が与えられた時、Group Lasso回帰では次の解を考える。ただし $\lambda \geq 0$

$\frac{1}{2} | Y - \sum_{j=1}^{J} X_j \beta_j |^2 + \lambda \sum_{j=1}^{J} |\beta_j|_{K_j}$

Bakin(1999)はこれをグループ変数によるLassoの拡張版として提案しました。

ここも特に気にする必要はなく、大切なのは正定値対称行列 $K$ により変数に「重み」が掛かっているところです。機械学習ではこのように変数に重みを加える動作をよくします。一例として
– マハラノビス距離(Mahalanobis distance)
があります。僕たちが無意識に距離として使っているものはユークリッドノルムで $K=I$ 、つまり単位行列の時です。

では先ほどの

「どうやってグループ単位で変数が潰れるのか？」

をみましょう。

下の図はグループが2つ(各係数はベクトルとスカラー)、つまり $\beta_1 =(\beta_{11}, \beta_{12})^t, \beta_2$ の場合を考え、 $K_j$ が単位行列の場合の正規化項を表しています（ラッソはダイヤ、リッジは円だったやつの三次元バージョン）

つまり、

Figure(a)は　　　 $||(\beta_1^t,\beta_2) ||_1^1 = 1$
Figure(e)は　　　 $||\beta_1||_K + ||\beta_2||_1^1 = 1$
Figure(i)は　　　 $||(\beta_1^t,\beta_2) ||_2^2 = 1$

もっと噛み砕くと

Figure(a)は　　Lassoの $L_1$ -norm
Figure(e)は　　 GroupLassoのnorm
Figure(i)は　　Ridgeの $L_2$ -norm

です。Figure(e)において $\beta_{11}- \beta_{12}$ 平面を見ると $\beta_1$ というグループ単位でゼロに潰れることがわかります。

だからGroup Lassoはグループ単位ででスパース性を与えているんです。

でわ。

Lassoを数式から実装まで(理論編)

2019年2月19日 / コメントする

こんにちは。

その１でラッソの概要と大きな特徴であるスパース性を確認しました。
今回からはラッソ実装に向けて数学を頑張りましょう。

Notationのについて

Obj (objective function)
OLS (Ordinary least squares) (回帰の残差平方和)
$m$ : データ数
$n$ : 次元(特徴)
$\theta_j$ : $j$ 番目の特徴 $\theta$
$x_j^{(i)}$ : $i$ 番目の観測データ $x$ のうち $j$ に対応するもの

目的関数について

今回、つまりラッソの目的関数は以下の通り

$Obj(\theta) & = OLS(\theta) + \lambda | \theta |_1^1$

$= \frac{1}{2} \sum_{i=1}^m \left[y^{(i)} - \sum_{j=0}^n \theta_j x_j^{(i)}\right]^2 + \lambda \sum_{j=0}^n |\theta_j|$

(OLS)の勾配について

上の目的関数は微分できません。なのでひとまず置いておいて
OLSをいつも通り微分して微分を計算します。
後々のため、今回は勾配でなく、 $\theta_j$ での微分を計算します。

$\frac{\partial }{\partial \theta_j} OLS(\theta) & = - \sum_{i=1}^m x_j^{(i)} \left[y^{(i)} - \sum_{j=0}^n \theta_j x_j^{(i)}\right]$

$= - \sum_{i=1}^m x_j^{(i)} \left[y^{(i)} - \sum_{k \neq j}^n \theta_k x_k^{(i)} - \theta_j x_j^{(i)}\right]$

$= - \sum_{i=1}^m x_j^{(i)} \left[y^{(i)} - \sum_{k \neq j}^n \theta_k x_k^{(i)} \right] + \theta_j \sum_{i=1}^m (x_j^{(i)})^2$

簡単のためこれを次のように表します。

$:= - \rho_j + \theta_j z_j$

L1-normの微分について

$L_1$ -normは微分できないとわかりました。なので
微分できるように、微分という概念を拡張します。そこで出てくるのが

　劣勾配(Subgradient)　

です。wikiのリンク

定義（劣微分）
　Convex $f: I \rightarrow \mathbb{R}$ の点 $x_0$ における劣微分は次の条件を満たす $c \in \mathbb{R}$ で定義する。

(簡単のため $\mathbb{R}^1$ についてはなし、微分の代わりに勾配という言葉を使っている)

$\frac{\partial f }{\partial x}\bigg |<em>{x=x_0} = { c \in \mathbb{R} : f(x) \geq f(x_0) + c(x - x_0) , \forall x \in I }$

書き換えると
$x_0$ における劣勾配は次の閉区間の集合である。 $a \leq b$

$a = \lim</em>{x \rightarrow -x_0} \frac{f((x) - f(x_0)}{x - x_0}$

$b = \lim_{x \rightarrow +x_0} \frac{f((x) - f(x_0)}{x - x_0}$

$\frac{\partial f }{\partial x}\bigg |_{x=x_0} = [a,b]$

これが勾配の拡張になっていることは微分可能な点においてその劣勾配が一点集合になることからわかる

実はこれは簡単で、場合分けして綺麗に微分できるところの間を１つの微分の集合と見るのです。
今回の $L_1$ ノルム、つまり

$f(x) = |x|$

を例にとると、

(1) $\begin{align*} \texttt{if } x > 0 ~~ , ~~ \frac{\partial f}{\partial x} = 1 \end{align*}$

(2) $\begin{align*} \texttt{if } x < 0 ~~ , ~~ \frac{\partial f}{\partial x} = -1 \end{align*}$

そしてこの間の集合を0での微分とするのが劣勾配

(3) $\begin{align*} \texttt{if } x = 0 ~~ , ~~ \frac{\partial f}{\partial x} = [-1,1] \end{align*}$

次のはイメージ図

ではこの劣勾配を使って今回の目的関数の第２項を微分してゼロとおきます。

(4) $\begin{align*} \frac{\partial }{\partial \theta_j} Obj(\theta) = \frac{\partial }{\partial \theta_j} OLS(\theta) + \frac{\partial }{\partial \theta_j} \lambda | \theta |_1^1 \end{align*}$

(5) $\begin{align*} 0 = -\rho_j + \theta_j z_j + \frac{\partial }{\partial \theta_j} \lambda | \theta_j |_1^1 \end{align*}$

ここで先ほどと同様に場合分けを行って

(6) $\begin{align*} \texttt{if } \theta_j < 0 ~~ , ~~ 0 = -\rho_j + \theta_j z_j - \lambda \end{align*}$

(7) $\begin{align*} \texttt{if } \theta_j > 0 ~~ , ~~ 0 = -\rho_j + \theta_j z_j + \lambda \end{align*}$

そして、

(8) $\begin{align*} \texttt{if } \theta_j = 0 ~~ , ~~ 0 = [-\rho_j - \lambda ,-\rho_j + \lambda ] \end{align*}$

ですね。文字が多くなってしまったので再確認ですが、目標は $\theta_j$ を求めることです。上の２つの式は簡単に求まりますが、３つ目の式については閉区間に0が含まれるという条件で考えましょう。なので

(9) $\begin{align*} 0 \in [-\rho_j - \lambda ,-\rho_j + \lambda ] \end{align*}$

(10) $\begin{align*} -\rho_j - \lambda \leq 0 \end{align*}$

(11) $\begin{align*} -\rho_j + \lambda \geq 0 \end{align*}$

なので

(12) $\begin{align*} - \lambda \leq \rho_j \leq \lambda \end{align*}$

まとめると

(13) $\begin{align*} \texttt{for } \rho_j < - \lambda ~~ , ~~ \theta_j = \frac{\rho_j + \lambda}{z_j} \end{align*}$

(14) $\begin{align*} \texttt{for } - \lambda \leq \rho_j \leq \lambda ~~ , ~~ \theta_j = 0 \end{align*}$

(15) $\begin{align*} \texttt{for } \rho_j > \lambda ~~ , ~~ \theta_j = \frac{\rho_j - \lambda}{z_j} \end{align*}$

となります。ちなみにこの場合分けからSoft thresholding functionという関数が定義されます。

最後に劣微分と一緒にplotを見てみましょう。

微分できなくてよく使う関数をたまたま思い出しました。ニューラルネットワークで使うRELU fucntionですね。
– 活性化関数ReLUとその類似関数
他にも色々あるようです。

色々計算づくしでしたが微分を計算できました。あとは最適化のみですが今までは勾配法を使っていましたが「微分不可能な場合は勾配法は無理」なので他の最適化法が必要です。そこでLassoの問題ではいくつかのメジャーな方法があります。

Coordinate descent（座標降下法）
ISTA (メジャライザー最小化)
FISTA (高速化)

です。もっとも簡単なのは座標降下法です。メジャライザーはテイラー展開を用いて上界を扱うものですがこれらについては次回解説します。

でわ。

Lassoはなんでスパース？

2019年2月15日 / コメントする

こんにちは。
素人にfish_shellは無理でした。kzです。

リッジが終わりましたね。ついに来ました

ラッソ

最近色々ラッソについて調べていたんですが、微分不可能な関数の最適化ってやっぱ難しいですね。機械学習において非常に重要な２つのキーワード
– ラグランジュの未定乗数法
– KKT条件
は別の記事でゆっくり解説します。では本題に入りましょう。

Lasso

過学習を考慮した回帰モデルの一つ
– L_1正則化項を使用した回帰model
– スパース性を考えるときに用いる（これについては次の記事で詳しく説明します。）

(1) $\begin{align*} \beta = \texttt{argmin}_b { |y-Xb |^2_2 + \lambda|b|_1^1 } \end{align*}$

リッジ回帰との唯一の違いは正規化項がL1(絶対値)であるということ。

微分できない？

ちょっと微分について復習しよう。おまけで複素解析も出てくるよ

L1とL2

これまでに何度か説明していると思いますがまずは
– 微分可能性
が大きな違いです。L1は尖っているので0で微分できませんね。

他の違いは？

リッジ回帰ではデカイパラメータ $b$ がでかくなりすぎないよう上限を設けましたね。
ラッソも上限はありますが、ゼロがKEYWORDです。なぜなら
– 無関係な特徴量はゼロで排除する
という特徴があります。ゆえに、モデルで初めに設定したものよりも少ない特徴量で済む可能性があります。これが先ほどのスパース性と言われる所以です。

パラメータが0を図から考察

では図を用いて直感的に理解しましょう。次の図はRidge,　Lassoを議論するときに使われるとてもポピュラーなものです。

「スパース」とは「スカスカ」なイメージ

つまり、ゼロがいっぱいあるイメージでいい。
しかし重要なのは常にスパースなわけじゃないこと。

スパースと言われるのはこの「尖ったポイント」に最適解がある時。

これはL1だからできますよね。L2は円だったのでスパース性はありません。

パラメータが0を数式から考察

次のように超シンプルな回帰モデルを考えます。(以下転置使ってますが無くてもいいです)

(2) $\begin{align*} y = bx + \epsilon \end{align*}$

次の最小化を考えます。

(3) $\begin{align*} \texttt{minimize } y^Ty - 2y^T xb + b x^Tx b + 2\lambda |b|_1^1 \end{align*}$

(この2は計算の都合上のもの)

ここで解が $b>0$ と仮定します。

すると、 $y = bx \Leftrightarrow <y,y> = b<y,x>$ より
$y^Tx >0$ と仮定するのと同値であることがわかります。さて、目的関数を $b$ に関して微分しましょう

仮定より今は微分できます

(4) $\begin{align*} \frac{\partial }{\partial b} \texttt{objective function} = -2y^Tx + 2x^Tx b + 2\lambda \end{align*}$

であり、これをゼロと置くことでこの解は次のものだとわかります

(5) $\begin{align*} b = \frac{y^Tx - \lambda}{x^Tx } \end{align*}$

ここで $\lambda$ を増やしていくと $\lambda = y^Tx$ でゼロになることがわかりますね。
この瞬間にLassoはスパース性を持ちます。

一方、リッジの場合は　

(6) $\begin{align*} \textrm{minimize } & y^Ty - 2y^T xb + b x^Tx b + 2\lambda |b|_2^2　 \end{align*}$

(7) $\begin{align*} b &= \frac{y^Tx}{x^Tx + \lambda} \end{align*}$

となるので、 $b \neq 0$ よりパラメータはゼロにはならないですね。

でわ

Ridgeを数式から実装まで(理論)

2019年2月10日 / コメントする

こんにちは。
料理にはまっているkzです。

本日はリッジ回帰行きましょう。次回実装します。🐷

Ridge Regression

過学習を防ごうということで導入された正則化回帰の一例。
– $L_2$ 正則化項を使用した回帰model
– パラメータに上限あり
たまにリッジ回帰はペナルティーありの線形回帰という記事を見かけますがそうではなく多項式回帰でもあります。つまりモデルが

(1) $\begin{align*} y = \beta_0 + \beta_1x_1 + \beta_2x_2 \cdots + \beta_ix_i + \epsilon \end{align*}$

だけでなく

(2) $\begin{align*} y = \beta_0 + \beta_1x + \beta_2x^2 \cdots + \beta_ix^i + \epsilon \end{align*}$

もOKということ。

Notation

Number of data: $N$
Dimension: $n$
Parameters: $b$

つまり各データは

(3) $\begin{align*} x_i := x^{(i)} = (x_1, x_2, \cdots, x_n)^T \end{align*}$

と表されます。 $b$ も同様。

Constrained minimizationは

(4) $\begin{align*} \texttt{minimize } \sum_{i=1}^N (y_i - x_i^Tb)^2 \end{align*}$

(5) $\begin{align*} \texttt{subject to } \sum_{i=1}^n b_i^2 \leq t ~~ , ~~ 0\leq t \end{align*}$

となりますが扱いずらいのでベクトル表示に変更し同時にKKT条件より最適化問題を同値なものに書き換えます。

(6) $\begin{align*} \texttt{minimize } |(y-Xb)|_2^2 + \lambda |b|_2^2 ~~ , ~~ 0\leq \lambda \end{align*}$

ではまず解を求めましょう。いつも通り微分をゼロとおいて解くだけ。今回のコスト関数は次の通り

(7) $\begin{align*} g: b \rightarrow |(y-Xb)|_2^2 + \lambda |b|_2^2 \end{align*}$

として $b$ で微分する。
– 「ベクトルで微分・行列で微分」公式まとめ

(8) $\begin{align*} \frac{\partial g}{\partial b} = -2X^Ty + 2bX^TX+2\lambda \end{align*}$

したがって

$b_{ridge} = (X^TX + \lambda I)^{-1}X^Ty$

$\lambda$ は正則をキープするためのものでしたね。

ここまでは前回の内容をちょっとだけ丁寧にやっただけのおさらいのようなものです。ここからが数学チック。

さて、この解が最適解かどうかを関数の凸性から確認してみましょう。凸性は前回やりましたがお椀型なら極小値が最小値になることでしたね。
– 第2章凸関数
– Convex functions

(9) $\begin{align*} f \in \mathcal{C}^\infty: \mathbb{R}^n \rightarrow \mathbb{R} \end{align*}$

に対して次は同値である。
– $f$ が凸である
– $\forall x\in \mathbb{R}^n$ に対してヘッセ行列 $\nabla^2 f(x)$ がpositive-semi-define

では

(10) $\begin{align*} \frac{\partial^2 g}{\partial b^2} = 2X^TX + 2I_n \end{align*}$

となりますね。ここでこれがpositive-semi-defineであることは次からわかります。 $\forall y\in \mathbb{R}^n$ に対して

(11) $\begin{align*} y^T X^TX y = |Xy|_2^2 \geq 0 \end{align*}$

とノルムの定義よりわかりました。

$g$ はほんまに凸か？

今見たのは

(12) $\begin{align*} H: b \rightarrow |(y-Xb)|_2^2 \end{align*}$

の凸性ですよね。そして

(13) $\begin{align*} K: b \rightarrow \lambda |b|_2^2 \end{align*}$

こちらは明らかに凸。では $g = H + K$ が本当に凸かどうか。これは直感的には明らかですが証明が簡単なのでやってみましょう。

集合 $X$ を定義域とする関数 $H,K$ がともに凸関数であると仮定する。
このとき $\forall x,y \in X ~~ , ~~ \forall \lambda \in [0,1]$ に対して $g$ が凸であることを示す。

(14) $\begin{align*} (H+K)(\lambda x + (1-\lambda)y) = H(\lambda x + (1-\lambda)y) + K(\lambda x + (1-\lambda)y) \end{align*}$

(15) $\begin{align*} \leq \lambda H(x) + (1-\lambda)H(y) + \lambda K(x) + (1-\lambda)K(y) \end{align*}$

(16) $\begin{align*} = \lambda (H+K)(x) + (1-\lambda)(H+K)(y) \end{align*}$

Q.E.D.

Lassoは結構めんどくさいんですよね。。
でわ

線形回帰とRidge回帰を数式からやる(理論編)

2019年2月3日 / コメントする

こんにちは。

正則化パート２ですね。前回の最後に

Ridge Regression (Tikhonov regularization)
Lasso Regression

を例に出しました。今回はこれらを用いて正則化について学びましょう。の前に線形回帰を思い出しましょう。その時に

$\texttt{minimize} | y - X\beta|^2$

を考えて

$X(X^TX)^{-1}X^T$

$X$ の列ベクトルが張る空間への射影行列をやりましたね。ここで $X$ が正則ということなのでもう一歩踏み込みQR分解を行います。つまり

$X = QR ~~ s.t. ~~ Q^TQ = I ~~ , ~~ R~ \textrm{is upper triangular matrix}$

を使うと

$(X^TX)^{-1}X^T = (R^TQ^TQR)^{-1}R^TQ^T = R^{-1}Q^T$

とできます。わかったことは

$X$ の列ベクトルが張る空間への射影行列 $X(X^TX)^{-1}X^T$ が $QQ^T$ となるということ。超シンプル。これが線形回帰。で、これに $L_2$ 正則化項を加えます。つまり $\lambda \in \mathbb{R}_+$ に対して

$\texttt{minimize} | y - X\beta|^2 + \lambda | \beta|_2^2$

を考えます。これをリッジ回帰というのでしたね。ラッソは $L_1$ 正則化項を使います。

違いは？

正則化項のノルムです。 $L_2$ と $L_1$ の違いは？もう少し親しみある言い方だと距離と絶対値の違いは？

Answer: 微分可能性

$L_1$ は微分できませんね。前回の図を思い出してください。 $L_1$ は尖っていました。

gists

勾配法は微分を使いましたね。なので一旦ラッソは無視してリッジについて学びましょう。

次の解を考えます。ただし $\lambda > 0$

$\beta = \texttt{argmin}_b { |y-Xb |^2_2 + \lambda|b|^2_2 }$

「アーグミン $b$ 」は右辺を最小化するときの $b$ ということ。これを展開し、微分すると

$|y-Xb |^2_2 + \lambda|b|^2_2 &= (y-Xb)^T(y-Xb) + \lambda b^Tb$

$= y^Ty -2y^TXb + b^TX^TXb + \lambda b^Tb$

$\frac{\partial}{\partial b} \left( |y-Xb |^2_2 + \lambda|b|^2_2 \right) = -2X^Ty + 2X^TXb + 2\lambda b = 0$

$\Leftrightarrow (X^TX + \lambda I)b = X^Ty$

$\Leftrightarrow \beta = (X^TX + \lambda I)^{-1}X^Ty$

これが欲しかった $\beta$ ですね。さっきの $X(X^TX)^{-1}X^T$ と似てませんか？邪魔なのは $\lambda$ の項。なのでここで $\lambda \downarrow 0$ とするとさっきのと同じになります。だって元々のコスト関数の正則化項が０になるわけですからね。

$\lambda I$ 、これ意味ある？

ここで $X^TX$ は半正定値対称行列であり、実数値からなる対称行列は実数の固有値を持ちます。よって $X^TX$ の固有値は0以上。行列式を計算して見よう。

$\det |X^TX + \lambda I - t I| = 0$

$\det |X^TX - (t - \lambda) I| = 0$

これを $(\lambda -t)$ について解くと $X^TX$ の固有値が求まりますね。そこでその固有値の一つを $e_1$ とするとリッジ回帰における固有値は $t = e_1 + \lambda$ となり、 $\lambda$ 分シフトされていますね。つまり、行列式をゼロから遠ざけており正則をキープしてくれています。