Top 5 This Week

関連記事

24. 高次元回帰:p≫n問題と理論的保証

- 本サイト運営者のサービスの紹介 -

高次元回帰の問題設定とOLSの破綻

統計的回帰分析において、観測数$n$に対して説明変数の数$p$が小さいという設定は古典的な前提です。しかし、ゲノム解析ではSNP(一塩基多型)の数が$p \approx 10^5$〜$10^6$に達する一方でサンプル数は$n \approx 10^3$〜$10^4$程度にとどまります。テキストマイニングでは文書数を大幅に超える語彙次元が生成され、金融リスクモデルでは資産数がサンプル期間を上回ることが通常の状況です。これらを高次元設定と総称し、$p \gg n$で特徴づけます。高次元設定では、低次元理論の基礎となってきた各種仮定が系統的に崩壊します。

通常最小二乗(OLS)推定量の閉形式解は

$$\hat{\beta}^{\mathrm{OLS}} = (X^TX)^{-1}X^Ty$$

と表されます。$X$を$n \times p$の計画行列とすると、$X^TX$は$p \times p$の半正定値行列であり、そのランクは$\min(n, p)$を超えません。$p \geq n$の場合、$X^TX$の行列ランクが高々$n$となり、$p \times p$正方行列として特異行列になります。逆行列は存在しないため、OLS推定量は一意に定義できません。解が存在しないか無数に存在するかのいずれかであり、どちらの場合も統計的推論として意味をなしません。

高次元漸近理論では、次元比$p/n$の挙動として$p = p_n \to \infty$かつ$p/n \to c > 0$(または$c \to \infty$)という枠組みが採用されます。$c > 1$の領域ではOLSが計算不能となり、$c \gg 1$の極限はGWAS(ゲノムワイド関連解析)のような超高次元設定を表します。$c \leq 1$であっても$p$が$n$に対して大きい場合には過適合が顕在化し、訓練誤差がゼロに近くなる一方でテスト誤差が増大します。OLSは$p \geq n$において一意解を持たず、高次元設定には直接適用できません。

高次元設定で統計的推定を成立させるには、パラメータの数を実質的に削減する正則化または次元削減が不可欠です。正則化アプローチはパラメータの大きさに対するペナルティを目的関数に加え、次元削減アプローチは説明変数を低次元の潜在空間に投影します。以降では、正則化によるスパース推定の理論的枠組みを体系的に整理します。

p/n比の増大に伴うOLSとLassoの予測誤差の変化

(Fig1. p/n比の増大に伴うOLSとLassoの予測誤差の変化)

スパース性仮定と制限固有値条件

高次元回帰を統計的に成立させるには、推定対象に構造的制約を課す必要があります。最も標準的な仮定がスパース性仮定であり、真のパラメータベクトル$\beta^* \in \mathbb{R}^p$のうち非ゼロ成分の個数が少ないことを要請します。$s$-スパースベクトルとは

$$\|\beta^*\|_0 \leq s$$

を満たすベクトルを指します。ここで$\|\cdot\|_0$は非ゼロ成分の個数です。$s \ll p$のとき、真の信号は少数の変数に集中し、残りの変数は応答変数に寄与しないという意味を持ちます。ゲノム解析の文脈では疾患関連SNPが全SNP数の中で少数にとどまるという生物学的根拠が、スパース性仮定を支持します。

スパース性仮定だけでは推定理論は構成されません。計画行列$X$に関する条件として、REC(制限固有値条件)が要請されます。非ゼロ支持集合$S \subseteq \{1,\ldots,p\}$($|S| \leq s$)に対してRECは

$$\kappa(s) = \min_{\substack{\|\delta_{S^c}\|_1 \leq 3\|\delta_S\|_1 \\ \delta \neq 0}} \frac{\|X\delta\|_2}{\sqrt{n}\|\delta\|_2} \geq \kappa > 0$$

と定義されます。$\delta_S$は$\delta$の支持集合$S$上の成分、$\delta_{S^c}$はその補集合上の成分を表します。この条件は、スパースな方向$\delta$に対して計画行列$X$が縮退しないこと、すなわち$X$のスパースな方向への射影が一定の伸展性を保つことを保証します。

互換性条件はRECの弱化版であり、Lassoの理論解析で用いられる最小限の条件です。支持集合サイズ$s$に対して

$$\phi^2(s) = \min_{\substack{\|\delta_{S^c}\|_1 \leq 3\|\delta_S\|_1 \\ \delta \neq 0}} \frac{s\|X\delta\|_2^2}{n\|\delta_S\|_1^2} \geq \phi_0^2 > 0$$

と定義されます。REC条件が成立すれば互換性条件も成立しますが逆は成立しません。ガウス分布に従う計画行列では、$n \geq Cs\log p$を満たす十分なサンプルサイズのもとで、これらの条件が高確率で成立することが確率論的に保証されます。

推定理論の成立には、真のパラメータベクトルが$s$-スパース($\|\beta^*\|_0 \leq s$)であること、および計画行列$X$がRECまたは互換性条件を満たすことが必要です。REC条件が満たされない場合(変数間に強い多重共線性が存在する場合など)、以降で述べる理論的保証は失われます。実際の応用では、$X$がRECを満たすかどうかを事前に確認することは一般に困難であり、問題構造から妥当性を判断するか、シミュレーション研究による経験的評価を行う必要があります。スパース性の程度$s$が未知の場合には、実践的な正則化パラメータの選択に困難が生じます。

Lasso推定量の理論的保証:オラクル不等式

Lassoは$\ell_1$ノルムによるペナルティを用いた推定法であり、目的関数は

$$\hat{\beta}^{\mathrm{Lasso}} = \arg\min_{\beta \in \mathbb{R}^p} \frac{1}{2n}\|y – X\beta\|_2^2 + \lambda\|\beta\|_1$$

と定義されます。$\lambda > 0$は正則化パラメータです。$\ell_1$ペナルティの非滑らかな幾何的構造($\ell_1$球の頂点が座標軸上に位置すること)により、最適解においてパラメータの一部が正確にゼロとなる性質があり、変数選択が自動的に行われます。

RECが成立するとき、Lasso推定量$\hat{\beta}^{\mathrm{Lasso}}$の予測誤差についてオラクル不等式

$$\frac{\|X\hat{\beta}^{\mathrm{Lasso}} – X\beta^*\|_2^2}{n} \leq C \cdot s\lambda^2$$

が成立します。定数$C$はREC定数$\kappa$のみに依存します。また推定誤差の上界として

$$\|\hat{\beta}^{\mathrm{Lasso}} – \beta^*\|_1 \leq C’ \cdot s\lambda$$

が得られます。誤差$\varepsilon$が部分ガウス性を持つとき、理論的に最適な正則化パラメータは

$$\lambda \asymp \sqrt{\frac{\log p}{n}}$$

であり、この選択のもとで予測誤差の収束速度は$s\log p / n$となります。$p$が収束速度に対数項でのみ寄与することが重要です。$p$が$n$に対して指数関数的に大きくなっても、$s\log p / n \to 0$を達成するサンプルサイズは$s\log p$のオーダーで十分であり、有効次元が$s\log p$で代表されることを意味します。

「オラクル不等式」という名称は、真のスパース支持集合$S$を事前に知っている仮想的なオラクルが構築するモデルに近い精度が、$\lambda$の適切な選択によって達成されることを指します。スパース回復(正しい変数選択の達成)と予測精度の確保は区別が必要であり、前者には追加の条件が要求されます。収束速度$s\log p / n$は真のスパース次元$s$と$\log p$の積で有効次元を表しており、変数選択後の$s$次元モデルにおける通常の収束速度$s/n$に対して$\log p$の余剰コストがかかることを示します。この余剰コストは$p$個の変数から$s$個を選び出すための探索コストに対応します。

仮定として、誤差$\varepsilon$の部分ガウス性、REC条件の成立、および$s\log p / n \to 0$(有効次元の漸近的小ささ)が要求されます。Lassoは一般に変数選択の符号一致性を保証しません。符号一致性を達成するには、非関連変数と関連変数間の相関を制限するIRREP条件が別途必要であり、これはREC条件より強い条件です。相関の強い変数が存在する場合には推定バイアスが残存します。

適応型Lasso:オラクル性の達成

適応型Lassoは各変数に異なるペナルティ重みを割り当てることで、Lassoが達成できなかった変数選択一致性(オラクル性)を実現します。目的関数は

$$\hat{\beta}^{\mathrm{AL}} = \arg\min_{\beta \in \mathbb{R}^p} \frac{1}{2n}\|y – X\beta\|_2^2 + \lambda\sum_{j=1}^p \hat{w}_j|\beta_j|$$

と定義されます。重みは初期推定量$\hat{\beta}^{\mathrm{init}}$を用いて$\hat{w}_j = 1/|\hat{\beta}_j^{\mathrm{init}}|^\gamma$($\gamma > 0$)と設定されます。初期推定量の絶対値が大きい変数(真の非ゼロ変数の候補)には小さい重みが与えられペナルティが相対的に軽減されます。一方、初期推定量がゼロに近い変数には大きい重みが与えられゼロへの縮小が促進されます。この非対称なペナルティ構造が、Lassoで残存する大係数への過剰バイアスを除去します。

オラクル性は次の2条件として定義されます。第一に、真のゼロ変数が確率1でゼロに推定される変数選択一致性。第二に、真の非ゼロ変数の推定量が$\sqrt{n}$一致性を達成し漸近正規性を持つ推定の漸近効率性です。これらの条件が成立するための正則化パラメータの漸近的条件は

$$\lambda\sqrt{n} \to \infty \quad \text{かつ} \quad \lambda n^{(1-\gamma/2)} \to 0$$

です。第1の条件は真のゼロ変数をゼロに縮小するために十分強いペナルティを要求し、第2の条件は真の非ゼロ変数の推定を妨げないようペナルティを抑制します。両条件を同時に満たす$\lambda$の存在が、適応型Lassoのオラクル性達成の核心です。

二段階の推定手順は、(1)初期推定量の計算($p < n$の場合はOLS、$p \geq n$の場合はRidgeを使用)、(2)重み$\hat{w}_j$の計算、(3)加重Lassoの適用、の順に実施されます。パラメータ$\gamma$は重みの感度を制御し、大きな$\gamma$ほど真のゼロ係数に対応する重みが相対的に増大します。

仮定として、初期推定量の$\sqrt{n}$一致性($p < n$の設定ではOLSがこれを満たし、$p \geq n$ではRidgeが代替として利用されます)と、真の非ゼロ係数の最小絶対値が正であるbeta-min条件が必要です。beta-min条件は検出可能な信号の強さに対する下限を要請しており、非ゼロ係数が実質的にゼロに近い場合には変数選択一致性が成立しません。$p \geq n$ではOLSを初期推定量に使用できないためRidge推定量が代替となりますが、Ridge推定量は真のスパース構造に対してバイアスを持つため二段階設計が複雑化します。高次元設定でのオラクル性の成立にはより強い仮定が要求されます。

サンプルサイズnの増大に伴うLassoと適応型Lassoの変数選択精度の比較

(Fig3. サンプルサイズnの増大に伴うLassoと適応型Lassoの変数選択精度の比較)

ノンコンケーブペナルティ:SCADとMCP

Lassoの$\ell_1$ペナルティは凸関数ですが、係数の絶対値が大きい場合にも一定の収縮バイアスが残存します。ノンコンケーブペナルティは凸性を緩和することで、大きな非ゼロ係数に対するバイアスを排除しながらゼロ係数の縮小を維持する設計が採用されています。

SCAD(Smoothly Clipped Absolute Deviation)のペナルティ導関数はFanとLiによって提案され、

$$p’_\lambda(t) = \lambda \cdot \mathbf{1}(t \leq \lambda) + \frac{(a\lambda – t)_+}{a – 1} \cdot \mathbf{1}(t > \lambda), \quad a > 2$$

と定義されます。$t \leq \lambda$の領域ではLassoと同等の収縮を示し、$\lambda < t < a\lambda$ではペナルティ導関数が線形に減少します。$t \geq a\lambda$でゼロとなり、係数がこの閾値を超えると収縮が完全に停止します。大係数に対してゼロバイアスが達成されることが、SCADの設計上の核心です。

MCP(Minimax Concave Penalty)のペナルティ導関数はZhangが提案し、

$$p’_\lambda(t) = \left(\lambda – \frac{t}{b}\right)_+, \quad b > 0$$

と定義されます。$t = 0$からの収縮量が$t$に対して線形に減少し、$t = b\lambda$でゼロに達します。SCADと異なり導関数が$t$に対して連続的に減衰する設計で、ペナルティの凹性を最小化するという意味でMinimaxと称されます。

SCADとMCPの設計上の違いは、閾値への到達方法にあります。SCADは$[\lambda, a\lambda]$の区間で段階的にペナルティ導関数を減衰させ、MCPは$[0, b\lambda]$の区間で最初から連続的に減衰させます。この設計の違いが収縮特性と局所解の性質に影響し、問題設定によって適性が異なります。

これらのノンコンケーブペナルティが引き起こす非凸最適化には、LQA(局所二次近似)が標準的な手法として用いられます。現在の推定値$\tilde{\beta}$の近傍でペナルティを二次関数で近似し、等価な加重Lasso問題を反復的に解くことでパラメータを更新します。適切な正則化パラメータのもとでは、局所極小解が真のスパース支持集合を回復する局所オラクル性が理論的に保証されます。

仮定として、正則化パラメータ$\lambda$と$a$(SCAD)または$b$(MCP)の適切な選択と、局所オラクル性の成立条件(局所極小解が真のスパース支持集合を回復する条件)の充足が必要です。非凸最適化であるため大域最適解の保証はなく、収束結果が初期値に依存します。実装はLassoより計算コストが高く、チューニングパラメータが多いという実践的課題があります。

L1(Lasso)・SCAD・MCPペナルティ関数の形状比較

(Fig2. L1(Lasso)・SCAD・MCPペナルティ関数の形状比較)

生物統計学への応用:ゲノムワイド関連解析

GWAS(ゲノムワイド関連解析)は高次元回帰の代表的な応用分野です。SNP(一塩基多型)の数$p \approx 10^5$〜$10^6$がサンプル数$n \approx 10^3$〜$10^4$を数十〜数百倍上回る超高次元設定において、疾患リスクに関連するスパースな遺伝子セットの特定が目標となります。疾患関連SNPの数が全SNP数に比べて十分少ないとみなせることが、スパース性仮定の生物学的根拠を与えます。

GWASに固有の困難は、LD(連鎖不平衡)と呼ばれる隣接SNP間の強い相関構造にあります。ゲノム上で物理的に近接するSNP同士は対立遺伝子の組み合わせが独立でなく、計画行列$X$の列間に強い相関が生じます。この多重共線性によりREC条件の成立が妨げられ、Lassoのオラクル不等式が保証する収束速度がGWAS実践では達成されにくくなります。LDブロックの範囲内では多数のSNPが同様の関連シグナルを示すため、スパース回復の精度も理論値を下回ることがあります。

疾患関連SNPの特定には多重検定の問題も伴います。$p \approx 10^6$の設定では偽陽性の制御が不可欠であり、スパース推定による変数選択は多重比較補正と補完的な役割を担います。ただし、スパース推定と多重検定補正は異なる統計的枠組みであり、一方が他方を代替するわけではありません。スパース推定が推定精度を目的とする一方、多重検定補正はタイプIエラーの制御を目的とします。

$\lambda$の実践的な選択にはクロスバリデーションとBIC(ベイズ情報量規準)の両方が利用されます。GWASのような設定ではサンプル分割による統計的検出力の低下が懸念されるため、BICが選好されることがあります。また、グループLassoなどの構造的スパース推定はSNPをゲノム領域でグループ化することでLD構造を明示的にモデルに組み込む拡張として活用されます。

仮定として、疾患関連SNPの数が全SNP数に比べて十分少ないスパース性仮定の成立が必要です。LDにより隣接SNP間に強い相関が生じるため、REC条件が成立しにくい状況が多く、スパース回復の理論的保証が実際のGWAS解析では理論通りに機能しないことがあります。同一LDブロック内の複数SNPが交換可能な推定解として現れる不識別性の問題も生じます。サンプルサイズの拡大がスパース回復精度に直結し、小規模研究では理論保証が実質的に機能しません。大規模コホート研究やメタ解析によるサンプルサイズの拡大が、GWASにおけるスパース推定の有効性を高める実践的な手段となります。

手法比較と実践的選択指針

高次元回帰における主要手法の理論的性質と実践的特性を以下の表に整理します。OLSは参照として掲載しており、$p < n$の設定にのみ適用可能です。

手法 スパース回復一致性 オラクル性 凸性 計算コスト 主な追加仮定
OLS($p < n$のみ) ×(変数選択機能なし) N/A 低(閉形式) $p < n$、計画行列フルランク
Lasso △(IRREP条件が別途必要) × 低〜中 REC条件、部分ガウス誤差
Adaptive Lasso 中(二段階推定) $\sqrt{n}$一致初期推定、beta-min条件
SCAD ○(局所オラクル性) ○(局所) × 高(非凸、LQA反復) 局所オラクル性条件、$\lambda$と$a$の選択
MCP ○(局所オラクル性) ○(局所) × 高(非凸、LQA反復) 局所オラクル性条件、$\lambda$と$b$の選択

$p/n$比と手法適性の観点では、$p < n$かつスパース性仮定が成立する場合はAdaptive LassoまたはSCAD・MCPが変数選択一致性の観点で優位です。Lassoは凸最適化であり計算コストが低いため、$p$が大きくスパース性の程度が不明な場合の探索的推定として有用です。GWASのように$p/n$が極端に大きい設定ではLassoの計算効率が際立ちますが、符号一致性は期待できず、得られた変数選択結果の解釈には注意が必要です。

スパース性仮定が成立しない設定では、これらの手法の理論的保証が適用外となります。変数間に広範な依存構造が存在する場合はリッジ回帰、主成分回帰、部分最小二乗法が代替として検討されます。Elastic Netは$\ell_1$と$\ell_2$ペナルティを組み合わせることで、スパース性仮定が弱い設定での推定安定性を高めます。解釈可能性(変数選択)を重視する場合はLasso系の手法が適し、予測精度を最優先する場合はリッジ回帰やElastic Netが有効な選択肢となります。

$\lambda$の選択はすべての手法に共通する実践的課題です。クロスバリデーションは予測誤差の最小化を目的とするため、変数選択一致性より予測性能が優先されます。変数選択の正確さを重視する場合はBICが選好されることがあります。BICは対数尤度とモデル複雑度の罰則項の和として定義され、BIC最小化はサンプルサイズが大きくなるにつれて真のスパースモデルを一致して選択する性質を持ちます。

これらの理論的保証はいずれも漸近的議論に基づいており、有限サンプルでの挙動が保証と乖離することがあります。特に$n$が小さい場合には、収束速度$s\log p / n$が理論上小さくなっていても実際の推定誤差が大きくなります。また、すべての高次元推定手法においてスパース性仮定が成立しない設定では、理論保証が適用外となる点に留意が必要です。

Popular Articles