13. Lasso回帰：L1正則化とスパース推定

Lasso回帰の定義と目的

Lasso回帰はL1正則化を利用した線形回帰の推定法であり、最小二乗目的関数に係数の絶対値の総和をペナルティとして加えることで、一部の係数を正確に0に縮小するスパース解を生成します。目的関数は次のように定義されます。

$$
\hat{\boldsymbol{\beta}}_{\text{Lasso}} = \arg\min_{\boldsymbol{\beta}} \left\{ \frac{1}{2n} \|\mathbf{y} – \mathbf{X}\boldsymbol{\beta}\|_2^2 + \lambda \|\boldsymbol{\beta}\|_1 \right\}
$$

ここで$\mathbf{y} \in \mathbb{R}^n$は応答変数ベクトル、$\mathbf{X} \in \mathbb{R}^{n \times p}$は設計行列、$\boldsymbol{\beta} \in \mathbb{R}^p$は回帰係数ベクトルです。L1ペナルティは次のように定義されます。

$$
\|\boldsymbol{\beta}\|_1 = \sum_{j=1}^p |\beta_j|
$$

正則化パラメータ$\lambda \geq 0$はペナルティの強さを制御します。$\lambda = 0$のとき通常の最小二乗法と一致し、$\lambda \to \infty$のとき全係数が0に収束します。Ridge回帰がL2ペナルティ$\lambda\sum_{j=1}^p \beta_j^2$を用いて係数を連続的に縮小するのに対し、L1ペナルティは係数を正確に0に縮小する点が異なり、変数選択の効果をもたらします。

前提として、応答変数と予測変数の関係が線形であること（線形性仮定）、および設計行列$\mathbf{X}$が固定された既知の行列であること（設計行列の固定性）が要求されます。

限界として、相関の強い変数が複数存在する場合、その中から任意の1変数を選択し残りを0にするため変数選択の結果が不安定になります。また$n < p$の状況では最大$n$個の変数しか選択できず、解の一意性が保証されません。

凸最適化の観点：KKT条件と部分勾配

Lassoの目的関数は$\beta_j = 0$において微分不可能ですが、L1ノルムは凸関数であり問題全体を凸最適化として扱えます。凸関数$f$の部分勾配の集合（劣微分）は次のように定義されます。

$$
\partial f(\mathbf{x}) = \bigl\{\mathbf{g} \mid f(\mathbf{z}) \geq f(\mathbf{x}) + \mathbf{g}^T(\mathbf{z} – \mathbf{x}),\ \forall \mathbf{z}\bigr\}
$$

絶対値関数の劣微分は次のとおりです。

$$
\partial |\beta_j| = \begin{cases} \{+1\} & (\beta_j > 0) \\ [-1,\, +1] & (\beta_j = 0) \\ \{-1\} & (\beta_j < 0) \end{cases} $$

最適解$\hat{\boldsymbol{\beta}}$における一階最適性条件（KKT条件）は次のとおりです。

$$
\frac{1}{n}\mathbf{X}^T\!\bigl(\mathbf{y} – \mathbf{X}\hat{\boldsymbol{\beta}}\bigr) = \lambda\,\hat{\mathbf{s}}, \quad \hat{s}_j \in \partial|\hat{\beta}_j|
$$

$\hat{\beta}_j \neq 0$の変数については残差との内積が$\pm\lambda$に等しく、$\hat{\beta}_j = 0$の変数については残差との内積の絶対値が$\lambda$以下であることを意味します。

Lagrange双対性の観点から双対問題を定式化します。制約形式の主問題は次のとおりです。

$$
\min_{\boldsymbol{\beta}} \frac{1}{2n}\|\mathbf{y} – \mathbf{X}\boldsymbol{\beta}\|_2^2 \quad \text{s.t.} \quad \|\boldsymbol{\beta}\|_1 \leq t
$$

双対変数$\boldsymbol{u} \in \mathbb{R}^n$を導入すると、Lasso最適化の双対問題は次の形式で表されます。

$$
\max_{\boldsymbol{u}} \left\{ \mathbf{y}^T\boldsymbol{u} – \frac{n}{2}\|\boldsymbol{u}\|_2^2 \right\} \quad \text{s.t.} \quad \|\mathbf{X}^T\boldsymbol{u}\|_\infty \leq \lambda
$$

目的関数が凸かつスレーターの条件が満たされるため強双対性が成立し、双対ギャップは0となります。双対問題の制約$\|\mathbf{X}^T\boldsymbol{u}\|_\infty \leq \lambda$は、各予測変数の残差との相関が$\lambda$以下であるというKKT条件に対応します。

凸性の仮定のもと、局所的最適解は大域的最適解と一致します。限界として、閉形式解が存在しないため反復的な数値計算が必要であり、$p$が大きい場合の計算複雑性が増大します。

スパース性の幾何学的解釈

L1ノルムとL2ノルムはそれぞれ次のように定義されます。

$$
\|\boldsymbol{\beta}\|_1 = \sum_{j=1}^p |\beta_j|, \qquad \|\boldsymbol{\beta}\|_2 = \sqrt{\sum_{j=1}^p \beta_j^2}
$$

制約形式Lassoの幾何学的解釈として、L1球$\{\boldsymbol{\beta} : \|\boldsymbol{\beta}\|_1 \leq t\}$は2次元では菱形、高次元では超多面体となります。最小二乗目的関数の等高線（楕円）がこのL1球に接する点が最適解です。

L1球の境界はコーナー（頂点）を持ち、これらは座標軸上に位置します。等高線楕円がL1球の平らな面よりもコーナーに接触する確率が高いという構造から、最適解が座標軸上（すなわちいくつかの係数が0）に位置しやすくなります。これがスパース解の幾何学的必然性です。L2球は滑らかな球面でありコーナーが存在しないため、同様の仕組みでスパース解は生まれません。

（Fig1. L1球・L2球と最小二乗解の交点：スパース解の幾何学的根拠）

限界として、3変数以上では幾何学的直感を視覚化することが困難になります。2次元での解釈は定性的理解の補助に留まり、高次元での厳密な性質はKKT条件と凸解析によって導出する必要があります。

推定アルゴリズム：座標降下法と軟閾値化

Lassoには閉形式解が存在しないため反復アルゴリズムが用いられます。座標降下法は1回の反復において1つの係数$\beta_j$のみを更新し、他の係数を固定した状態で最適化を繰り返します。$j$番目の係数の更新式は次のとおりです。

$$
\hat{\beta}_j \leftarrow S\!\left(\frac{1}{n}\sum_{i=1}^n x_{ij}\Bigl(y_i – \sum_{k \neq j} x_{ik}\hat{\beta}_k\Bigr),\ \lambda\right)
$$

ここで$S(z,\lambda)$は軟閾値化演算子であり、次のように定義されます。

$$
S(z,\,\lambda) = \mathrm{sign}(z)\cdot\max\!\bigl(|z| – \lambda,\; 0\bigr)
$$

$|z| \leq \lambda$のとき出力は正確に0となり、$|z| > \lambda$のとき入力を$\lambda$だけ0方向に縮小します。この演算子は絶対値関数の劣微分を用いた近接勾配法の解として導出されます。

収束判定条件として、係数の変化量の最大値が閾値$\epsilon$以下になった時点でアルゴリズムを停止します。

$$
\max_j \bigl|\hat{\beta}_j^{(t+1)} – \hat{\beta}_j^{(t)}\bigr| < \epsilon $$

（Fig2. 正則化パス：λに対する係数の軌跡）

前提として、理論的収束性を保証するためにはデータが標準化されていること（各変数の標本平均0・標本標準偏差1）が推奨されます。限界として、初期値の選択によって収束速度が異なります。また$\beta_j = 0$近傍では目的関数が微分不可能であるため、非微分点での数値的な取り扱いに注意が必要です。

正則化パラメータの選択法

$\lambda$はデータから選択します。$K$分割交差検証における検証誤差は次のように定義されます。

$$
\mathrm{CV}(\lambda) = \frac{1}{K}\sum_{k=1}^K \frac{1}{n_k}\left\|\mathbf{y}_{(k)} – \mathbf{X}_{(k)}\hat{\boldsymbol{\beta}}^{(-k)}_\lambda\right\|_2^2
$$

ここで$\hat{\boldsymbol{\beta}}^{(-k)}_\lambda$は$k$番目のフォールドを除いたデータによる推定係数、$\mathbf{y}_{(k)}$と$\mathbf{X}_{(k)}$は$k$番目のフォールドです。一般化情報量規準（GIC）は次のように定義されます。

$$
\mathrm{GIC}(\lambda) = \frac{1}{n}\|\mathbf{y} – \mathbf{X}\hat{\boldsymbol{\beta}}_\lambda\|_2^2 + \frac{a_n}{n}\,|\hat{S}_\lambda|
$$

ここで$|\hat{S}_\lambda|$は選択された変数の個数、$a_n$はサンプルサイズ依存のペナルティ因子です。$a_n = \log n$のとき BIC（ベイズ情報量規準）と対応します。

1標準誤差ルールは、交差検証誤差が最小となる$\lambda_{\min}$の1標準誤差以内で最大の$\lambda$を選択する規則です。

$$
\lambda_{1\mathrm{se}} = \max\bigl\{\lambda : \mathrm{CV}(\lambda) \leq \mathrm{CV}(\lambda_{\min}) + \hat{\sigma}_{\mathrm{CV}}\bigr\}
$$

$\hat{\sigma}_{\mathrm{CV}}$は交差検証誤差の標準誤差です。$\lambda_{1\mathrm{se}}$は$\lambda_{\min}$より大きく、より少ない変数を選択するスパースなモデルを返します。

（Fig3. 交差検証曲線と最適正則化パラメータの選択）

前提として、データが$K$分割可能な十分なサンプルサイズを持つことが必要です。限界として、計算量はグリッドサイズと交差検証分割数の積に比例して増大します。小サンプルでは各フォールドのデータが少なくなり、交差検証誤差の推定が不安定になります。

生物統計での実践応用：遺伝子選択

マイクロアレイや次世代シーケンサーによる遺伝子発現データは、患者数（サンプル数）$n$が数十から数百程度であるのに対し、遺伝子数$p$が数千から数万に達する高次元データです。この$p \gg n$状況においてOLSは適用不可能ですが、Lassoは最大$n$個の変数を選択するスパース解を与えます。

疾患の予後予測や診断モデルの構築では、多数の遺伝子の中から疾患と関連する少数の遺伝子バイオマーカーを特定することが求められます。Lassoによる特徴選択は、そのようなスパースな生物学的シグナル構造に対応する推定を可能にします。

前提として、疾患に関連する遺伝子の数が全体の中でごく少数であるというスパース性の仮定が設定されます。疾患が多数の遺伝子の微小な寄与の集積によって生じる場合、この仮定は成立しません。

限界として、Lassoによる変数選択は正則化の枠組みに基づくものであり、古典的な有意性検定の枠組みとは相容れません。選択された遺伝子に対して従来のp値解釈を直接適用することはできず、二者の融合は現在も方法論的に困難な課題です。またサンプルサイズが不足する状況では選択される変数が不安定になり、異なるサンプルで異なる遺伝子が選択される可能性があります。

他の回帰手法との比較

OLS、Ridge回帰、Elastic Net、前進的変数選択とLassoを比較します。L0正則化は変数の個数そのものをペナルティとする非凸問題であり、L1正則化はそのL0正則化の凸緩和として位置付けられます。

（Fig4. Lasso・Ridge・OLSのバイアス・分散曲線）

正則化パラメータを増加させると、いずれの手法でもバイアスが増加して分散が低下します。Lassoはある閾値を超えた係数を段階的に0に縮小するため、Ridgeよりも急峻なバイアス増加と分散低下が生じます。バイアス・分散トレードオフの観点から、最適な$\lambda$は両者のバランス点として$\mathrm{MSE} = \mathrm{Bias}^2 + \mathrm{Var}$を最小化する値となります。

Ridge回帰はスパース性を持たず全変数を保持するため、変数解釈が必要な場面では不向きです。LassoはRidgeに比べて相関変数への振る舞いが不安定であり、相関の強い変数群から選択される変数がサンプルに依存します。Elastic NetはL1とL2ペナルティを組み合わせることでこの不安定性を緩和しグループとして変数を選択できますが、正則化パラメータが2つになることで最適化の探索空間が増大します。