木材近赤外スペクトルで読む PRML 3.5-3.6 解説ノート

まずコンペとして何を解くのか

このノートでは、PRML 3.5-3.6 の数式を、木材の近赤外スペクトルから含水率を予測するコンペを主な題材として読む。最初に実問題を押さえ、そのあとで「なぜエビデンス近似が必要になるのか」「なぜ固定基底関数には限界があるのか」を PRML の式に接続する。

この題材は、近赤外スペクトル（近赤外光を当て、波数ごとの吸光度を並べた数値ベクトル）から、木材の含水率（木材に含まれる水の量を表す割合）を予測する回帰課題（連続値を予測する課題）である。評価指標は RMSE（予測値と正解値のずれを二乗して平均し、その平方根を取った指標）であり、小さいほどよい。

扱っている対象は、次のような 19 樹種の木材である。

木材は樹種（木の種類）によって色、木目、密度（同じ体積あたりの重さ）、組織構造（木材内部の細胞や繊維の並び）が大きく異なる。そのため、木材の含水率をスペクトルから予測するには、水に由来する情報と、樹種や組織構造に由来する情報を区別する必要がある。

測定では、木材を 20 mm x 20 mm x 12 mm の試料に切り出し、飽水状態（水を十分に含ませた状態）から室温で乾燥させながら、重量と近赤外スペクトルを繰り返し測る。

含水率は次の式で定義される。

$$\text{含水率}(\%) = \frac{\text{水の質量}}{\text{全乾状態の木材質量}} \times 100$$

予測入力は画像ではなく、CSV（表形式のデータを保存するテキストファイル）に入った波数ごとの吸光度である。水は近赤外領域に特徴的な吸収を持つため、スペクトルには含水率に関する情報が含まれる。

データを見て分かる難しさ

EDA（探索的データ分析。モデルを作る前に、分布、外れ値、変数間の関係を確認する作業）から、この課題の構造が見える。基本構造は次の通りである。

train（学習用データ）と test（評価用データ）のスペクトル列は一致しており、波数はおよそ 9994 cm^-1 から 4000 cm^-1 まで降順に並ぶ。つまり、入力 $\mathbf{x}$ は 1555 次元のスペクトルベクトルである。ここから、まず 1 つの試料を表す入力が高次元であるという難しさが分かる。

目的変数である含水率は、右側に長い裾を持つ偏った分布である。

中央値は 29.05%、平均は 49.94%、最大値は 298.58% である。75%点は 70.26%、90%点は 124.89% であり、高含水率側の値がかなり大きい。RMSE は大きな誤差に強く反応するので、高含水率側を外すとスコアに大きく効く。

含水率の分布は樹種ごとにも大きく異なる。

train は 13 樹種、test は 6 樹種で構成され、両者の樹種は重なっていない。したがって、この課題は単に「見たことのある樹種の測定値を補間する問題」ではなく、未知樹種に対して含水率を予測する問題である。

この事実は、ランダムにデータを分ける検証が危険であることも示す。ランダム分割では、同じ樹種の似た乾燥系列が訓練側と検証側に入る可能性がある。その場合、検証スコアは良く見えるが、未知樹種に対する性能を過大評価しやすい。したがって、この課題では LOGO-CV（Leave-One-Species-Out Cross Validation。1 つの樹種を丸ごと検証用にし、残りの樹種で学習する交差検証）が重要になる。

スペクトルを見ると、低含水率群と高含水率群では平均スペクトルに明確な差がある。

特に水分に関係する波数帯（ある範囲の波数）では差が大きい。train 内で含水率との相関（片方が増えるともう片方も増減しやすい関係）が最も高い波数は 7618 cm^-1 付近で、相関は約 0.779 である。

固定した水分帯特徴でも、mean_7500_4800、mean_7200_6200、mean_5400_4900 などは含水率と強い相関を持つ。これは、近赤外スペクトルに含水率の情報が確かに入っていることを示す。ただし、相関が高い波数を単純に選ぶだけでは、樹種差や散乱の影響まで拾う可能性がある。

つまり、1555 次元のスペクトルには水分の情報が含まれている一方で、水分だけでなく、樹種差、木材組織、散乱（光が木材内部や表面で方向を変える現象）、表面状態、測定位置の違いも混ざっている。この多数の変動の中から、未知樹種にも通用する、含水率に本当に関係する安定した変動を取り出すことが、このデータを見て分かる中心的な難しさである。

PCA（主成分分析。高次元データを、変動が大きい少数の方向に要約する方法）の結果も重要である。raw スペクトル（前処理していない生のスペクトル）を標準化（平均やばらつきをそろえる処理）して PCA すると、PC1（第 1 主成分、最も大きな変動方向）だけで分散の 94.65% を説明する。これは、raw スペクトルの大部分が全体強度やベースライン（スペクトル全体の底上げや傾き）のような大きな変動に支配されている可能性を示す。SNV 前処理（各スペクトルごとに平均とばらつきをそろえる前処理）後では、PC1 は 51.08%、PC2 は 24.01%、PC3 は 17.91% となり、構造が複数の方向に分かれる。

この図は、train と test の分布が完全には同じでないことを示す確認用の図である。ただし、test 全体の分布を使って補正や特徴量作成を行うと、評価データへの依存が生じる。したがって、test 分布は sanity check（明らかにおかしな点がないかを確認すること）として見るだけにし、前処理やモデル選択は train 内で完結させる必要がある。

EDA では、train のみを使った LOGO-CV の簡易ベースラインとして PLS（Partial Least Squares。目的変数と関係する方向を優先して次元を減らす回帰法）も計算されている。上位は次の通りである。

pooled RMSE は全検証サンプルをまとめて計算した RMSE であり、worst 樹種 RMSE は最も悪かった樹種の RMSE である。この結果から、SNV + 低成分 PLS は最初の安定した候補になる。一方で、raw + PLS も近い性能を示しており、raw の強度情報を完全に捨てるべきとは限らない。また、worst 樹種 RMSE が大きいため、平均スコアだけで判断せず、どの樹種で大きく外すかを見る必要がある。

ここまでで見えている問題は、PRML 3.5 の問題意識とよく対応する。モデルは訓練データに合わせるだけでは不十分であり、どの程度の複雑さを許すか、どの方向の重みを信用するか、ノイズや偶然の変動にどれだけ合わせないかを決める必要がある。

検証用語の整理

この題材では、単に「CV スコアが良い」と言うだけでは不十分である。train と test で樹種が重ならないため、検証では「見たことのない樹種に対してどれくらい外さないか」を確認する必要がある。ここで使う用語を先に整理しておく。

CV（Cross Validation、交差検証）は、train データをいくつかの分割に分け、一部を検証用、残りを学習用にして性能を測る方法である。1 回の分割を fold と呼ぶ。重要なのは、検証用 fold のデータを学習に使わないことである。これにより、手元の train データだけで未知データへの性能を近似的に見る。

OOF（Out-of-Fold prediction）は、各サンプルを「そのサンプルを学習に使っていないモデル」で予測した値である。たとえば LOGO-CV である樹種を検証用にした fold では、その樹種の全サンプルに対する予測は OOF 予測になる。全 fold の OOF 予測をつなげると、train 全体に対して「学習に使わずに予測した値」がそろう。これを使って計算する RMSE が pooled OOF RMSE である。

LOGO-CV（Leave-One-Species-Out Cross Validation）は、1 つの樹種を丸ごと検証用にし、残りの樹種で学習する検証である。13 樹種の train なら 13 fold になる。未知樹種への外挿を直接見るため、この題材の基本検証である。

LKSO（Leave-K-Species-Out）は、1 樹種ではなく K 個の樹種を同時に検証用にする方法である。LOGO-CV より厳しく、複数の未知樹種が同時に来たときにモデルが崩れないかを見る stress test として使う。K が大きいほど、学習側に残る樹種が少なくなり、検証は難しくなる。

Species Cluster CV は、train 内の樹種を、含水率分布や平均スペクトルが似ているもの同士でクラスタにまとめ、そのクラスタ単位で検証用にする方法である。ランダムに樹種を抜くのではなく、似た性質の樹種群をまとめて抜くことで、分布シフトに対してモデルが頑健かを見る。

Species GroupKFold は、樹種をグループとして扱い、同じ樹種が学習側と検証側に同時に入らないようにする K-fold 検証である。LOGO-CV や LKSO の結論が、特定の分割方法だけに依存していないかを確認する補助的な検証として使う。

LOGO と Strict LOGO の違いも重要である。LOGO は「1 樹種を抜く」という分割方法そのものを指す。一方、Strict LOGO は、その LOGO 結果を使った採用判定である。平均スコアだけでなく、既存の基準モデルに対して何 fold 改善したか、悪化した fold の悪化幅が大きすぎないか、bootstrap CI（fold 差分を再標本化して見積もる信頼区間）の上側が改善側にあるか、などを確認する。つまり、LOGO は測り方、Strict LOGO は「局所的に壊れていないか」を見る厳しめの gate である。

Multi-CV gate は、LOGO、LKSO、Species GroupKFold、Species Cluster CV など複数の検証をまとめた採用判定である。複数の CV を単純に平均して 1 つの点数にするのではなく、どれか 1 つで大きく壊れていないかを確認する。Public LB は見える一部評価なので、最終汎化性能そのものとは限らない。したがって、提出候補は Public LB だけで選ばず、Multi-CV gate を通るかを重視する。

このノートで出てくる主な検証指標は次の通りである。

PRMLの記号に翻訳する

この題材と PRML の記号は、次のように対応させて読むと分かりやすい。

この題材ではスペクトル列が 1555 個ある。しかし、近い波数の吸光度は強く相関し、水分に関係する波数帯も限られる。そのため、形式的な入力次元は大きくても、含水率予測に実際に使える独立した情報方向はずっと少ない可能性がある。この「見かけのパラメータ数」と「実際にデータで決まる自由度」の違いを表すのが、PRML 3.5.3 の $\gamma$ である。

この節で何をしたいのか

線形基底関数モデル（入力を基底関数で特徴量に変換し、その特徴量の重み付き和で予測するモデル）では、入力 $\mathbf{x}$ を基底関数 $\boldsymbol{\phi}(\mathbf{x})$で特徴量に変換し、それらの線形結合で予測を行う。この題材で言えば、raw スペクトル、SNV 後スペクトル、水分帯特徴、PLS 成分などを $\boldsymbol{\phi}(\mathbf{x})$ とみなし、それらの重み付き和で含水率を予測するという見方である。

このモデルには、少なくとも次の 2 つの重要なハイパーパラメータ（学習で直接求める重みとは別に、モデルの性質を決める設定値）がある。

• $\alpha$ : 重み $\mathbf{w}$ の事前分布（データを見る前に重みがどの範囲にありそうかを表す分布）の精度（分散の逆数）。大きいほど重みを 0 付近に強く押し込み、モデルを単純にする。
• $\beta$ : 観測ノイズ（測定値や目標値に含まれる誤差）の精度。大きいほどノイズ分散は小さいとみなす。

通常の正則化付き最小二乗では、正則化の強さを交差検証などで決めることが多い。この節では、ベイズ的な考え方（不確実性を確率分布で表す考え方）を使い、訓練データだけから $\alpha$ と $\beta$ を決める方法を考える。

その中心になる量が、周辺尤度 $p(\mathbf{t}|\alpha,\beta)$ である。これは、ある $\alpha,\beta$ のもとで、観測された訓練データ $\mathbf{t}$ がどれくらい自然に説明されるかを表す。この周辺尤度をエビデンスと呼ぶ。

正則化付き最小二乗との対応で見ると、比 $\alpha/\beta$ は正則化パラメータ（正則化の強さを決める値）に似た役割を持つ。したがって、この節は「木材近赤外スペクトルのような高次元スペクトル回帰で、正則化の強さをデータからベイズ的に決める方法」として読むことができる。

記号の整理

上の対応表に出てこなかった、式を追うための補助記号を整理する。

• $N$ : 訓練データ点の数。
• $M$ : 重みベクトル $\mathbf{w}$ の次元、つまりパラメータ数。
• $\mathbf{m}_N$ : 事後分布（データを見たあとに更新された分布）における重みの平均、またはモード（分布が最も高い位置）。
• $\mathbf{S}_N$ : 事後分布における重みの共分散行列（複数の重みが一緒にどう変動するかを表す行列）。
• $\mathbf{A}$ : 誤差関数のヘッセ行列（2 階微分を並べた行列）。
• $\lambda_i$ : データが第 $i$ 方向の重みをどれくらい強く制約しているかを表す固有値。
• $E_D(\mathbf{w})$ : データへの当てはまりの悪さを表す誤差項。
• $E_W(\mathbf{w})$ : 重みの大きさを表す正則化項。

精度は分散の逆数である。したがって $\beta$ が大きいほどノイズ分散は小さく、 $\alpha$ が大きいほど重みの事前分布は狭くなる。

この図は、$\alpha$ と $\beta$ がそれぞれ何を制御しているかを示している。$\alpha$ は重みの大きさをどれくらい抑えるかを決める量であり、大きいほど重みは 0 付近に集まり、モデルは単純になる。この題材では、樹種差やノイズに由来する不安定なスペクトル方向へ大きな係数を与えすぎないための調整に対応する。$\beta$ は観測ノイズをどれくらい小さいとみなすかを決める量であり、大きいほどデータ点のまわりのばらつきは小さいと考える。この 2 つをデータから決めることが、この節の中心的な目的である。

3.5 エビデンス近似

この題材で特徴量や PLS 成分数を増やすと、訓練データへの当てはまりは改善しやすい。しかし、それが未知樹種でも効くとは限らない。この節の問いは、重みの不確実性と正則化の強さを確率的に扱い、訓練データが本当に支持する複雑さを選べないかというものである。

この図は、完全なベイズ推論からエビデンス近似へ進む流れを先取りして示している。左側は $\mathbf{w},\alpha,\beta$ をすべて平均する考え方、右側は $\alpha,\beta$ をデータが最も支持する値に固定し、$\mathbf{w}$ だけを平均する考え方である。以下で、この違いを式で確認していく。

完全なベイズ推論は何をするか

完全にベイズ的に扱うなら、未知の重み $\mathbf{w}$ だけでなく、ハイパーパラメータ $\alpha$ と $\beta$ についても不確実性を持たせる。つまり、ある 1 つの値に固定するのではなく、それらの可能な値すべてを平均する。

そのとき予測分布は次のようになる。

$$p(t|\mathbf{t}) = \iiint p(t|\mathbf{w}, \beta)p(\mathbf{w}|\mathbf{t}, \alpha, \beta)p(\alpha, \beta|\mathbf{t}) \,d\mathbf{w} \,d\alpha \,d\beta$$

(3.74)

この式は、次の 3 つを掛けて、 $\mathbf{w}$ 、 $\alpha$ 、 $\beta$ について積分している。

• $p(t|\mathbf{w}, \beta)$ : 重みとノイズ精度が決まったとき、新しい目標値 $t$ がどう分布するか。
• $p(\mathbf{w}|\mathbf{t}, \alpha, \beta)$ : 訓練データとハイパーパラメータが決まったとき、重みがどう分布するか。
• $p(\alpha, \beta|\mathbf{t})$ : 訓練データを見たあと、ハイパーパラメータがどう分布するか。

直感的には、「重みも、正則化の強さも、ノイズの大きさも、本当は分からないのだから、全部の可能性を平均しよう」という式である。

ただし、この完全な積分は一般には解析的に解けない。そこで、 $\alpha$ と $\beta$ の事後分布がある値の周辺で鋭く集中していると仮定する。その代表値を $\widehat{\alpha}$ と $\widehat{\beta}$ とおく。

すると、$\alpha$ と $\beta$ は代表値に固定し、$\mathbf{w}$ だけを周辺化（積分して平均することで、その変数を明示的に扱わなくすること）すればよい、という近似になる。

$$p(t|\mathbf{t}) \simeq p(t|\mathbf{t}, \widehat{\alpha}, \widehat{\beta}) = \int p(t|\mathbf{w}, \widehat{\beta}) p(\mathbf{w}|\mathbf{t}, \widehat{\alpha}, \widehat{\beta}) \, d\mathbf{w}.$$

(3.75)

この近似がエビデンス近似の考え方である。

この題材で読むと、完全なベイズ推論は「スペクトルのどの方向をどれだけ使うか」という重み $\mathbf{w}$ だけでなく、「どれくらい重みを抑えるか」という $\alpha$ と、「残差をどれくらいノイズとみなすか」という $\beta$ まで平均する考え方である。一方、エビデンス近似では、$\alpha,\beta$ は訓練データが最も支持する代表値に固定し、そのうえで $\mathbf{w}$ の不確実性を残す。これは、コンペで言えば、正則化の強さやノイズ水準を毎回手で決め打ちするのではなく、含水率データ自身から自然な複雑さを選ばせる発想である。

どの $\alpha,\beta$ を選ぶか

ベイズの定理（条件付き確率を逆向きに計算する公式）から、ハイパーパラメータの事後分布は次の比例関係（定数倍を無視すれば等しい関係）で書ける。

$$p(\alpha, \beta | \mathbf{t}) \propto p(\mathbf{t} | \alpha, \beta) p(\alpha, \beta).$$

(3.76)

ここで重要なのは $p(\mathbf{t} | \alpha, \beta)$ である。これは、$\alpha,\beta$ を固定したときに訓練データ $\mathbf{t}$ がどれくらい起こりやすいかを表す。この量を周辺尤度、またはエビデンスと呼ぶ。

もし $p(\alpha,\beta)$ が比較的平坦な事前分布なら、$p(\alpha, \beta | \mathbf{t})$ を最大にすることは、ほぼ $p(\mathbf{t} | \alpha, \beta)$ を最大にすることと同じになる。つまり、データ自身が最もよく支持する $\alpha,\beta$ を選ぶ。

この方法は、統計学では経験ベイズ（ハイパーパラメータをデータから推定するベイズ手法）、タイプ 2 最尤法（重みを積分した後の尤度を最大化する方法）、一般化最尤法などと呼ばれ、機械学習ではエビデンス近似と呼ばれる。

ただし、ここでの「データ自身」は train の含水率データである。このコンペでは test 樹種が未知なので、エビデンス最大化だけで未知樹種への汎化を保証できるわけではない。したがって、エビデンスはモデル内部の正則化を決める視点として役立ち、LOGO-CV や LKSO は未知樹種への外挿を確認する視点として役立つ、と分けて考える必要がある。

なぜ完全な周辺化ではなくエビデンス最大化を使うのか

$\alpha$ と $\beta$ にガンマ事前分布（正の値だけを取るパラメータに使いやすい事前分布）を置けば、ハイパーパラメータについての周辺化の一部は解析的（数値計算ではなく式変形で扱えること）にできる。その結果、重み $\mathbf{w}$ 上に Student の t 分布（ガウス分布より裾が厚い分布）が出てくる。しかし、そこからさらに $\mathbf{w}$ について積分する部分は、一般には解析的に扱えない。

ラプラス近似（分布の山の頂点付近をガウス分布で近似する方法）を使えばよさそうに見えるかもしれないが、実際には被積分関数（積分される中身の関数）が強く歪んだ形になりやすく、モード周辺だけをガウス分布で近似すると、確率質量（分布が持つ確率の量）の大部分をうまく捉えられないことがある。そのため、この節ではエビデンスを最大化する方法を採用する。

エビデンス最大化には大きく 2 つの方法がある。

• エビデンス関数を解析的に書き下し、微分して $\alpha$ と $\beta$ の再推定式を得る方法。
• EM アルゴリズム（隠れた変数があるモデルでよく使われる反復最適化法）を使う方法。

この節では前者を扱う。

3.5.1 エビデンス関数の評価

ここでは、候補となる重み $\mathbf{w}$ を 1 つに決め打ちせず、あり得る重み全体を積分してモデルを評価する。この題材で言えば、あるスペクトル方向に大きな係数を置いた 1 つのモデルだけを見るのではなく、その正則化設定のもとで自然に許される重み全体として、含水率データをどれだけ説明できるかを見る考え方である。

周辺尤度は重みを積分して得る

周辺尤度 $p(\mathbf{t}|\alpha,\beta)$ は、重み $\mathbf{w}$ を消去することで得られる。

$$p(\mathbf{t}|\alpha,\beta) = \int p(\mathbf{t}|\mathbf{w},\beta)p(\mathbf{w}|\alpha)\,\mathrm{d}\mathbf{w}.$$

(3.77)

これは、「重み $\mathbf{w}$ のどれか 1 つを選ぶ」のではなく、**「その $\alpha$ のもとであり得るすべての重みを考慮して、データ全体の説明力を評価する」**という意味である。

この式の中では、尤度（パラメータを固定したときに観測データがどれくらい起こりやすいか） $p(\mathbf{t}|\mathbf{w},\beta)$ と事前分布 $p(\mathbf{w}|\alpha)$ の両方がガウス分布（正規分布）である。そのため、指数部（$\exp\{\cdot\}$ の中身）を平方完成（2 次式を「中心からのずれの二乗」の形に直すこと）すれば積分できる。

この題材では、$p(\mathbf{t}|\mathbf{w},\beta)$ は「この回帰係数で各試料の含水率をどれくらい説明できるか」を表し、$p(\mathbf{w}|\alpha)$ は「波数方向や PLS 成分に対する係数が大きくなりすぎないか」を表す。つまり (3.77) は、当てはまりの良さと係数の控えめさを同時に見て、その正則化設定全体を評価している。

指数部の形

(3.11)、(3.12)、(3.52) から、エビデンス関数は次のように書ける。

$$p(\mathbf{t}|\alpha,\beta) = \left(\frac{\beta}{2\pi}\right)^{N/2} \left(\frac{\alpha}{2\pi}\right)^{M/2} \int \exp\left\{-E(\mathbf{w})\right\} d\mathbf{w}$$

(3.78)

ここで $E(\mathbf{w})$ は次の量である。

$$E(\mathbf{w}) = \beta E_D(\mathbf{w}) + \alpha E_W(\mathbf{w}) = \frac{\beta}{2} \|\mathbf{t} - \mathbf{\Phi}\mathbf{w}\|^2 + \frac{\alpha}{2} \mathbf{w}^{\mathrm{T}} \mathbf{w}.$$

(3.79)

この式の第 1 項は、モデルの予測 $\mathbf{\Phi}\mathbf{w}$ と実際の目標値 $\mathbf{t}$ のずれを測る。第 2 項は、重みが大きくなりすぎることを罰する。

したがって $E(\mathbf{w})$ は、正則化付き二乗和誤差と本質的に同じものである。$\beta$ はデータへの当てはまりをどれだけ重視するか、$\alpha$ は重みを小さく保つことをどれだけ重視するかを決める。

この題材に当てはめると、第 1 項は含水率の予測誤差である。高含水率側を大きく外すと二乗誤差が大きくなるため、この項は強く悪化する。第 2 項は、特定の波数や成分に極端な係数を置くことへの罰である。$\alpha$ が小さすぎると、モデルは水分帯だけでなく、樹種固有のベースライン、散乱、偶然のノイズにも重みを置きやすくなる。$\alpha$ は、そのような未知樹種で壊れやすい重みの使い方を抑える役割を持つ。

平方完成と事後分布の中心

積分を行うために、 $E(\mathbf{w})$ を $\mathbf{w}$ について平方完成する。

$$E(\mathbf{w}) = E(\mathbf{m}_N) + \frac{1}{2}(\mathbf{w} - \mathbf{m}_N)^{\mathrm{T}} \mathbf{A}(\mathbf{w} - \mathbf{m}_N)$$

(3.80)

これは、$E(\mathbf{w})$ が $\mathbf{m}_N$ のところで最小になり、その周りでは 2 次形式で増えることを表している。

ここで行列 $\mathbf{A}$ は次のように定義される。

$$\mathbf{A} = \alpha \mathbf{I} + \beta \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}$$

(3.81)

また、中心での誤差は次のようになる。

$$E(\mathbf{m}_N) = \frac{\beta}{2} \|\mathbf{t} - \mathbf{\Phi} \mathbf{m}_N\|^2 + \frac{\alpha}{2} \mathbf{m}_N^{\mathrm{T}} \mathbf{m}_N.$$

(3.82)

$\mathbf{A}$ は、誤差関数の 2 階微分の行列である。

$$\mathbf{A} = \nabla \nabla E(\mathbf{w})$$

(3.83)

このような 2 階微分の行列をヘッセ行列と呼ぶ。ヘッセ行列は、各方向にどれくらい急に誤差が増えるか、つまり曲率を表す。

この題材では、曲率が大きい方向は、係数を少し動かすだけで含水率の説明力が大きく変わるスペクトル方向である。たとえば水分吸収帯と強く結びついた安定した方向は、データによって強く制約されやすい。一方、曲率が小さい方向は、係数を動かしても train の誤差があまり変わらない方向であり、樹種差やノイズに近い不安定な方向である可能性がある。

$\mathbf{m}_N$ は次のように与えられる。

$$\mathbf{m}_N = \beta \mathbf{A}^{-1} \mathbf{\Phi}^{\mathrm{T}} \mathbf{t}.$$

(3.84)

これは、重みの事後分布の平均である。ガウス分布の場合、平均とモードは一致するので、最もありそうな重みでもある。

したがって $\mathbf{m}_N$ は、この題材では「正則化を考慮したうえで、各スペクトル特徴をどれだけ含水率予測に使うか」を表す係数ベクトルである。単なる最小二乗解ではなく、事前分布によって大きすぎる係数が抑えられた解になっている点が重要である。

ガウス積分を実行する

平方完成できたので、多変量ガウス分布（複数の変数を同時に扱う正規分布）の正規化係数（確率分布全体の面積を 1 にするための係数）を使って積分できる。

$$\int \exp\left\{-E(\mathbf{w})\right\} d\mathbf{w}$$

$$= \exp\left\{-E(\mathbf{m}_N)\right\} \int \exp\left\{-\frac{1}{2}(\mathbf{w} - \mathbf{m}_N)^{\mathrm{T}} \mathbf{A} (\mathbf{w} - \mathbf{m}_N)\right\} d\mathbf{w}$$

$$= \exp\left\{-E(\mathbf{m}_N)\right\} (2\pi)^{M/2} |\mathbf{A}|^{-1/2}.$$

(3.85)

この式の意味は次の通り。

• $\exp\{-E(\mathbf{m}_N)\}$ : 最もよい重みでの当てはまりの良さを表す。
• $|\mathbf{A}|^{-1/2}$ : その周りに、どれくらい広い重みの範囲が許されるかを表す。ここで $|\mathbf{A}|$ は行列式（行列が空間の体積を何倍にするかを表す量）である。

誤差の最小値だけでなく、良い重みがどれくらい広く存在するかもエビデンスに入っている点が重要である。複雑すぎるモデルは、訓練データにはよく合っても、許される重みの領域が狭くなり、エビデンス上は不利になることがある。

これはコンペの感覚とも合う。あるモデルが train や特定の CV 分割でだけよく当たり、係数を少し変えるとすぐ悪化するなら、そのモデルは不安定である可能性が高い。この題材では未知樹種に外挿する必要があるため、狭い条件でだけ成立する重みよりも、多少の不確実性を含んでも自然に説明できる重みの範囲が広いモデルの方が望ましい。

対数エビデンス

(3.78) と (3.85) から、対数周辺尤度は次の形になる。

$$\ln p(\mathbf{t}|\alpha,\beta) = \frac{M}{2} \ln \alpha + \frac{N}{2} \ln \beta - E(\mathbf{m}_N) - \frac{1}{2} \ln |\mathbf{A}| - \frac{N}{2} \ln(2\pi)$$

(3.86)

この式はエビデンス関数の中心である。各項はおおまかに次のように読める。

• $\frac{M}{2} \ln \alpha$ : 重みの事前分布の正規化に由来する項。
• $\frac{N}{2} \ln \beta$ : ノイズモデルの正規化に由来する項。
• $- E(\mathbf{m}_N)$ : 最良の重みでのデータへの当てはまりと重みの大きさへの罰。
• $-\frac{1}{2} \ln |\mathbf{A}|$ : モデルの複雑さに関係する項。
• $-\frac{N}{2} \ln(2\pi)$ : 定数項。

エビデンスは「データへの当てはまり」と「モデルの複雑さ」の両方を自動的に考慮する。このため、モデル選択にも使える。

この題材で言えば、対数エビデンスは「含水率をよく当てること」と「波数方向を使いすぎないこと」のバランスを取る基準である。これは、PLS 成分数を増やすか、水分帯特徴を追加するか、スケーリングを強くするかといった判断と同じ問題意識を持っている。ただし、実際のコンペでは未知樹種への分布シフトがあるため、エビデンス的な複雑度の見方を、LOGO-CV や Multi-CV gate と組み合わせて使うのが自然である。

この図は、エビデンスが単に訓練データへの当てはまりだけを見ているわけではないことを表している。左側の「当てはまり」が良いほどエビデンスは上がりやすいが、右側の「複雑さ」が大きすぎると、良い重みの範囲が狭くなり、エビデンス上は不利になる。したがって、エビデンスはよく当たるが、必要以上に複雑ではないモデルを選ぶための基準として働く。

多項式回帰でのモデル選択

多項式回帰（入力の 2 乗、3 乗などを使って曲線を当てはめる回帰モデル）では、多項式の次数 $M$ を変えることでモデルの複雑さが変わる。次数が低すぎるとデータに合わない。次数が高すぎると柔軟すぎて、不要な複雑さを持つ。

図 3.14 は、多項式の次数に対するモデルエビデンスを示している。

ここでの $M$ は多項式の次数を表しており、前節までの「重みベクトルの次元」としての $M$ とは文脈が少し違う。PRML では文脈に応じて同じ文字が使われるため、ここでは「モデルの複雑さを表す量」として読むのがよい。

この例では、事前分布として (1.65) の形を仮定し、 $\alpha=5\times 10^{-3}$ に固定している。したがって、ここで比較しているのは主に多項式の次数 $M$ の違いである。

この図では、エビデンスは $M=3$ のモデルを最も支持している。重要なのは、訓練誤差だけを見れば次数を上げるほどよくなりがちだが、エビデンスは複雑さへのペナルティも含むため、単に高次数のモデルを選ぶわけではないという点である。

$M=0$ ではデータへの当てはまりが悪いため、エビデンスは低くなる。 $M=1$ では当てはまりが大きく改善し、エビデンスが上がる。 $M=2$ では当てはまりがわずかしか改善しない一方、複雑さは増えるため、エビデンスは下がる。これは、元の人工データが正弦波に由来し、偶関数的な 2 次項を追加しても改善が小さいためである。 $M=3$ では当てはまりが大きく改善するため、最も高いエビデンスになる。それ以上の次数では、当てはまりの改善が小さいのに複雑さが増えるため、エビデンスは下がる。

この例から、エビデンスは「観測データをよく説明する、できるだけ単純なモデル」を好むことが分かる。

この題材での対応物は、多項式の次数ではなく、たとえば PLS 成分数、使う波数帯の数、前処理後に残すスペクトル方向の数である。PLS 成分を増やすと train 内の説明力は上がりやすいが、未知樹種の樹種差まで拾うと Private で崩れやすい。図 3.14 の教訓は、コンペでは「成分数や特徴量を増やせばよい」ではなく、改善に見合うだけの安定した説明力があるかを見るという形で使える。

3.5.2 エビデンス関数の最大化

次に、エビデンスを最大にする $\alpha$ と $\beta$ を求める。このコンペの実験ではハイパーパラメータを LOGO-CV などで比較するが、PRML ではまず、モデル内部の確率計算だけから正則化の強さとノイズの大きさを再推定する道筋を示す。

ここからは、対数エビデンス (3.86) を $\alpha$ と $\beta$ について最大化する。直接の目的は、エビデンスを最大にする $\widehat{\alpha}$ と $\widehat{\beta}$ を求めることである。

$\alpha$ の最大化

まず $\alpha$ について最大化する。そのために、次の固有ベクトル方程式を考える。

$$(\beta \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}) \mathbf{u}_i = \lambda_i \mathbf{u}_i.$$

(3.87)

この式は、行列 $\beta \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}$ の固有値（その方向が何倍に伸びるかを表す値）と固有ベクトル（行列を掛けても向きが変わらない方向）を定義している。固有ベクトル $\mathbf{u}_i$ はパラメータ空間の特別な方向を表し、固有値 $\lambda_i$ はその方向でデータがどれくらい重みを強く制約しているかを表す。

この題材では、$\mathbf{u}_i$ は「複数の波数や特徴量を組み合わせた 1 つのスペクトル方向」と考えられる。$\lambda_i$ が大きい方向は、含水率データがその方向の係数を強く決めている方向である。水分吸収帯と安定して結びつく方向はこの候補になる。一方、$\lambda_i$ が小さい方向は、train だけでは係数をあまり決められない方向であり、樹種差や偶然の変動に近い可能性がある。

(3.81) より、 $\mathbf{A}$ の固有値は $\alpha + \lambda_i$ になる。したがって、 $\ln |\mathbf{A}|$ の $\alpha$ 微分は次のように計算できる。

$$\frac{d}{d\alpha}\ln|\mathbf{A}| = \frac{d}{d\alpha}\ln\prod_{i}(\lambda_{i} + \alpha) = \frac{d}{d\alpha}\sum_{i}\ln(\lambda_{i} + \alpha) = \sum_{i}\frac{1}{\lambda_{i} + \alpha}.$$

(3.88)

行列式は固有値の積なので、対数を取ると固有値の対数の和になる。この変形により、微分が扱いやすくなっている。

対数エビデンス (3.86) を $\alpha$ で微分して 0 とおくと、停留条件（関数が最大または最小になる候補を表す条件）は次のようになる。

$$0 = \frac{M}{2\alpha} - \frac{1}{2}\mathbf{m}_N^{\mathrm{T}}\mathbf{m}_N - \frac{1}{2}\sum_{i}\frac{1}{\lambda_i + \alpha}.$$

(3.89)

これを整理すると、次の式が得られる。

$$\alpha \mathbf{m}_N^{\mathrm{T}} \mathbf{m}_N = M - \alpha \sum_i \frac{1}{\lambda_i + \alpha} = \gamma.$$

(3.90)

ここで新しく $\gamma$ が導入される。 $\gamma$ は次のようにも書ける。

$$\gamma = \sum_{i} \frac{\lambda_i}{\alpha + \lambda_i}.$$

(3.91)

この $\gamma$ は、後で有効パラメータ数と解釈される。

(3.90) を $\alpha$ について解くと、次の再推定式になる。

$$\alpha = \frac{\gamma}{\mathbf{m}_N^{\mathrm{T}} \mathbf{m}_N}.$$

(3.92)

ただし、これは閉じた形で一発で解ける式ではない。なぜなら、右辺の $\gamma$ も $\mathbf{m}_N$ も $\alpha$ に依存するからである。したがって、次のような反復計算を行う。

1. $\alpha$ の初期値を選ぶ。
2. その $\alpha$ を使って $\mathbf{m}_N$ を計算する。
3. 固有値から $\gamma$ を計算する。
4. (3.92) で $\alpha$ を更新する。
5. 収束するまで繰り返す。

このとき $\mathbf{\Phi}^{T}\mathbf{\Phi}$ はデータから決まっていて固定なので、その固有値は最初に一度計算しておけばよく、 $\lambda_i$ はそれに $\beta$ を掛けることで得られる。

ここで重要なのは、$\alpha$ を訓練データだけから決めていることである。最尤法（データを最も起こりやすくするパラメータを選ぶ方法）で複雑さを調整する場合のように、独立した検証データは必要ない。

コンペ上の感覚では、$\alpha$ の更新は「今の特徴量と前処理のもとで、どれくらい係数を抑えるべきか」を train から見積もる処理である。もし含水率を説明するために大きな係数が本当に必要なら、$\mathbf{m}_N^{\mathrm{T}}\mathbf{m}_N$ が大きくなり、$\alpha$ は小さくなる。逆に、小さい係数で十分なら、$\alpha$ は大きくなり、より強く正則化される。

$\beta$ の最大化

次に $\beta$ について対数エビデンスを最大化する。固有値 $\lambda_i$ は $\beta$ に比例するため、 $d\lambda_i/d\beta = \lambda_i/\beta$ である。したがって、

$$\frac{d}{d\beta}\ln|\mathbf{A}| = \frac{d}{d\beta}\sum_{i}\ln(\lambda_{i} + \alpha) = \frac{1}{\beta}\sum_{i}\frac{\lambda_{i}}{\lambda_{i} + \alpha} = \frac{\gamma}{\beta}.$$

(3.93)

となる。

対数エビデンスを $\beta$ で微分して 0 とおくと、停留条件は次のようになる。

$$0 = \frac{N}{2\beta} - \frac{1}{2} \sum_{n=1}^{N} \left\{ t_n - \mathbf{m}_N^{\mathrm{T}} \boldsymbol{\phi}(\mathbf{x}_n) \right\}^2 - \frac{\gamma}{2\beta}$$

(3.94)

これを整理すると、

$$\frac{1}{\beta} = \frac{1}{N - \gamma} \sum_{n=1}^{N} \left\{ t_n - \mathbf{m}_N^{\mathrm{T}} \boldsymbol{\phi}(\mathbf{x}_n) \right\}^2.$$

(3.95)

が得られる。

この式は、ノイズ分散 $1/\beta$ を、予測誤差の二乗和から推定している式である。ただし、分母が $N$ ではなく $N-\gamma$ になっている点が重要である。モデルがデータに合わせるために実質的に $\gamma$ 個の自由度（モデルがデータに合わせて自由に調整できる実質的な数）を使ったので、その分を差し引いていると理解できる。

$\beta$ も $\alpha$ と同様、右辺が $\beta$ に依存するため反復的に解く。 $\alpha$ と $\beta$ の両方をデータから決める場合は、 $\gamma$ を更新するたびに両方を更新する。

この題材では、残差には測定ノイズだけでなく、モデルが表現しきれない樹種差、乾燥挙動の違い、高含水率側の非線形性も含まれる。したがって $\beta$ は、単なる装置ノイズの精度ではなく、今の特徴量とモデルで説明しきれずに残るばらつきの小ささとして読むとよい。$\gamma$ を差し引くことは、モデルが train に合わせるために使った自由度を考慮し、残差を過小評価しないための補正である。

再推定の実用的な流れ

実装上は、概念的に次の流れになる。

1. $\alpha$ と $\beta$ の初期値を置く。
2. $\mathbf{A} = \alpha \mathbf{I} + \beta \mathbf{\Phi}^{\mathrm{T}} \mathbf{\Phi}$ を作る。
3. $\mathbf{m}_N$ を計算する。
4. 固有値 $\lambda_i$ から $\gamma$ を計算する。
5. (3.92) で $\alpha$ を更新する。
6. (3.95) で $\beta$ を更新する。
7. 値がほとんど変わらなくなるまで繰り返す。

この反復によって、交差検証を使わずに、エビデンスが大きくなるようなハイパーパラメータを得る。

ただし、ここで得られるのは、あくまで現在の特徴量とモデル内部で自然な $\alpha,\beta$ である。最終採用では、未知樹種への安定性を別途検証する必要がある。

図 3.15 は、尤度と事前分布がパラメータ空間でどのような関係にあるかを示している。

赤い等高線は尤度を表し、緑の等高線は事前分布を表す。$\alpha=0$ なら事前分布による抑制がないため、最尤解（尤度を最大にする重み） $\mathbf{w}_{\mathrm{ML}}$ が中心になる。$\alpha$ が 0 でないときは、事前分布が重みを 0 に近づけるため、MAP 解（事後分布を最大にする重み） $\mathbf{w}_{\mathrm{MAP}}=\mathbf{m}_N$ は最尤解よりも原点寄りになる。

ここで重要なのは、すべての方向で同じように重みが決まるわけではないという点である。データが強く制約している方向では、事前分布の影響は小さく、重みは最尤値に近くなる。データがほとんど情報を持っていない方向では、事前分布の影響が強く、重みは 0 に近くなる。

この題材では、水分に対応する安定した方向では、データが係数を比較的はっきり決める。一方で、特定樹種だけに現れるスペクトルの癖や、測定位置に由来する揺らぎの方向では、データの支持が弱く、事前分布によって係数が抑えられる。この図は、使うべき方向と抑えるべき方向を、データの制約の強さで分けるイメージとして読める。

3.5.3 有効パラメータ数

このようなスペクトル回帰では、入力列数が多いことと、実際に使われる情報方向が多いことは同じではない。$\gamma$ は、この違いを数式として表すための量である。

$\gamma$ は何を表すか

$\gamma$ は次の式で定義されていた。

$$\gamma = \sum_{i} \frac{\lambda_i}{\alpha + \lambda_i}.$$

(3.91)

各項 $\lambda_i/(\alpha+\lambda_i)$ は、0 と 1 の間の値を取る。なぜなら $\lambda_i$ は非負であり、 $\alpha$ は正だからである。

この比の意味は次の通り。

• $\lambda_i \gg \alpha$ のとき、 $\lambda_i/(\alpha+\lambda_i)$ は 1 に近くなる。この方向のパラメータはデータによってよく決まっている。
• $\lambda_i \ll \alpha$ のとき、 $\lambda_i/(\alpha+\lambda_i)$ は 0 に近くなる。この方向のパラメータはデータではあまり決まらず、事前分布によって 0 付近に抑えられる。

したがって、$\gamma$ は「データによって実質的に使われているパラメータ数」を表す。各項が 0 と 1 の間にあるので、全体として $0 \le \gamma \le M$ である。形式上のパラメータ数は $M$ だが、正則化が強かったり、データが少なかったりすると、実際に自由に動けるパラメータ数は $M$ より小さくなる。

このため、$\gamma$ は有効パラメータ数と呼ばれる。

この題材で重要なのは、$\gamma$ が「選んだ波数の本数」や「PLS 成分数」そのものではない点である。1555 個の波数があっても、近い波数同士は強く相関している。PLS 成分を 4 個使っても、その全てが同じ強さで含水率を説明しているとは限らない。$\gamma$ は、正則化を受けたあとに、実際にどれだけの方向がデータによって使われているかを読む量である。

$\alpha$ の再推定式の直感

$\alpha$ の再推定式は次の通りであった。

$$\alpha = \frac{\gamma}{\mathbf{m}_N^{\mathrm{T}} \mathbf{m}_N}.$$

(3.92)

分母 $\mathbf{m}_N^{\mathrm{T}} \mathbf{m}_N$ は、事後平均の重みベクトルの大きさの二乗である。重みが大きいなら、事前分布の精度 $\alpha$ は小さくなる。つまり、重みを強く 0 に押し込めるべきではない、という判断になる。

逆に、重みが小さいなら $\alpha$ は大きくなる。これは、データが大きな重みを必要としていないので、より強い正則化をかけてもよいという判断である。

ただし分子は単なる $M$ ではなく $\gamma$ である。実質的に使われているパラメータ数だけを考慮して、重みの大きさと正則化の強さを釣り合わせている。

この式は、この題材の実験で「特徴量を増やしたが、実際には使われていない方向が多い」という状況を考えると分かりやすい。形式的な特徴量数 $M$ が大きくても、$\gamma$ が小さければ、モデルはそれほど多くの自由度を使っていない。その場合、正則化の強さは、見かけの特徴量数ではなく、実際に使われた方向の数に合わせて調整される。

$\beta$ の再推定式と不偏分散

$\beta$ の再推定式は次の通りであった。

$$\frac{1}{\beta} = \frac{1}{N - \gamma} \sum_{n=1}^{N} \left\{ t_n - \mathbf{m}_N^{\mathrm{T}} \boldsymbol{\phi}(\mathbf{x}_n) \right\}^2.$$

(3.95)

これは、ノイズ分散を残差平方和（予測誤差を二乗して足し合わせた量）から推定する式である。分母が $N-\gamma$ であることは、単純なガウス分布の分散推定における $N-1$ と似ている。

単一変数のガウス分布で、平均もデータから推定する場合、分散の最尤推定（尤度を最大にする推定）は次の形である。

$$\sigma_{\rm ML}^2 = \frac{1}{N} \sum_{n=1}^{N} (x_n - \mu_{\rm ML})^2$$

(3.96)

しかしこれは、平均 $\mu_{\rm ML}$ がデータに合わせて推定されているため、分散を小さく見積もりがちである。平均の推定に 1 つの自由度を使ったと考えると、不偏推定（平均的に真の値からずれないよう補正した推定）では分母を $N-1$ にする。

$$\sigma_{\text{MAP}}^2 = \frac{1}{N-1} \sum_{n=1}^{N} (x_n - \mu_{\text{ML}})^2.$$

(3.97)

線形回帰では、平均関数 $\mathbf{w}^{\mathrm{T}} \phi(\mathbf{x})$ を推定するために複数のパラメータを使う。ただし、その全てが完全にデータに合わせて自由に動くわけではなく、実際に使われた自由度は $\gamma$ である。したがって分母は $N-\gamma$ になる。

この見方をすると、エビデンス近似は、ノイズ分散の推定においてもモデルの実効的な複雑さを自然に補正していることが分かる。

この題材では、高含水率側や特定樹種で大きな残差が出ることがある。このとき、単に残差平方和を $N$ で割ると、モデルが train に合わせるために使った自由度を無視してしまう。$N-\gamma$ で割ることは、モデルがどれだけデータに合わせ込んだかを差し引いて、残りのばらつきを見るという意味を持つ。

正弦波データでの $\alpha$ の決定

図 3.16 は、正弦波の人工データ集合に対して、エビデンス手続きがどのように $\alpha$ を決めるかを示している。この例では、ガウス基底関数が 9 個あり、バイアス（入力によらず足される切片項）を含めるとパラメータ数は $M=10$ である。説明を簡単にするため、 $\beta$ は真の値 11.1 に固定されている。

左図では、$\gamma$ と $2\alpha E_W(\mathbf{m}_N)$ の交点が、(3.92) を満たす $\alpha$ に対応する。右図では、その同じ点で対数エビデンスが最大になっていることが示されている。さらに、テスト誤差の最小付近とも近いことが示されている。

これは、エビデンス最大化が訓練データだけを使っているにもかかわらず、汎化性能（未知データに対する予測性能）のよいハイパーパラメータを選びやすいことを示している。

この題材でこの図を読むなら、横軸の $\alpha$ は「スペクトル方向への重みをどれくらい抑えるか」であり、テスト誤差に近い考え方は未知樹種への性能に対応する。ただし、図 3.16 の人工データと違い、このコンペでは樹種分布が train/test で変わるため、この図は正則化選択の直感として読む必要がある。

$\alpha$ と重みの大きさ

図 3.17 は、 $\alpha$ を変化させたときに、各重みがどのように変わるかを、有効パラメータ数 $\gamma$ に対して示している。

$\alpha$ が大きいと、事前分布による抑制が強いため、多くの重みは 0 近くになる。このとき $\gamma$ は小さく、実質的に使われるパラメータ数は少ない。

$\alpha$ が小さいと、重みへの抑制が弱くなり、より多くのパラメータがデータに合わせて動けるようになる。このとき $\gamma$ は $M$ に近づく。

この題材では、$\alpha$ が大きいモデルは、水分帯のような強い共通パターンだけを使う保守的なモデルに近い。$\alpha$ が小さいモデルは、より多くの波数方向を使えるため train には合いやすいが、樹種固有の癖や散乱まで拾う危険がある。したがって $\alpha$ は、未知樹種に持ち出せる情報だけを使うためのブレーキとしてイメージできる。

データ数が非常に多い場合

データ数 $N$ がパラメータ数 $M$ に比べて非常に大きい場合、すべてのパラメータはデータによってよく決まるようになる。これは $\mathbf{\Phi}^{T}\mathbf{\Phi}$ がデータ点についての和を含むため、データが増えるほど固有値 $\lambda_i$ が大きくなるからである。

この極限では $\gamma = M$ となり、再推定式は次のように簡単になる。

$$\alpha = \frac{M}{2E_W(\mathbf{m}_N)}$$

$$\beta = \frac{N}{2E_D(\mathbf{m}_N)}$$

(3.98)

ここで $E_W$ は重みの大きさに関する項、 $E_D$ はデータ誤差に関する項である。この近似は、ヘッセ行列の固有値スペクトル（固有値を小さい順または大きい順に並べたもの）を評価しなくてよいので、計算が簡単である。

このコンペの train は 1322 件であり、スペクトル列数 1555 と同程度の規模である。さらに樹種単位の分布差が大きいので、単純に「データ数が十分大きいから $\gamma=M$ とみなす」とは考えにくい。むしろ、正則化や PLS のような次元削減を使い、実質的に使う方向を抑えることが重要である。

3.6 固定基底関数の限界

ここまでの議論は、前処理や特徴量 $\boldsymbol{\phi}(\mathbf{x})$ が固定されていることを前提にしていた。この題材で言えば、水分帯平均や 7618 cm^-1 付近の特徴のように手で定めた特徴、あるいは PLS のように train 内でいったん低次元方向を決め、その後はその score 上で線形重みを学習する表現が対応する。この節では、そのような固定的な表現の限界を確認する。

ここまでの線形基底関数モデルの利点

本章で扱ってきたモデルは、固定された非線形基底関数（入力に対して直線ではない変換をする基底関数）の線形結合である。入力 $\mathbf{x}$ に対して非線形な基底関数を使うので、入力と目標値の関係としては非線形なものを表せる。一方で、重み $\mathbf{w}$ には線形なので、最小二乗解やベイズ的な事後分布を解析的に扱いやすいという利点がある。

つまり、固定基底関数モデルは次のような良い性質を持つ。

• 最小二乗問題に閉形式解（反復計算せず、式で直接書ける解）がある。
• ガウス事前分布とガウス尤度を使うと、ベイズ推論が扱いやすい。
• 基底関数をうまく選べば、入力から目標値への複雑な非線形写像（直線では表せない変換）も表せる。

このため、一見するとパターン認識（データから規則性を見つけて分類や予測を行う分野）の一般的な枠組みに見える。

この題材で PLS や固定水分帯特徴が使いやすいのは、この利点があるからである。厳密には PLS の方向は train データから推定されるが、いったん score 表現を作れば、その後は少数の固定された特徴上の線形回帰として扱える。前処理と特徴量を固定すれば、回帰係数の意味や水分吸収帯の利用有無を追いやすく、コンペ資料としても予測の根拠を説明しやすい。

問題は基底関数を事前に固定すること

しかし固定基底関数モデルには重大な限界がある。最大の問題は、基底関数 $\phi_j(\mathbf{x})$ を訓練データを見る前に固定していることである。

入力空間（入力が取り得る値全体の空間）の次元を $D$ とすると、空間全体を十分に覆うために必要な基底関数の数は、$D$ が増えるにつれて急速に増える。多くの場合、指数的に増える。これが次元の呪い（次元が増えると必要なデータや基底関数が急増する問題）である。

低次元では、基底関数を格子状に置いて空間を覆うことができるかもしれない。しかし高次元では、空間の体積が非常に大きくなり、ほとんどの領域にはデータが存在しない。それにもかかわらず全空間を覆うように基底関数を置くと、膨大な数の基底関数が必要になる。

この題材の 1555 次元スペクトル空間でも同じ問題が起きる。すべての波数の組み合わせを細かく覆うような特徴を手で用意することは現実的ではない。さらに test は未知樹種なので、train に存在しない領域へ外挿する必要がある。固定した特徴量だけで空間全体を覆う発想は、この問題に対して弱い。

実データが持つ助けになる性質

幸い、実データにはこの問題を和らげる性質がある。

第一に、データは入力空間全体に一様に散らばっているわけではなく、より低次元の非線形多様体（高次元空間の中にある、少ない自由度で表せる曲がった面のような構造）の近くに存在することが多い。たとえば手書き数字画像は、画素数だけ見れば高次元だが、実際に自然な数字画像が占める領域はその中のごく限られた部分である。

この題材でも、1555 次元すべてが独立に動くわけではない。近い波数同士は連続したスペクトルとして一緒に変化し、水分量、散乱、ベースライン、樹種差といった少数の要因が大きな変動を作る。EDA の PCA で少数の主成分が大きな分散を説明していたことは、この見方と対応する。

この場合、局所的な基底関数（入力空間の一部だけに強く反応する基底関数）を使うなら、データがある領域にだけ基底関数を配置すればよく、空間全体を覆う必要はない。この考え方は、動径基底関数ネットワーク（局所的な山型の基底関数を使うモデル）、サポートベクトルマシン（境界や予測に効く少数のデータ点を重視するモデル）、関連ベクトルマシン（少数の代表的なデータ点を残すベイズ的なスパースモデル）などで使われる。

第二に、目標変数（予測したい値）は、データ多様体のすべての方向に依存しているわけではなく、少数の重要な方向にだけ強く依存していることがある。ニューラルネットワーク（層状に並んだ多数の変換をデータから学習するモデル）は、基底関数を固定するのではなく、パラメータを学習することで、どの方向に反応すべきかをデータから適応的に選べる。

この題材では、含水率は水分吸収帯に強く関係するが、スペクトル上のすべての変動が含水率に必要なわけではない。固定特徴では「どの方向が重要か」を人間が事前に決める必要がある。適応的なモデルは、その重要方向をデータから学べる可能性がある一方、データ数や樹種分布の制約があるため、過学習への注意が必要である。

固定基底から適応基底へ

固定基底関数モデルでは、基底関数の位置や形はあらかじめ決まっている。そのため、データの分布や目標値との関係に合わせて基底関数自体を変えることができない。

一方、ニューラルネットワークのようなモデルでは、基底関数に相当する中間表現が学習によって変わる。これにより、データが存在する多様体や、目標値に関係する方向に合わせて、モデルが表現を調整できる。

このため、線形基底関数モデルは理論的に重要で扱いやすい出発点だが、高次元の実問題では、固定基底を入力空間全体にあらかじめ敷き詰める発想だけでは足りないことが多い。データがある領域を中心に基底を置く RBF、SVM、RVM 的な方法や、基底に相当する表現を学習するニューラルネットワークが、後続の候補になる。

この題材で見ると、固定水分帯特徴や PLS は、少ないデータでも安定しやすいベースラインである。一方、特徴量や前処理を少し変えるだけで、Public LB と樹種単位 CV の評価がずれることがある。これは、固定特徴だけでは未知樹種への頑健性を完全には制御しきれないことを示している。PRML 3.6 は、その限界を理論的に説明している。

この図は、3.5 と 3.6 の全体像をまとめている。前半では、エビデンス近似によって $\alpha$ と $\beta$ をデータから決め、$\gamma$ によって実質的に使われているパラメータ数を読む。後半では、固定基底関数モデルが高次元では限界を持つため、データのある領域や重要な方向に合わせて表現を変えられるモデルが必要になる、という流れである。

実例: 採用モデルは何をしているか

ここまでのエビデンス近似、有効パラメータ数、固定基底関数モデルの議論を、採用モデルに当てはめる。

2026-05-28 時点の単一 champion は、次のブレンドモデルである。

step094 = 0.500 * step081 + 0.500 * region-wise PLS

ここで step081 は、エビデンス近似を使う ARD 型の基礎モデルである。構成は次の通りである。

SG window 47 -> stable bands -> Pareto scaling -> PLS 4-score basis -> ARDRegression -> clip0

一方、region-wise PLS は、水分に関係しやすい複数の波数帯を領域ごとに分け、それぞれで PLS score を作ってから回帰するモデルである。step094 は、エビデンス近似を使う step081 と、領域別に水分帯情報を拾う region-wise PLS を半分ずつ混ぜたモデルである。

step081 の Public LB は 14.566924156579136 であり、step094 の Public LB は 13.359541375095974 である。つまり、region-wise PLS 方向を加えることで Public LB は大きく改善した。

ただし、この Public LB は最終汎化性能そのものではない。最終順位は Private（コンペ終了後に使われる非公開評価部分）で決まり、Public は test の一部だけで計算される暫定スコアである。この題材では、Public と Private の上位層に大きな乖離が見られる。

この表は、Public で極端に良いモデルが、Private でも同じように良いとは限らないことを示している。Public 上位10名の Public 平均は 2.6721 と非常に低いが、同じ参加者群の Private 平均は 18.8475 まで悪化している。一方、Private 上位10名の Public 平均は 16.4490 であり、Public だけを見れば必ずしも目立たない。したがって、Public LB は補助的な外部観測として記録し、採用判断の中心は train 内で完結する LOGO-CV、LKSO、Species Cluster CV に置く。

PRML の言葉で言えば、訓練誤差だけで複雑なモデルを選ぶと過学習するのと同じように、Public LB だけで前処理、特徴量、モデル構造を選び続けると、Public 部分にだけ適合したモデルを選んでしまう。エビデンス近似は「train データが自然に支持する複雑さ」を見る道具であり、Public LB はその外側にある、偏りうる観測値として扱う必要がある。

step081: エビデンス近似を使う基礎モデル

step081 のポイントは、最後の回帰器が単なる PLS 回帰ではなく、PLS で作った 4 次元の基底に ARD をかけていることである。実装上は pls_pareto_ard で、ARDRegression(compute_score=True) を使っている。

compute_score=True は、反復中の対数周辺尤度、つまりエビデンスの値を記録するための設定である。ARDRegression 自体は、重みの精度とノイズ精度をエビデンス最大化の反復で推定する。scikit-learn の属性名ではノイズ精度が alpha_、重み精度が lambda_ と呼ばれるため、ここでは混乱を避けるために、PRML の意味に合わせて「重みを抑える精度」を $\alpha_j$ と書く。

PRML の記号に対応させると、次のように読める。

SG 平滑化の window は、1 点の値をならすときに周囲の何点を使うかを表す。window が小さいと細かな形を残しやすく、window が大きいとノイズを強くならす一方で、細かなピークも弱める。ここでは、スペクトルの局所的なノイズや測定ゆらぎを抑え、PLS が安定した方向を拾いやすくするために使っている。

PLS は PCA と似て低次元表現を作るが、目的が違う。PCA は入力 $\mathbf{x}$ の分散が大きい方向を探す。一方、PLS は目標値 $t$、つまり含水率と関係しやすい方向を探す。したがって、単にスペクトル上で大きく変動する方向ではなく、含水率予測に使いやすい方向を基底として作れる。ここでいう 4-score basis は、「元の 1555 次元を、含水率に関係しやすい 4 個の要約値へ変換したもの」と考えるとよい。

PLS 方向は train fold の中で推定されるので、PRML 3.6 が想定する「訓練データを見る前から完全に固定された基底」とは厳密には同じではない。ただし、ARDRegression から見ると、入力はすでに 4 個の PLS score に変換されており、その score を $\boldsymbol{\phi}(\mathbf{x})$ とする線形基底関数モデルとして扱える。

つまり、step081 では、元の 1555 次元のスペクトルに直接重みを置くのではなく、まず PLS によって「含水率と関係しやすい 4 つの方向」を作っている。この 4 つの PLS score が、PRML の $\boldsymbol{\phi}(\mathbf{x})$ に相当する。最終的な予測は、おおまかには

$$y(\mathbf{x}, \mathbf{w}) = \mathbf{w}^{\mathrm{T}}\boldsymbol{\phi}(\mathbf{x})$$

という線形基底関数モデルとして読める。

ただし、4 つの PLS score がすべて同じくらい信用できるとは限らない。ある成分は水分に対応する安定した方向かもしれないが、別の成分は特定樹種だけの散乱、ベースライン、測定位置の癖を含んでいるかもしれない。そこで ARD では、各係数 $w_j$ に対して別々の重み精度を持たせる。

PRML 本文の単純な形では、重み全体に 1 つの $\alpha$ を置いていた。

$$p(\mathbf{w}|\alpha) = \mathcal{N}(\mathbf{w}|0,\alpha^{-1}\mathbf{I})$$

ARD ではこれを、成分ごとの精度を持つ形に拡張していると考えられる。

$$p(\mathbf{w}|\boldsymbol{\alpha}) = \prod_j \mathcal{N}(w_j|0,\alpha_j^{-1})$$

$\alpha_j$ が大きい成分では、対応する $w_j$ は 0 に強く近づけられる。逆に、データがその成分を強く支持していれば、$\alpha_j$ は相対的に小さくなり、係数は残りやすい。これは、PLS score の中で、どの方向を使い、どの方向を抑えるかをエビデンスから決める操作である。

この意味で、step081 はエビデンス近似の実用例になっている。人間が「4 つの PLS 成分のうち、何番目の係数をどれくらい縮めるか」を手で決めているのではない。訓練データに対する周辺尤度、つまりエビデンスが大きくなるように、ARD が係数ごとの正則化の強さを推定している。

重要なのは、ここでエビデンス近似が解いている問題と、コンペの最終目的を混同しないことである。ARD のエビデンス最大化は、現在の特徴量と train データの中で自然な重みの抑え方を決める。一方、この題材の本質的な難しさは未知樹種への外挿である。したがって、エビデンス近似で作った候補でも、LOGO-CV、LKSO、Multi-CV gate で別途検証する必要がある。

step081 が基礎モデルとして採用された理由は、エビデンス近似を使っているだけではなく、未知樹種を想定した検証でも崩れにくかったからである。

step081 の構成では、PLS 成分数、安定波数帯、Pareto scaling、ARD、clip0 を固定し、前処理の平滑化幅も含めて train 内の検証で破綻しないことを確認している。したがって、step081 は「エビデンス近似を含む ARD 型の安定した基底回帰を、未知樹種想定の検証と Public LB の両方で確認した基礎モデル」と解釈できる。

step094: 安定した基礎モデルに別の水分帯情報を足す

step094 は、step081 を捨てたモデルではない。step081 の予測を半分残し、残り半分を region-wise PLS に置き換えたモデルである。

region-wise PLS は、1555 次元のスペクトル全体を一度に扱うのではなく、主に水分と関係する複数の波数帯を領域ごとに分ける。そして各領域で PLS score を作り、それらを少数の特徴量として回帰する。この処理は、PRML の $\boldsymbol{\phi}(\mathbf{x})$ を「水分帯ごとの PLS score」として作り直している、と読める。

この追加は、3.6 の固定基底関数の限界とも関係する。step081 では、安定波数帯をまとめて 4 個の PLS score に圧縮していた。一方、region-wise PLS では、水分帯の中でも領域ごとに異なる情報を別々に取り出す。つまり、同じ水分情報でも、1 つの固定された要約だけに任せず、領域ごとの基底表現を追加する操作である。

step094 の主な指標は次の通りである。

ただし、ここで重要なのは「Public LB が良かったから無条件で採用する」ことではない。step094 周辺では、w=0.475 と w=0.500 の LOGO primary 差は 0.000648 と非常に小さい。つまり、w=0.500 は鋭い一点だけが偶然よかったというより、w=0.45 から w=0.525 付近の平坦な領域にある。

一方、追加の train-only stress CV（train だけを使い、未知樹種セットをより厳しく作る検証）では、w=0.450 から w=0.475 側が平均的にはやや安定していた。したがって、w=0.500 は Public LB では最良だが、Private 安定性だけを見ると少し攻めた位置である。

この判断は、PRML 3.5 のエビデンスの考え方とよく対応する。エビデンスは、訓練データへの当てはまりだけでなく、モデルがどれくらい狭い条件に依存しているかも見る。コンペでも同じで、特定の Public LB だけに合わせて w=0.525、w=0.550 と右側へ追い続けると、Public 部分にだけ過適合する危険がある。したがって、step094 は採用するが、さらに右側の重みを追わないという判断になる。

normal_pass だけでは十分ではない

この実験では、step095 として 0.700 * step094 + 0.300 * LGBM segment-stat という別候補も試した。これは normal_pass（通常提出候補として用意した strict scorecard を通ること）を満たしていた。つまり、LOGO 上の平均性能、複数 fold での改善、bootstrap CI などは一見よく見えていた。

しかし、step095 の Public LB は 15.03225537584012 まで悪化した。これは、normal_pass していても、test 予測の移動量が大きすぎる候補は危険であることを示す。PRML 的に言えば、訓練データ内の基準でよく見える複雑な補正でも、未知樹種に対して同じ構造が成り立つとは限らない。

したがって、この題材での採用判断は次のように整理できる。

• step081 は、ARD によるエビデンス近似を使った安定した基礎モデルである。
• step094 は、step081 に region-wise PLS の水分帯情報を加えた現在の単一 champion である。
• step094 は Public LB と LOGO 近傍の両方で妥当だが、右側へさらに重みを増やす根拠は弱い。
• normal_pass しても、test 予測が大きく動く候補は Public/Private の分布差で崩れる可能性がある。

逆に、エビデンス近似を使っていれば常に勝てるわけではない。似た構成でも、特定の LKSO subset で崩れる候補や、Public LB で悪化する候補はある。エビデンスはモデル内部の重みをどう抑えるかには強いが、未知樹種の分布ずれを完全に知っているわけではない。この題材では、エビデンス近似は step081 の中核であり、step094 でもその予測を半分残している。ただし、それを採用可能にしているのは、樹種単位の CV gate と組み合わせている点である。

まとめ

この節の要点は次の通り。

• エビデンス近似では、ハイパーパラメータ $\alpha,\beta$ を周辺尤度 $p(\mathbf{t}|\alpha,\beta)$ が最大になるように決める。
• 周辺尤度は、重み $\mathbf{w}$ を積分して消去することで得られる。
• 対数エビデンスには、データへの当てはまりとモデル複雑度の両方が入っている。
• $\alpha$ と $\beta$ は反復的な再推定式 (3.92) と (3.95) によって求められる。
• $\gamma$ は、形式的なパラメータ数ではなく、データによって実質的に決定された有効パラメータ数を表す。
• $\beta$ の推定式で分母が $N-\gamma$ になるのは、モデルがデータに合わせるために使った自由度を補正しているから。
• 固定基底関数モデルは解析的に扱いやすいが、高次元では基底関数の数が急増するという限界がある。
• 実データは低次元多様体上にあることが多く、局所基底、カーネル法、ニューラルネットワークなどはこの性質を利用できる。
• このような高次元スペクトル回帰では、エビデンス近似は正則化と有効パラメータ数を理解する実践的な視点を与える。
• step081 は、SG 平滑化と PLS で 4 次元の基底を作り、その PLS score に ARDRegression をかけることで、エビデンス近似を実際の正則化として使っている。
• 現在の単一 champion である step094 は、step081 に region-wise PLS を半分ブレンドし、水分帯ごとの安定した変動を追加で利用している。
• w=0.500 は Public LB では強いが、train-only stress CV では w=0.450 から w=0.475 側がやや安定していたため、w=0.525 以上を追い続けるべきではない。
• ただし未知樹種への汎化性能は、LOGO-CV、LKSO、Species Cluster CV のような樹種単位の検証で別途確認する必要がある。
• Public LB と Private は大きく乖離しうるため、Public の低スコアだけを最終汎化性能の証拠として扱わない。Public は補助情報として記録し、モデル選択は train 内で完結する複数の検証と組み合わせて行う。