はじめに

X線天文学のスペクトル解析では、各エネルギービンの信号を「カウント」として扱います。そのため、揺らぎの記述には Poisson 分布が自然に現れ、Cash 統計（C-stat; Cash 1979） が広く用いられています。

一方で、実際の解析では

なぜ Cash 統計が最尤推定に対応しているのか
なぜ $\chi^2$ と似た扱いができると言われるのか
XSPEC で使われている C-stat は、Cash (1979) の定義と何が同じで何が違うのか
$\Delta C$ は、どのような量として解釈すればよいのか

といった点が、明示されないまま使われていることも少なくありません。

本記事では、「Poisson 分布に基づく最尤推定を行う」という問題設定を出発点に、
Cash 統計がどのような発想で導入され、
その後 「フィットの良さ」や「モデル間の差」を評価しやすい形へと整理されてきたのかを、数式を追いながら整理します。

特に、

Cash (1979) における元々の立場
観測を完全に再現する仮想モデルを基準にした尤度比という考え方
現在 XSPEC で用いられている C-stat の位置づけ

の関係を、ひとつの流れとして理解することを目標とします。

1. Poisson 尤度

各ビン $i$ の観測 $k_i$ が Poisson($\mu_i$) に従うとすると、尤度は

\mathcal{L}
= \prod_{i=1}^N
\frac{e^{-\mu_i} \mu_i^{k_i}}{k_i!}

対数を取ると

\ln \mathcal{L}
= \sum_i
\left(
- \mu_i + k_i \ln \mu_i - \ln k_i!
\right)\tag{1}

ここで $k_i$ は観測なので固定です。また $\ln k_i!$ は $\mu_i$ に依存しないため、最尤推定の位置そのものには影響しません。

2. Cash (1979) における統計量の定義

式(1)で与えられる対数尤度のうち、Cash (1979) はモデルパラメータ $\mu_i$ に依存する部分のみを取り出し、それに $-2$ を掛けた統計量として

\tilde{C}\equiv 2\sum_i(\mu_i - k_i\ln\mu_i)\tag{2}

という形を用いました。
これは $\tilde C$ を最小化することは $\ln\mathcal{L}$ を最大化することと同値です。

3. 観測を完全に再現する仮想モデル

次に、「観測を完全に再現する仮想モデル」を考えます。これは各ビンの平均 $\mu_i$ を 独立に自由に動かせるモデルです。
このとき、式(1)を各 $\mu_i$ で最大化すると

\frac{\partial \ln \mathcal{L}}{\partial \mu_i}
= -1 + \frac{k_i}{\mu_i}=0
\quad\Rightarrow\quad
\mu_i = k_i\tag{3}

が得られます。

つまり、各ビンだけを見れば「平均を観測値に一致させた Poisson」が最も尤もらしい、ということになります。これが $\mu_i=k_i$ が特別になる理由です。

4. Castor statistic（XSPECで現在使われている形式）

ここで、モデルと完全再現モデルの尤度比を

C
\equiv -2\ln
\frac{\mathcal{L}(\text{model})}
{\mathcal{L}(\text{perfect})}\tag{4}

と定義します。
perfect 側には $\mu_i=k_i$ を代入します。

式(1)を用いて差を取ると、$\ln k_i!$ は打ち消し合い、

C
= 2\sum_i\left(
\mu_i - k_i + k_i\ln\frac{k_i}{\mu_i}
\right)\tag{5}

が得られます。
この $C$ は、Wilks の定理により漸近的には $\chi^2$ 分布に従うことが知られていますが、$\mu_i$ が極端に小さい場合などでは注意が必要です（e.g., Kaastra 2017, Bonamente 2020）。

XSPEC への組み込み

XSPEC では 1998 年 10 月以降、この形の $C$（いわゆる Castor statistic）が
C-stat として実装されています。この統計量は、
John Castor による private communication に基づくものとして知られています。

出典：https://wwwastro.msfc.nasa.gov/qdp/help/fit.html

5. なぜベストフィットは変わらないのか

$\tilde C$ と $C$ の違いは、観測データ $k_i$ のみで決まる項のみです。つまり、

C = \tilde C + A(\{k_i\})

と書けます。ここで $A(\{k_i\})$ は観測が固定なら定数なので、最小化のパラメータには影響しません。
そのため、$\tilde C$ を最小化しても $C$ を最小化しても、得られるベストフィットは同じです。

6. Cash 統計の思想と現在の立場（私見）

Cash (1979) が当初から $-2$ を掛けた形を用いたのは、$\tilde C$ の絶対値そのものによってモデルの良し悪しを判断することを目的としていたからではなく、最尤推定の結果をもとにパラメータの不確かさを評価することを重視していたためだと考えられます。実務的には、あるモデルが観測をどの程度説明できているかは、残差の分布やスペクトルの見え方を確認することで、ある程度は把握できます。それ以上に重要になるのは、ベストフィットの周辺でパラメータをどこまで動かせるのか、すなわちパラメータ誤差や相関をどのように評価するか、という点です。

この文脈で、モデル A と B に対して

\Delta \tilde C = \tilde C(\text{A}) - \tilde C(\text{B})

を考えると、観測データ $k_i$ のみに依存する項は完全に打ち消し合い、$\Delta \tilde C$ は純粋な尤度比（$-2\log$ 尤度比）として解釈できます。この構造は、$\chi^2$ フィットにおける $\Delta\chi^2$ と本質的に同じであり、十分なデータ数がある極限では、Wilks の定理により、$\Delta \tilde C$ は 両モデルのパラメータ数の差を自由度とする $\chi^2$ 分布に従うと期待されます。

一方で、実際の解析ではパラメータの許容範囲だけでなく、モデル全体として観測とどの程度整合しているかを把握したい場面もあります。そのため、完全再現モデルを基準にした尤度比として整理された現在の $C$ が、実用的な統計量として XSPEC で用いられているのだと考えられます。

まとめ

本記事では、Cash 統計が Poisson 分布に基づく最尤推定のために導入された統計量であることを整理しました。Cash (1979) の本来の立場では、モデルパラメータの推定そのものが主眼であり、$\Delta\tilde{C}$ を用いて $\Delta\chi^2$ と対応づけるという考え方に重点が置かれていました。

その後、フィットの妥当性評価やモデル間比較の必要性から、完全再現モデルを基準にした尤度比という形に整理された $C$ が XSPEC で用いられるようになったと考えられます。

本記事の整理が、XSPEC における C-stat フィットの背景や設計思想を理解する一助となれば幸いです。

参考