Elastic Netの理論と正則化パスの導出

遺伝子発現データを使って疾患の有無を予測するモデルを作りたいとしましょう。測定された遺伝子の数は数万に及ぶ一方、サンプル数はせいぜい数百です。このような $p \gg n$(変数の数がサンプル数をはるかに超える)状況で、LASSOを適用すると困った問題が生じます。LASSOは最大でも $n$ 個の変数しか選択できず、しかも互いに強く相関する遺伝子群がある場合、グループの中から1つだけを恣意的に選び、残りを捨ててしまうのです。

Elastic Netは、リッジ回帰のL2ペナルティとLASSOのL1ペナルティを同時に課すことで、両者の長所を組み合わせた正則化手法です。2005年にZouとHastieによって提案されました。L1ペナルティがスパース性(変数選択能力)を担保し、L2ペナルティが強く相関する変数群を一括で扱う「グループ効果」を提供します。

Elastic Netを理解すると、以下のような幅広い応用が開けます。

  • ゲノミクス・バイオインフォマティクス: $p \gg n$ の超高次元設定で、互いに相関する遺伝子群を同時に選択するバイオマーカー探索
  • テキストマイニング・自然言語処理: TF-IDFベクトルのような高次元スパース特徴量で、共起する単語群を同時に選択する文書分類
  • 金融工学: 相関のあるファクター群のグループ選択やポートフォリオの正則化
  • 気象学・環境科学: 空間的に相関する観測地点の変数を同時に選択する予測モデル

本記事では、LASSOの理論的限界を明確にした上で、Elastic Netの最適化問題を定式化し、ナイーブElastic Netとスケーリング補正の関係、グループ効果の理論的保証を厳密に導出します。さらに座標降下法によるアルゴリズムを解説し、Pythonで正則化パスを可視化してLASSOとの違いを確認します。

本記事の内容

  • LASSOの2つの理論的限界 — $p > n$ 問題とグループ効果の欠如
  • Elastic Netの数学的定式化と制約領域の幾何学
  • ナイーブElastic Netとスケーリング補正の導出
  • グループ効果の理論的保証の証明
  • 座標降下法によるElastic Netの解法
  • Pythonでの実装と正則化パスの可視化

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

LASSOの理論的限界

$p > n$ 問題 — 選択変数数の上限

LASSOの最適化問題を復習しましょう。

$$ \begin{equation} \hat{\bm{\beta}}^{\text{lasso}} = \arg\min_{\bm{\beta}}\left\{\frac{1}{2n}\|\bm{y} – \bm{X}\bm{\beta}\|^2 + \lambda\|\bm{\beta}\|_1\right\} \end{equation} $$

LASSOの解は $\bm{X}^T(\bm{y} – \bm{X}\hat{\bm{\beta}}) = n\lambda\bm{s}$ という最適性条件(KKT条件)を満たします。ここで $\bm{s}$ は劣勾配ベクトルです。非ゼロの係数に対応する列を $\bm{X}_A$ とすると、$\bm{X}_A$ は $n$ 行の行列です。$\bm{X}_A^T\bm{X}_A$ が正則であるためには $\bm{X}_A$ の列数が行数 $n$ を超えてはなりません。

すなわち、LASSOが選択できる変数の数は最大でも $n$ 個です。$p = 20{,}000$ で $n = 200$ のゲノムデータでは、200個を超える遺伝子は選択できません。真に重要な遺伝子が200個を超える場合、LASSOは必ず重要な変数を見逃します。

グループ効果の欠如 — 相関変数の取り扱い

LASSOのもう一つの深刻な問題は、強く相関する変数群の取り扱いです。

具体例で考えましょう。$\bm{x}_1, \bm{x}_2, \bm{x}_3$ が互いに相関係数0.99で相関しており、真のモデルが $y = 2x_1 + 2x_2 + 2x_3 + \varepsilon$ であるとします。LASSOはこの3つの変数のうち1つだけを選び(例えば $\hat{\beta}_1 \approx 6, \hat{\beta}_2 = 0, \hat{\beta}_3 = 0$)、残り2つの係数をゼロにする傾向があります。

これは解釈上も予測上も問題です。解釈の観点では、3つの変数は等しく重要であるにもかかわらず、1つだけが選ばれるのは誤解を招きます。予測の観点では、選ばれた1つの変数にノイズが乗ると予測が大きく不安定になります。

一方、リッジ回帰は相関のある変数群に似た係数を与える「グループ効果」を持っていますが、全ての係数を非ゼロのままにするため変数選択ができません。

Elastic Netは、L1の変数選択能力とL2のグループ効果を統合することで、この二律背反を解消します。具体的にどのような定式化でこれを実現するのか、次のセクションで見ていきましょう。

Elastic Netの数学的定式化

最適化問題

Elastic Netは次の最適化問題の解として定義されます。

$$ \begin{equation} \hat{\bm{\beta}}^{\text{enet}} = \arg\min_{\bm{\beta}}\left\{\frac{1}{2n}\|\bm{y} – \bm{X}\bm{\beta}\|^2 + \lambda P_\alpha(\bm{\beta})\right\} \end{equation} $$

ここでペナルティ関数 $P_\alpha$ はL1ノルムとL2ノルムの凸結合です。

$$ \begin{equation} P_\alpha(\bm{\beta}) = \alpha\|\bm{\beta}\|_1 + \frac{1-\alpha}{2}\|\bm{\beta}\|_2^2 \end{equation} $$

2つのハイパーパラメータが存在します。

  • $\lambda \geq 0$: 正則化の強さを制御するパラメータ。大きいほど係数が強く縮小されます
  • $\alpha \in [0, 1]$: L1とL2の混合比を制御するパラメータ

$\alpha$ の値に応じて、Elastic Netは以下の手法に特殊化されます。

$\alpha$ の値 手法 特徴
$\alpha = 1$ LASSO 純粋なL1ペナルティ。スパース解を与えるがグループ効果なし
$\alpha = 0$ リッジ回帰 純粋なL2ペナルティ。グループ効果はあるが変数選択不可
$0 < \alpha < 1$ Elastic Net L1とL2の混合。変数選択能力とグループ効果を両立

$\alpha = 0.5$ が典型的な選択ですが、データの性質に応じて交差検証で最適な値を選ぶことが推奨されます。

ペナルティの展開形

個々の係数について展開すると、Elastic Netのペナルティは次のように書けます。

$$ \begin{equation} P_\alpha(\bm{\beta}) = \sum_{j=1}^{p}\left(\alpha|\beta_j| + \frac{1-\alpha}{2}\beta_j^2\right) \end{equation} $$

各 $\beta_j$ に対するペナルティ関数 $\alpha|t| + (1-\alpha)t^2/2$ は、原点で角を持つ(L1由来)が、$|t|$ が大きい領域では二次的に増加する(L2由来)ハイブリッドな関数です。L1の角がスパース解を生み出し、L2の二次項が極端に大きい係数を抑制します。

制約領域の幾何学

ペナルティ付き最適化問題は、等価な制約付き最適化問題として表現できます。

$$ \begin{equation} \min_{\bm{\beta}} \frac{1}{2n}\|\bm{y} – \bm{X}\bm{\beta}\|^2 \quad \text{subject to} \quad \alpha\|\bm{\beta}\|_1 + \frac{1-\alpha}{2}\|\bm{\beta}\|_2^2 \leq t \end{equation} $$

この制約領域 $\{\bm{\beta} : \alpha\|\bm{\beta}\|_1 + (1-\alpha)\|\bm{\beta}\|_2^2/2 \leq t\}$ の形状は $\alpha$ に依存します。

  • $\alpha = 1$ のとき、制約領域はL1ボールのダイヤモンド形(正方形を45度回転した形)です。角があるため、等高線がその角に接触しやすく、スパース解(座標軸上の解)が生じます
  • $\alpha = 0$ のとき、制約領域はL2ボールの円形(球形)です。角がないため、解が座標軸上に落ちにくく、全変数が非ゼロになります
  • $0 < \alpha < 1$ のとき、制約領域は丸みを帯びたダイヤモンドの形状になります。L1の角は残っているため座標軸上に解が落ちうる(スパース性)一方、全体の形状は丸みを帯びているため安定性も確保されます

2次元の場合を考えると、制約領域の境界は次の曲線で与えられます。

$$ \begin{equation} \alpha(|\beta_1| + |\beta_2|) + \frac{1-\alpha}{2}(\beta_1^2 + \beta_2^2) = t \end{equation} $$

$\alpha$ が1に近いほど角が尖り、0に近いほど丸くなります。この中間的な形状が、Elastic Netの「適度なスパース性と安定性の両立」を幾何学的に表現しています。

制約領域の幾何学的理解はElastic Netの振る舞いを直感的に把握する上で重要ですが、実際の理論的保証はより厳密な分析が必要です。次に、Elastic Netの解の性質を詳しく調べるために、ナイーブElastic Netとスケーリング補正について解説します。

ナイーブElastic Netとスケーリング補正

問題提起 — 二重の縮小

前節の最適化問題の解をそのまま使う方法をナイーブElastic Netと呼びます。しかし、この解には望ましくない性質があります。

リッジ回帰の解は $\hat{\bm{\beta}}^{\text{ridge}} = (\bm{X}^T\bm{X} + n\lambda\bm{I})^{-1}\bm{X}^T\bm{y}$ であり、L2ペナルティが係数を原点に向けて縮小します。LASSOの解もL1ペナルティにより係数を縮小します。ナイーブElastic Netでは、L1とL2の両方のペナルティによって係数が縮小されるため、二重の縮小(double shrinkage)が起こります。

具体的には、直交デザイン($\bm{X}^T\bm{X}/n = \bm{I}$)の場合、ナイーブElastic Netの解は:

$$ \begin{equation} \hat{\beta}_j^{\text{naive}} = \frac{1}{1 + \lambda(1-\alpha)} S(\hat{\beta}_j^{\text{ols}}, \lambda\alpha) \end{equation} $$

ここで $S(z, \gamma) = \text{sign}(z)(|z| – \gamma)_+$ はソフト閾値関数です。$\hat{\beta}_j^{\text{ols}}$ は最小二乗推定量です。

この式から、ナイーブElastic Netの解は: 1. まずL1ペナルティによりソフト閾値処理(LASSOと同じ縮小)が行われ 2. さらにL2ペナルティにより $1/(1 + \lambda(1-\alpha))$ 倍に縮小されます

LASSOだけでも既にバイアス(縮小バイアス)があるのに、L2ペナルティがさらにバイアスを追加してしまうのです。

スケーリング補正

ZouとHastieは、この二重縮小を補正するために、ナイーブElastic Netの解にスケーリング係数を掛けることを提案しました。

$$ \begin{equation} \hat{\bm{\beta}}^{\text{enet}} = (1 + \lambda(1-\alpha))\hat{\bm{\beta}}^{\text{naive}} \end{equation} $$

この補正を施すと、直交デザインの場合の解は:

$$ \begin{equation} \hat{\beta}_j^{\text{enet}} = S(\hat{\beta}_j^{\text{ols}}, \lambda\alpha) \end{equation} $$

となり、L2による余計な縮小が打ち消されます。補正後の推定量はLASSOと同じ縮小度を持ちつつ、非直交の場合にはグループ効果(後述)を維持するという理想的な性質を持ちます。

補正の導出

なぜ $(1 + \lambda(1-\alpha))$ という係数が適切なのかを見ましょう。ナイーブElastic Netの最適化問題を変形すると:

$$ \frac{1}{2n}\|\bm{y} – \bm{X}\bm{\beta}\|^2 + \lambda\alpha\|\bm{\beta}\|_1 + \frac{\lambda(1-\alpha)}{2}\|\bm{\beta}\|_2^2 $$

L2ペナルティ項を二乗誤差項に吸収させることを考えます。$\tilde{\bm{X}} = \bm{X}/\sqrt{1+\lambda(1-\alpha)}$, $\tilde{\bm{\beta}} = \sqrt{1+\lambda(1-\alpha)}\bm{\beta}$ と変数変換すると:

$$ \frac{1}{2n}\|\bm{y} – \tilde{\bm{X}}\tilde{\bm{\beta}}\|^2 + \frac{\lambda\alpha}{\sqrt{1+\lambda(1-\alpha)}}\|\tilde{\bm{\beta}}\|_1 + \text{const} $$

これは $\tilde{\bm{\beta}}$ に対するLASSO問題に帰着します。$\tilde{\bm{\beta}}$ のLASSO解を $\hat{\tilde{\bm{\beta}}}$ とすると、元の $\bm{\beta}$ のナイーブ解は $\hat{\bm{\beta}}^{\text{naive}} = \hat{\tilde{\bm{\beta}}}/\sqrt{1+\lambda(1-\alpha)}$ です。

LASSOの解 $\hat{\tilde{\bm{\beta}}}$ が適切な縮小度を持つと仮定すれば、元の尺度に戻すためのスケーリング係数は $(1+\lambda(1-\alpha))$ の平方根の逆数ではなく $(1+\lambda(1-\alpha))$ そのもので近似的に補正できます。ZouとHastieはシミュレーションとオラクル不等式の分析から、$(1+\lambda(1-\alpha))$ 倍のスケーリングが実用上最も良い性能を示すことを確認しました。

ここまでの議論で、Elastic Netの定式化と解の修正方法がわかりました。次に、Elastic Netの最も重要な理論的特徴であるグループ効果を厳密に証明しましょう。

グループ効果の理論的保証

グループ効果とは

Elastic Netの最も重要な特性は、相関の高い変数群に似た係数を与えるという「グループ効果」です。直感的には、2つの変数 $\bm{x}_i$ と $\bm{x}_j$ の相関が高ければ高いほど、Elastic Netは $\hat{\beta}_i$ と $\hat{\beta}_j$ に近い値を割り当てます。

定理の述べ方

定理(ZouとHastie, 2005): データは標準化されている($\|\bm{x}_j\| = 1$ for all $j$)と仮定する。ナイーブElastic Netの解 $\hat{\bm{\beta}}^{\text{naive}}$ において、$\hat{\beta}_i \cdot \hat{\beta}_j > 0$(同符号)のとき:

$$ \begin{equation} \frac{1}{\|\bm{y}\|}|\hat{\beta}_i – \hat{\beta}_j| \leq \frac{1}{\lambda(1-\alpha)}\sqrt{2(1 – r_{ij})} \end{equation} $$

ここで $r_{ij} = \bm{x}_i^T\bm{x}_j$ は $\bm{x}_i$ と $\bm{x}_j$ の標本相関係数です。

証明

Elastic Netの最適性条件(KKT条件)から出発します。目的関数を $L(\bm{\beta})$ とすると:

$$ L(\bm{\beta}) = \frac{1}{2n}\|\bm{y} – \bm{X}\bm{\beta}\|^2 + \lambda\alpha\|\bm{\beta}\|_1 + \frac{\lambda(1-\alpha)}{2}\|\bm{\beta}\|_2^2 $$

$\hat{\beta}_i \neq 0$ かつ $\hat{\beta}_j \neq 0$ のとき、KKT条件は:

$$ -\frac{1}{n}\bm{x}_i^T(\bm{y} – \bm{X}\hat{\bm{\beta}}) + \lambda\alpha\text{sign}(\hat{\beta}_i) + \lambda(1-\alpha)\hat{\beta}_i = 0 $$

$$ -\frac{1}{n}\bm{x}_j^T(\bm{y} – \bm{X}\hat{\bm{\beta}}) + \lambda\alpha\text{sign}(\hat{\beta}_j) + \lambda(1-\alpha)\hat{\beta}_j = 0 $$

残差ベクトルを $\bm{r} = \bm{y} – \bm{X}\hat{\bm{\beta}}$ とおきます。上の2式の差を取ると:

$$ -\frac{1}{n}(\bm{x}_i – \bm{x}_j)^T\bm{r} + \lambda\alpha(\text{sign}(\hat{\beta}_i) – \text{sign}(\hat{\beta}_j)) + \lambda(1-\alpha)(\hat{\beta}_i – \hat{\beta}_j) = 0 $$

$\hat{\beta}_i$ と $\hat{\beta}_j$ が同符号のとき、$\text{sign}(\hat{\beta}_i) = \text{sign}(\hat{\beta}_j)$ なので:

$$ \lambda(1-\alpha)(\hat{\beta}_i – \hat{\beta}_j) = \frac{1}{n}(\bm{x}_i – \bm{x}_j)^T\bm{r} $$

両辺の絶対値を取り、コーシー・シュワルツの不等式を適用します。

$$ \lambda(1-\alpha)|\hat{\beta}_i – \hat{\beta}_j| = \frac{1}{n}|(\bm{x}_i – \bm{x}_j)^T\bm{r}| \leq \frac{1}{n}\|\bm{x}_i – \bm{x}_j\|\|\bm{r}\| $$

$\|\bm{x}_i\| = \|\bm{x}_j\| = 1$(標準化の仮定)より:

$$ \|\bm{x}_i – \bm{x}_j\|^2 = \|\bm{x}_i\|^2 – 2\bm{x}_i^T\bm{x}_j + \|\bm{x}_j\|^2 = 2(1 – r_{ij}) $$

したがって $\|\bm{x}_i – \bm{x}_j\| = \sqrt{2(1-r_{ij})}$ です。また $\|\bm{r}\| \leq \|\bm{y}\|$(正則化により残差ノルムはデータのノルム以下)を用いると:

$$ \lambda(1-\alpha)|\hat{\beta}_i – \hat{\beta}_j| \leq \frac{1}{n}\sqrt{2(1-r_{ij})}\|\bm{y}\| $$

$n$ で標準化した形に直すと:

$$ \begin{equation} \frac{1}{\|\bm{y}\|}|\hat{\beta}_i – \hat{\beta}_j| \leq \frac{\sqrt{2(1-r_{ij})}}{n\lambda(1-\alpha)} \end{equation} $$

(注: 元論文では $\|\bm{y}\|$ のスケーリングが若干異なりますが、本質的な結論は同じです。)

定理の解釈

この不等式の意味を考えましょう。

右辺は $r_{ij}$ が1に近づくとゼロに近づきます。つまり、2つの変数の相関が完全に近いほど、Elastic Netの係数の差は小さくなります。$r_{ij} = 1$(完全相関)のとき、$|\hat{\beta}_i – \hat{\beta}_j| = 0$、すなわち2つの係数は完全に等しくなります。

右辺は $\lambda(1-\alpha)$ に反比例します。$\alpha$ が0に近いほど(L2成分が強いほど)、上界が小さくなりグループ効果が強まります。$\alpha = 0$(純粋なリッジ回帰)のとき最もグループ効果が強く、$\alpha = 1$(LASSO)のとき分母がゼロになり上界が発散する、すなわちグループ効果が全くありません。

この定理は、Elastic NetがLASSOとリッジの「いいとこ取り」をする理由を数学的に保証しています。L1成分が変数選択を行い、L2成分が相関のある変数群をまとめて扱うのです。

グループ効果の理論的保証が確認できたので、次にElastic Netを実際に計算するためのアルゴリズムを見ていきましょう。

座標降下法による解法

アルゴリズムの概要

Elastic Netの最適化問題は凸問題であるため、大域最適解が存在し、効率的なアルゴリズムで求めることができます。最も広く使われているのが座標降下法(Coordinate Descent)です。

座標降下法の基本アイデアは単純です。全ての係数を同時に最適化するのではなく、1つの係数 $\beta_j$ のみを変化させて残りを固定し、目的関数を最小化します。これを全ての係数について順番に繰り返し、収束するまで続けます。

1変数の部分問題

$\beta_j$ 以外の係数を固定したとき、$\beta_j$ に関する部分問題は:

$$ \min_{\beta_j}\left\{\frac{1}{2n}\sum_{i=1}^{n}\left(r_{ij} – x_{ij}\beta_j\right)^2 + \lambda\alpha|\beta_j| + \frac{\lambda(1-\alpha)}{2}\beta_j^2\right\} $$

ここで $r_{ij} = y_i – \sum_{k \neq j}x_{ik}\beta_k$ は $\beta_j$ を除いた残差の第 $i$ 成分です。$\tilde{r}_j = \frac{1}{n}\bm{x}_j^T\bm{r}_j = \frac{1}{n}\sum_i x_{ij}r_{ij}$($\bm{x}_j$ が標準化されている場合)とおくと、この部分問題の解は閉じた形で書けます。

$$ \begin{equation} \hat{\beta}_j = \frac{S(\tilde{r}_j, \lambda\alpha)}{1 + \lambda(1-\alpha)} \end{equation} $$

ここで $S(z, \gamma) = \text{sign}(z)(|z| – \gamma)_+$ はソフト閾値関数(soft thresholding operator)です。

この更新式の構造を見ると、Elastic Netの仕組みが明確にわかります。

  • 分子のソフト閾値処理 $S(\tilde{r}_j, \lambda\alpha)$: L1ペナルティの効果。部分残差の相関 $|\tilde{r}_j|$ が閾値 $\lambda\alpha$ を超えなければ $\hat{\beta}_j = 0$ となります(スパース性)
  • 分母の $(1 + \lambda(1-\alpha))$: L2ペナルティの効果。非ゼロの係数を原点に向けてさらに縮小します(安定性)

座標降下法のアルゴリズム

以上をまとめると、Elastic Netの座標降下法は次のアルゴリズムで表されます。

  1. 初期化: $\bm{\beta} = \bm{0}$
  2. 収束するまで以下を繰り返す: – $j = 1, 2, \dots, p$ の各座標について:
    • 部分残差を計算: $\tilde{r}_j = \frac{1}{n}\bm{x}_j^T(\bm{y} – \bm{X}\bm{\beta} + \bm{x}_j\beta_j)$
    • 係数を更新: $\beta_j \leftarrow S(\tilde{r}_j, \lambda\alpha) / (1 + \lambda(1-\alpha))$
  3. 収束判定: $\|\bm{\beta}^{(t+1)} – \bm{\beta}^{(t)}\| < \varepsilon$

このアルゴリズムの計算量は1反復あたり $O(np)$ であり、大規模データにも適用可能です。また、アクティブセット戦略(ゼロでない係数のみを更新する)を組み合わせることで、スパースな場合にはさらに高速化できます。

正則化パスの効率的な計算

実用上は単一の $\lambda$ ではなく、$\lambda$ の値の系列 $\lambda_1 > \lambda_2 > \dots > \lambda_K$ に対して解を計算する「正則化パス」が必要です。ウォームスタート(前の $\lambda$ の解を次の $\lambda$ の初期値として使う)により、パス全体を効率的に計算できます。

最大の $\lambda$ 値(全ての係数がゼロになる値)は次のように計算できます。

$$ \begin{equation} \lambda_{\max} = \frac{\max_j |(\bm{X}^T\bm{y})_j|}{n\alpha} \end{equation} $$

$\lambda_{\max}$ より大きい $\lambda$ では全ての係数がゼロになるため、$\lambda_{\max}$ から始めて徐々に小さくしていけばよいのです。

アルゴリズムの理論的な裏付けが得られたので、Pythonで実装して動作を確認しましょう。

Pythonによる実装と可視化

スクラッチ実装 — 座標降下法

まず、Elastic Netを座標降下法でスクラッチ実装します。

import numpy as np
import matplotlib.pyplot as plt

def soft_threshold(z, gamma):
    """ソフト閾値関数 S(z, gamma) = sign(z)(|z| - gamma)+"""
    return np.sign(z) * np.maximum(np.abs(z) - gamma, 0.0)

def elastic_net_cd(X, y, lam, alpha, max_iter=1000, tol=1e-6):
    """座標降下法によるElastic Netの解法"""
    n, p = X.shape
    beta = np.zeros(p)

    for iteration in range(max_iter):
        beta_old = beta.copy()
        for j in range(p):
            # 部分残差
            r_j = y - X @ beta + X[:, j] * beta[j]
            # x_j との相関
            rho_j = X[:, j] @ r_j / n
            # 座標更新
            beta[j] = soft_threshold(rho_j, lam * alpha) / (1 + lam * (1 - alpha))

        # 収束判定
        if np.max(np.abs(beta - beta_old)) < tol:
            break

    return beta

# --- データ生成: 相関のある変数群 ---
np.random.seed(42)
n = 150
z1 = np.random.randn(n)
z2 = np.random.randn(n)
noise = 0.05  # 相関の強さを制御

# 相関グループ1 (x1, x2, x3): z1 に基づく
# 相関グループ2 (x4, x5): z2 に基づく
# x6-x12: 独立なノイズ変数
X = np.column_stack([
    z1 + noise * np.random.randn(n),
    z1 + noise * np.random.randn(n),
    z1 + noise * np.random.randn(n),
    z2 + noise * np.random.randn(n),
    z2 + noise * np.random.randn(n),
] + [np.random.randn(n) for _ in range(7)])

# 標準化
X = (X - X.mean(axis=0)) / X.std(axis=0)
p = X.shape[1]

# 真の係数
beta_true = np.zeros(p)
beta_true[:3] = [2.0, 2.0, 2.0]
beta_true[3:5] = [-1.5, -1.5]
y = X @ beta_true + 0.5 * np.random.randn(n)
y = y - y.mean()

print("変数間の相関係数(グループ1):")
print(f"  r(x1,x2) = {np.corrcoef(X[:,0], X[:,1])[0,1]:.4f}")
print(f"  r(x1,x3) = {np.corrcoef(X[:,0], X[:,2])[0,1]:.4f}")
print(f"  r(x2,x3) = {np.corrcoef(X[:,1], X[:,2])[0,1]:.4f}")

相関係数の出力から、グループ1の3変数($x_1, x_2, x_3$)が互いに0.99以上の高い相関を持つことが確認できます。この設定で、LASSOとElastic Netの振る舞いがどう異なるかを正則化パスで比較します。

正則化パスの比較

import numpy as np
import matplotlib.pyplot as plt

def soft_threshold(z, gamma):
    return np.sign(z) * np.maximum(np.abs(z) - gamma, 0.0)

def elastic_net_cd(X, y, lam, alpha, max_iter=2000, tol=1e-7):
    n, p = X.shape
    beta = np.zeros(p)
    for iteration in range(max_iter):
        beta_old = beta.copy()
        for j in range(p):
            r_j = y - X @ beta + X[:, j] * beta[j]
            rho_j = X[:, j] @ r_j / n
            beta[j] = soft_threshold(rho_j, lam * alpha) / (1 + lam * (1 - alpha))
        if np.max(np.abs(beta - beta_old)) < tol:
            break
    return beta

# データ生成(前のブロックと同じ)
np.random.seed(42)
n = 150
z1, z2 = np.random.randn(n), np.random.randn(n)
noise = 0.05
X = np.column_stack([
    z1 + noise*np.random.randn(n), z1 + noise*np.random.randn(n),
    z1 + noise*np.random.randn(n), z2 + noise*np.random.randn(n),
    z2 + noise*np.random.randn(n),
] + [np.random.randn(n) for _ in range(7)])
X = (X - X.mean(0)) / X.std(0)
p = X.shape[1]
beta_true = np.zeros(p)
beta_true[:3], beta_true[3:5] = 2.0, -1.5
y = X @ beta_true + 0.5*np.random.randn(n)
y = y - y.mean()

# 正則化パスの計算
lambdas = np.logspace(1, -2.5, 80)
configs = [
    (1.0, "LASSO ($\\alpha=1.0$)"),
    (0.5, "Elastic Net ($\\alpha=0.5$)"),
    (0.1, "Elastic Net ($\\alpha=0.1$)"),
]

fig, axes = plt.subplots(1, 3, figsize=(17, 5))
for idx, (alpha, title) in enumerate(configs):
    coefs = np.array([elastic_net_cd(X, y, lam, alpha) for lam in lambdas])
    ax = axes[idx]
    colors = ['#e41a1c', '#ff7f00', '#984ea3', '#377eb8', '#4daf4a']
    labels = [f'$x_{j+1}$' for j in range(5)]
    for j in range(5):
        ax.plot(np.log10(lambdas), coefs[:, j], '-', color=colors[j],
                linewidth=2, label=labels[j])
    for j in range(5, p):
        ax.plot(np.log10(lambdas), coefs[:, j], '--',
                color='gray', linewidth=0.5, alpha=0.5)
    ax.axhline(0, color='k', linewidth=0.5)
    ax.set_xlabel('$\\log_{10}(\\lambda)$', fontsize=11)
    ax.set_ylabel('Coefficient', fontsize=11)
    ax.set_title(title, fontsize=12)
    ax.grid(True, alpha=0.3)
    if idx == 0:
        ax.legend(fontsize=8, loc='upper left')

plt.tight_layout()
plt.show()

正則化パスの比較から、LASSOとElastic Netの根本的な違いが視覚的に確認できます。

  1. LASSO($\alpha = 1.0$、左図): 相関の高い変数群($x_1, x_2, x_3$)のうち、1つの変数の係数だけが大きく成長し、残りの2つは早い段階でゼロになるか小さな値に留まっています。LASSOは3つの等しく重要な変数から1つだけを「代表」として選んでいます。これがグループ効果の欠如です

  2. Elastic Net $\alpha = 0.5$(中央図): 3つの相関変数の係数パスが互いに近い値を保ちながら同時に成長しています。L2成分のグループ効果により、相関の高い変数に似た係数が割り当てられています。同時に、ノイズ変数(灰色の破線)はゼロに留まっており、L1成分による変数選択も機能しています

  3. Elastic Net $\alpha = 0.1$(右図): L2成分がさらに強くなり、3つの相関変数の係数はほぼ完全に一致しています。グループ効果が最も強く現れていますが、$\alpha$ が小さいためスパース性は弱まり、ノイズ変数の中にも若干の非ゼロ係数が残る場合があります

制約領域の可視化

Elastic Netの制約領域の幾何学を2次元で可視化します。

import numpy as np
import matplotlib.pyplot as plt

fig, axes = plt.subplots(1, 3, figsize=(15, 5))
alphas = [1.0, 0.5, 0.1]
titles = ['$\\alpha=1$ (LASSO)', '$\\alpha=0.5$ (Elastic Net)',
          '$\\alpha=0.1$ (near Ridge)']

for idx, (alpha, title) in enumerate(zip(alphas, titles)):
    ax = axes[idx]
    beta1 = np.linspace(-1.5, 1.5, 500)
    beta2 = np.linspace(-1.5, 1.5, 500)
    B1, B2 = np.meshgrid(beta1, beta2)
    # ペナルティ: alpha * |b1| + |b2|) + (1-alpha)/2 * (b1^2 + b2^2)
    penalty = alpha * (np.abs(B1) + np.abs(B2)) + (1-alpha)/2 * (B1**2 + B2**2)
    ax.contour(B1, B2, penalty, levels=[0.5], colors='blue', linewidths=2)
    ax.contourf(B1, B2, penalty, levels=[0, 0.5], colors=['lightblue'], alpha=0.3)
    ax.axhline(0, color='k', linewidth=0.3)
    ax.axvline(0, color='k', linewidth=0.3)
    ax.set_xlabel('$\\beta_1$', fontsize=12)
    ax.set_ylabel('$\\beta_2$', fontsize=12)
    ax.set_title(title, fontsize=12)
    ax.set_aspect('equal')
    ax.set_xlim(-1.5, 1.5)
    ax.set_ylim(-1.5, 1.5)
    ax.grid(True, alpha=0.2)

plt.tight_layout()
plt.show()

制約領域の形状変化から、$\alpha$ パラメータの役割が幾何学的に理解できます。

  1. $\alpha = 1$(LASSO): 制約領域は正方形を45度回転させたダイヤモンド形です。4つの角が座標軸上にあり、等高線がこの角に接触しやすいため、解が座標軸上(少なくとも1つの係数がゼロ)に落ちやすくなります

  2. $\alpha = 0.5$(Elastic Net): 制約領域は丸みを帯びたダイヤモンドです。角はまだ存在するためスパース解が可能ですが、LASSOほど尖っていないため、すべての係数がゼロでなくなる確率も高くなります

  3. $\alpha = 0.1$(リッジに近い): 制約領域はほぼ円形に近づいています。角が非常に丸いため、解が座標軸上に落ちにくく、多くの係数が非ゼロになります

グループ効果の数値検証

理論で導出したグループ効果の上界を数値的に検証します。

import numpy as np

# 同じデータを使用
np.random.seed(42)
n = 150
z1, z2 = np.random.randn(n), np.random.randn(n)
noise = 0.05
X = np.column_stack([
    z1 + noise*np.random.randn(n), z1 + noise*np.random.randn(n),
    z1 + noise*np.random.randn(n), z2 + noise*np.random.randn(n),
    z2 + noise*np.random.randn(n),
] + [np.random.randn(n) for _ in range(7)])
X = (X - X.mean(0)) / X.std(0)
p = X.shape[1]
beta_true = np.zeros(p)
beta_true[:3], beta_true[3:5] = 2.0, -1.5
y = X @ beta_true + 0.5*np.random.randn(n)
y -= y.mean()

# alpha を変えてグループ効果を検証
print("=== グループ効果の検証 ===")
print(f"r(x1,x2) = {np.corrcoef(X[:,0], X[:,1])[0,1]:.6f}")
print()

lam = 0.1
for alpha in [0.9, 0.5, 0.2]:
    beta = elastic_net_cd(X, y, lam, alpha)
    diff_12 = abs(beta[0] - beta[1])
    r_12 = np.corrcoef(X[:,0], X[:,1])[0,1]
    upper_bound = np.linalg.norm(y) * np.sqrt(2*(1-r_12)) / (n * lam * (1-alpha))

    print(f"alpha={alpha:.1f}: |beta1-beta2| = {diff_12:.6f}, "
          f"upper bound = {upper_bound:.6f}, "
          f"beta1={beta[0]:.4f}, beta2={beta[1]:.4f}, beta3={beta[2]:.4f}")

この数値検証から、以下のことが確認できます。まず、$|\\hat{\\beta}_1 – \\hat{\\beta}_2|$ は理論的上界を常に下回っており、定理の正しさが裏付けられています。次に、$\\alpha$ が小さくなるほど(L2成分が強くなるほど)、$|\\hat{\\beta}_1 – \\hat{\\beta}_2|$ が小さくなり、3つの相関変数の係数がより近い値に収束しています。これがグループ効果の実際の動作です。

まとめ

本記事では、Elastic Netの理論を体系的に解説しました。

  • 動機: LASSOは $p > n$ で選択変数数に制限があり、相関変数群のグループ効果がないという2つの限界を持つ。Elastic NetはL1とL2の混合ペナルティでこれらを克服する
  • 定式化: ペナルティ $P_\alpha(\bm{\beta}) = \alpha\|\bm{\beta}\|_1 + (1-\alpha)\|\bm{\beta}\|_2^2/2$ で、$\alpha$ がL1/L2の混合比を制御。$\alpha = 1$ でLASSO、$\alpha = 0$ でリッジに特殊化される
  • スケーリング補正: ナイーブElastic Netの二重縮小を $(1 + \lambda(1-\alpha))$ 倍のスケーリングで補正する
  • グループ効果: 相関の高い変数に似た係数を与えることが理論的に保証されている。上界は $\sqrt{2(1-r_{ij})}/(\lambda(1-\alpha))$ に比例する
  • 座標降下法: ソフト閾値関数 $S(\tilde{r}_j, \lambda\alpha)/(1+\lambda(1-\alpha))$ による効率的な更新式で、大規模データにも適用可能

Elastic Netはscikit-learnの ElasticNet クラスやglmnetライブラリとして広く実装されており、実務でもすぐに利用できます。ハイパーパラメータ $\lambda$ と $\alpha$ の選択には交差検証(ElasticNetCV)を使うことが推奨されます。

次のステップとして、以下の記事も参考にしてください。