9.4 MARS - 多変量適応回帰スプライン

データから学習するモデルを作るとき、私たちは常に「単純すぎる」と「複雑すぎる」の ジレンマに直面します。直線(線形回帰)は単純で解釈しやすいが、現実の非線形なデータには 当てはまりません。かといって高次多項式は訓練データに過剰に適合してしまい(過学習)、 新しいデータでは失敗します。

この章では、そのジレンマを解決する巧みな手法 MARS(Multivariate Adaptive Regression Splines、 多変量適応回帰スプライン)を紹介します。データが自ら「どこで折れ曲がるべきか」を 教えてくれるという発想が核心です。ヒンジ関数という単純な部品から出発し、 複雑な非線形関係を自動的に発見するアルゴリズムを一緒に見ていきましょう。

ヒンジ関数 - 折れ曲がる直線

データを予測するとき、私たちはしばしばこんな現象に出会います。 「ある閾値を超えると急に影響が出てくる」という現象です。 たとえば、気温が20度を超えると冷房の使用量が急増する、 年収が一定以上になると税率が変わる、などです。

こういった「閾値効果」を表現するのに最適な道具がヒンジ関数(hinge function)です。 まず、アニメーションで見てみましょう。

ヒンジ関数の形状 - ノット位置を境に0から一次関数へと切り替わる
ノット位置(黄色の破線)を境に、左はゼロ、右は直線として立ち上がる

このアニメーションを見てください。ある点(ノット)を超えるまでは 0 のまま、 超えた瞬間から一次関数として立ち上がります。 まるでドアのヒンジ(蝶番)のように、ある点を境に動き出すから「ヒンジ関数」と呼びます。

$$(x - t)_+ = \begin{cases} x - t & \text{if } x > t \\ 0 & \text{otherwise} \end{cases}$$

ここで t がノット(データの折れ曲がり点)です。 この式は次のようにも書けます。

$$(x - t)_+ = \max(0, x - t)$$

$(\cdot)_+$」という記法は「正の部分」を意味します。 つまり、中身が正なら中身のまま、負なら 0 にしてしまいます。

一方、反対向きのヒンジ関数も存在します。 こちらは t を下回るときだけ活性化する、鏡写しの関数です。

$$(t - x)_+ = \begin{cases} t - x & \text{if } x < t \\ 0 & \text{otherwise} \end{cases}$$

この 2 つのヒンジ関数を合わせて「リフレクテッドペア(reflected pair)」と呼びます。 なぜペアにするのか、次のセクションで一緒に考えてみましょう。

リフレクテッドペアの意味

なぜ 2 つのヒンジ関数をペアにするのでしょうか?

考えてみてください。ノット位置 t を境に「左側での効果」と「右側での効果」を 両方捉えるためです。

リフレクテッドペア - 青(右上がり)と赤(左上がり)の2つのヒンジ関数
同じノット位置(黄色)に対して、右向きと左向きの2つのヒンジ関数を用意する

アニメーションを見てください。青い曲線(右上がり)と赤い曲線(左上がり)が、 中央のノット位置を基準に鏡写しになっています。

$(x - t)_+$ はノットの右側での変化を捉え、$(t - x)_+$ はノットの左側での変化を捉えます。

この 2 つを組み合わせることで、モデルはデータを見ながら 「このノット位置では、左に向かう効果が大きいか、右に向かう効果が大きいか」を 自動的に判断できます。

MARS では、すべての入力変数と全観測値に対してリフレクテッドペアを作ります。 これが「候補関数の集合 $\mathcal{C}$」です。

$$\mathcal{C} = \{(X_j - t)_+,\, (t - X_j)_+\}_{t \in \{x_{1j}, \ldots, x_{Nj}\},\, j = 1, \ldots, p}$$

ここで $j$ は変数のインデックス(全 $p$ 変数)、$N$ は観測値の数です。 各変数 $X_j$ の各観測値をノット位置として、 2 方向のヒンジ関数を作ると、合計 2Np 個の候補関数が生まれます。 これがモデルを構築するための「部品箱」です。

MARSモデルの構造

MARSのモデルは次の形をしています。

$$f(X) = \beta_0 + \sum_{m=1}^{M} \beta_m h_m(X)$$

各記号の意味を整理しましょう。

一見すると線形回帰 $f(X) = \beta_0 + \beta_1 X_1 + \beta_2 X_2 + \cdots$ と 同じ形に見えるかもしれません。しかし違いは $h_m(X)$ にあります。線形回帰では $h_m(X) = X_j$ という単純な変数ですが、 MARSでは $h_m(X) = (X_j - t)_+$ のようなヒンジ関数になっています。

複数のヒンジ関数(青・緑・赤)を加算して白い合計曲線(MARSモデル)を作る
個々のヒンジ関数(色付きの曲線)を重み付きで足し合わせると、複雑な形状が生まれる

アニメーションで確認してみましょう。 青・緑・赤の 3 つのヒンジ関数が順に現れ、それらを足し合わせると白い曲線が生まれます。 この白い曲線こそが MARS モデルの予測関数です。

さらに、2 つのヒンジ関数のも使えます。

$$h_m(X) = (X_j - t_1)_+ \cdot (X_k - t_2)_+$$

これは「$X_j$$t_1$ を超え、 かつ $X_k$$t_2$ を超えたときだけ活性化する」 という交互作用を表現できます。詳しくはSection 6 で見ます。

つまり MARS は、折れ曲がる直線(ヒンジ関数)を組み合わせることで、 あらゆる非線形関係を近似できる柔軟なモデルなのです。

前進ステップワイズ - モデルを育てる

MARSはどのようにして最適な項を選ぶのでしょうか?

アルゴリズムは 2 段階で進みます。 まず前進ステップワイズ(Forward Stepwise)でモデルを「育て」、 次に後退削除(Backward Pruning)でモデルを「刈り込む」。

前進ステップワイズでモデルが段階的に成長する様子 - 白い点群に青い曲線がフィットしていく
最初は水平線(定数モデル)から始まり、項を追加するごとにデータに近づいていく

前進ステップワイズの手順:

  1. 最初は定数 $h_0(X) = 1$ だけのモデルから始まる
  2. 候補集合 $\mathcal{C}$ の中から、モデルに追加したとき最も誤差を減らす項を 1 つ選ぶ
  3. 選ばれた項(とその反対向きペア)をモデルに追加する
  4. モデルが十分大きくなるまで繰り返す

アニメーションで見ると、曲線が段階的にデータ点へと近づいていく様子がわかります。 最初は水平な線(平均)から始まり、項が追加されるごとに形が変わっていきます。

この階層的な構造により、本当に重要な項だけが選ばれていきます。 前進ステップワイズは積極的に項を追加するため、 最終的には過学習気味の大きなモデルになります。 これは意図的です。次のステップで「刈り込む」ために。

後退削除とGCV - モデルを刈り込む

前進ステップワイズで大きく育ったモデルを、今度は「刈り込み」ます。

後退削除(Backward Pruning)では:

  1. モデルから 1 項ずつ削除してみる
  2. 削除してもモデルの性能が最も落ちない項を実際に削除する
  3. これを繰り返し、様々なサイズ(項数)のモデルを生成する

でも「どのサイズのモデルが最適か」はどう判断すれば良いでしょうか?

ここで登場するのが GCV(一般化交差検証、Generalized Cross-Validation)です。

$$\text{GCV}(\lambda) = \frac{\sum_{i=1}^{N} (y_i - \hat{f}_\lambda(x_i))^2}{\left(1 - \frac{M(\lambda)}{N}\right)^2}$$

各記号の意味:

分子は予測誤差の二乗和(小さいほど良い)、分母は「モデルの複雑さ」に対するペナルティです。 GCV が最小になるモデルサイズを選びます。

GCV値のU字型カーブ - 最小点(赤い点)で最適モデルサイズが決まる
モデルが単純すぎると予測誤差が大きく、複雑すぎても過学習でGCVが増える

アニメーションを見てください。GCV 値(縦軸)はモデルサイズ(横軸)に対して U 字型を描きます。 最小点(赤い点と黄色の垂直線)が最適なモデルサイズを示しています。

なお、$M(\lambda)$ には係数の数だけでなく選んだノット数にも 3 倍のペナルティがかかります。

$$M(\lambda) = r + 3K$$

ここで $r$ は係数の数、$K$ は選んだノット数です。 ノット位置を「データから探す」という行為は係数を推定するより自由度が多いため、 シミュレーション研究から経験的に 3 倍のペナルティが有効とわかっています。

GCV は通常の交差検証と似た性質を持ちながら、 計算が格段に速いため、大規模なデータでも使いやすい手法です。

交互作用項 - 変数を掛け合わせる

MARS の強力な特徴の一つが、変数間の交互作用を自動検出できることです。

単純な加算モデルでは、各変数の効果が独立に加算されます。

$$f(X) = f_1(X_1) + f_2(X_2) + \cdots$$

しかし現実のデータでは、「$X_1$ が大きいときかつ $X_2$ も大きいとき」だけ 特別な効果が生じることが多いです。これが交互作用です。

2次元格子で右上象限(X1とX2が両方ノットを超える領域)が光るアニメーション
縦・横の破線(ノット位置)より右上の領域だけが活性化する - これが交互作用項

アニメーションで見てください。縦の破線が $X_1$ のノット位置、 横の破線が $X_2$ のノット位置です。 両方の条件を満たす「右上」の領域だけが光ります。これが交互作用項の意味です。

MARS では 2 つのヒンジ関数の積で交互作用を表現します。

$$h_m(X) = (X_j - t_1)_+ \cdot (X_k - t_2)_+$$

この関数は「$X_j > t_1$ かつ $X_k > t_2$ のときだけ活性化する」 という 2 次元的な領域を定義します。

ただし、次数を制限することも重要です。次数を 2(2 変数の積まで)に制限すると、 計算コストを抑えながら主要な交互作用を捉えられます。 次数 1(積なし)に制限すると、各変数が独立に働く加法モデルになります。

複雑な問題では高次の交互作用が必要かもしれませんが、 多くの実用的な問題では次数 2 で十分なことが多いです。

MARSと回帰木の比較

MARS が交互作用まで含む柔軟なモデルを自動構築できることがわかりました。 では、同じく「データ構造を自動発見する」手法として有名な回帰木(CART)と比べるとどうでしょうか?

左:階段状の赤い曲線(CART)と右:滑らかな折れ線の青い曲線(MARS)の比較
同じデータに対して、CARTは「段差」、MARSは「折れ線」でフィットする

アニメーションで比べてみましょう。左側(赤)が CART(回帰木)、右側(青)が MARS です。 同じデータに対して、フィットの仕方が全く異なります。

回帰木は、データ空間を矩形領域に分割し、各領域で定数値を予測します。 境界は「段差」で表現されます。

$$h_m(X) = I(X_j > t) \quad \text{(ステップ関数)}$$

MARSは、折れ曲がる直線でデータを近似します。

$$h_m(X) = (X_j - t)_+ \quad \text{(連続なヒンジ関数)}$$

この違いは重要です。回帰木は境界で「ジャンプ」するため、 滑らかな関係を表現するには多くのノードが必要になります。 MARS は滑らかな曲線を少ない項で近似できます。

面白いことに、MARS のヒンジ関数をステップ関数に置き換えると、 実質的に回帰木のアルゴリズムと同じになることが証明されています。 MARS は回帰木を一般化したものとも言えます。

また、MARS は加算的な構造を保てるため、「どの変数が重要か」の解釈がしやすいという 利点もあります。回帰木の複雑な分岐構造と比べて、MARS のモデルは人間が読みやすいです。

まとめ - MARSの4つのステップ

MARS は単純な部品(ヒンジ関数)から始まり、 データ自身が「どこで折れ曲がるべきか」を教えてくれるアルゴリズムです。

MARSの4ステップ(候補集合→前進→後退→GCV選択)がフロー図で順次光るアニメーション
4つのステップが順番に実行され、最終的な最適モデルが選ばれる

MARSの流れをまとめると:

  1. 候補集合の構築: 全ての変数・全ての観測値に対してヒンジ関数を作る(2Np 個)
  2. 前進ステップワイズ: 誤差を最も減らす項を順次追加してモデルを育てる
  3. 後退削除: 重要度の低い項を順次削除してモデルを刈り込む
  4. GCV で選択: 最適なモデルサイズを選ぶ

MARSの強み:

次の章では、「専門家の混合モデル(Hierarchical Mixtures of Experts)」を学びます。 MARS が「ヒンジ関数でデータ空間を分割する」のに対し、 こちらは「複数の専門家がデータを取り合う」という発想です。 入力によって担当する専門家が変わり、それぞれの専門家が自分の得意な領域を学習します。 MARS が「硬い境界」で空間を分けるのと比べて、専門家の混合モデルは「柔らかい境界」で データを分けるという違いがあります。