二値アウトカムのCVでAUCだけを報告すると、リスクの高い人を高い順に並べられたかどうかは分かりますが、予測した確率そのものが正しいかは分かりません。30%と表示された患者の実際の発生割合が15%でも、順序が保たれていればAUCは変わらないためです。
生存時間解析では、さらに打ち切りが加わります。誰がイベントを起こしたかが追跡期間の途中までしか分からないので、評価指標も打ち切りを扱う形に替える必要があります。
この記事では、AUCと較正を分けて読む理由、生存時間解析でのCV、そしてfold間の標準偏差から信頼区間を作れるかという問題を扱います。分割の単位については、施設・時系列・イベント数が少ない場合で扱いました。
AUCと較正を別に読む理由
AUCは、イベント例の予測スコアが、非イベント例より高くなる程度を評価します。予測確率の絶対的な水準よりも、順序に注目する指標です。
ここで、イベント例の集合を$\mathcal P$、非イベント例の集合を$\mathcal N$、その人数を$n_1,n_0$、予測スコアを$s_i$とすると、同点を半分として扱う経験的AUCは次の式です。
$$\widehat{\mathrm{AUC}} =\frac{1}{n_1n_0} \sum_{i\in\mathcal P}\sum_{j\in\mathcal N} \left\{\mathbf1(s_i>s_j)+\frac12\mathbf1(s_i=s_j)\right\} \tag{1}$$
式(1)は、イベント例と非イベント例の全ての組を比べ、イベント例のスコアが高い組の割合を求めています。記号$\mathbf1(\cdot)$は、括弧内の条件が成り立つ場合に1、それ以外に0を取ります。
fold平均と全OOFのAUC
全OOF予測を一つにまとめたAUCには、異なるfoldの対象者同士の比較が含まれます。一方、各foldのAUCを平均した場合には、その比較は含まれません。
そのため、両者が一致するとは限らず、全OOFのAUCが常に最も正しい集約法という説明は避けます。特に、foldごとにモデルのスコアの位置や尺度が大きく変わる場合には、別モデルのスコアを直接比べる意味を確認します。予測値を集めた後に全体の正解を使ってスコアを調整し、そのままAUCを再評価する方法も、独立した評価にはなりません。
予測モデルの検証で扱うように、識別性能が良くても、予測確率と実際の発生割合が対応しているとは限りません。Van Calsterらは、この較正の評価が予測研究で見落とされやすいことを論じています。
較正の評価と更新の区別
そこで、二値アウトカムでは、予測確率のlogitを説明変数とするモデルから、較正の切片$a$と傾き$b$を調べる方法があります。ここでlogitは、確率$p$を$\log\{p/(1-p)\}$へ変換する関数です。
$$\operatorname{logit}\Pr(Y=1\mid\widehat p) =a+b\,\operatorname{logit}(\widehat p) \tag{2}$$
式(2)は、予測確率と実際の発生割合の対応を、切片と傾きで要約しています。この形で理想的な対応が成り立つ場合には$a=0$、$b=1$ですが、その値に近いことだけで、あらゆる範囲の較正が良いと確認できるわけではありません。
また、傾きと切片を同時に推定する較正モデルの切片と、傾きを1に固定して求めるcalibration-in-the-largeは区別します。名称が似ていても、同じ計算ではありません。
なお、評価集合のアウトカムを使って較正のズレを測ること自体は、性能評価の一部です。しかし、そこから推定した切片や傾きで予測を修正し、同じ集合で修正後の性能を最終評価する場合には、更新に使ったデータで採点しています。閾値選択や確率の再較正まで含めた予測手順を評価するなら、その更新過程も学習側へ入れます。
同じAUCで較正だけが崩れる例
2,000人について、一つのマーカーからイベント確率を生成し、正しい確率を返すモデルと、そこへ一定のズレを加えたモデルを比べます。後者は全員の予測確率を高い側へ動かしているだけなので、患者の並び順は前者と完全に同じです。
set.seed(42)
marker <- rnorm(2000)
event <- rbinom(2000, size = 1, prob = plogis(-1 + 1.2 * marker))
calibrated <- plogis(-1 + 1.2 * marker)
shifted <- plogis(-1 + 1.2 * marker + 0.8)
auc <- function(outcome, score) {
ranks <- rank(score); positives <- sum(outcome == 1)
(sum(ranks[outcome == 1]) - positives * (positives + 1) / 2) /
(positives * sum(outcome == 0))
}
brier <- function(p) mean((event - p)^2)
round(c(auc_calibrated = auc(event, calibrated), auc_shifted = auc(event, shifted),
brier_calibrated = brier(calibrated), brier_shifted = brier(shifted)), 4)#> auc_calibrated auc_shifted brier_calibrated brier_shifted
#> 0.7778 0.7778 0.1710 0.1917結果として、AUCはどちらも0.7778で、小数第4位まで一致しました。一方、Brierスコアは0.1710と0.1917で差が出ています。順序だけを見る指標は、予測確率がそろって高すぎることを検出しません。
治療の開始基準を予測確率で決める場面では、この差がそのまま判断へ影響します。並び順が同じでも、20%を閾値にしたときに選ばれる患者が変わるためです。
生存時間解析のCV
生存時間解析では、追跡期間とイベント指標を用いて予測モデルを学習します。ところが、打ち切り例の観測期間を、その人のイベント発生時間とみなして通常のRMSEを計算することはできません。
たとえば、ある時点のイベントリスクを評価する場合、その時点より前に打ち切られた対象者については、評価時点までのイベントの有無が未確定です。したがって、データを分ける方法に加えて、打ち切りの下で性能をどう測るかが必要になります。
時点と予測量の固定
そこで、Coxモデルから将来の生存確率を予測する場合には、係数だけでなく、学習集合で推定したベースライン生存関数も使います。評価集合のイベント情報からベースラインを推定し直して、その集合の予測確率を作る方法は、元のモデルを持ち出した評価と異なります。
また、異なるfoldのCoxモデルの線形予測子には、中心化の基準が異なる場合があります。それらをそのまま全体のリスク順位へまとめる際には、モデル間の尺度を確認します。固定した時点の絶対リスクを使う場合にも、推定対象の時点を揃える必要があります。
さらに、Kaplan-Meier推定量の右端で観測情報が少なくなる問題は、予測性能の評価にも関わります。全期間を通じた一つの成績だけで、情報の乏しい遠い将来まで評価できたとは記述しません。
打ち切りを反映した評価指標
この課題について、Unoらは、打ち切り分布への依存を抑えて、リスクスコアの識別性能を評価するC統計量を提案しています。また、Grafらは、打ち切りを伴う生存確率の予測に対する誤差評価を扱っています。いずれも、打ち切りに関する仮定や、評価できる時間範囲が必要です。
そして、逆確率重み付けで性能を推定する場合には、打ち切りモデルをどのデータから推定するかも、評価手順の定義に含めます。評価指標を推定するために評価集合の打ち切り情報を使うことと、その情報で予測モデルを作り直すことは区別します。
また、競合イベントを含む研究で、何の発生確率を予測するかによって、正解の定義と性能指標が変わります。通常の生存確率、原因別ハザード、累積発生関数を、同じ予測量として採点しないようにします。
fold間の標準偏差と信頼区間
5個のfoldの成績が得られると、その標準偏差を計算し、5の平方根で割って平均の標準誤差を作りたくなる場合があります。しかし、CVのfold間では学習データが重なっているため、独立した5標本の平均と同じ分散計算にはなりません。
ここで、各foldの誤差を$\widehat R_k$とすると、その単純平均の分散は次の形です。
$$\operatorname{Var}\!\left( \frac{1}{K}\sum_{k=1}^{K}\widehat R_k \right) = \frac{1}{K^2} \left\{ \sum_{k=1}^{K}\operatorname{Var}(\widehat R_k) + 2\sum_{k<\ell} \operatorname{Cov}(\widehat R_k,\widehat R_\ell) \right\} \tag{3}$$
式(3)は、平均の分散に、各foldの分散だけでなく、fold同士の共分散も入ることを表しています。学習集合が重なる構成では、この共分散を0とみなす根拠がありません。
この問題について、BengioとGrandvaletは、K-fold CVの分散について、あらゆる分布で通用する普遍的な不偏推定量が存在しないことを証明しています。この結果は、どのような条件でも不確実性を評価できないという意味ではなく、単純な公式を万能の方法として使えないことを表しています。
不確実性を評価する対象
予測誤差の信頼区間を作る前に、学習済みの一つのモデルを固定した性能を知りたいのか、データを集め直して開発し直す手順の性能を知りたいのかを分けます。ここを混同すると、再標本化している対象と、結論で述べる対象が異なってしまいます。
たとえば、完全に独立した評価データで、学習済みモデルを固定して評価する場合には、評価対象者の標本変動に対応する不確実性を扱えます。ただし、その後に評価データを学習へ加えてモデルを更新すると、元の区間はその更新後モデルを直接評価したものではなくなります。
CVを包む再標本化
学習手順の変動まで含めたい場合には、元の対象者の単位で再標本化し、その内側で前処理、チューニング、学習、評価をやり直す方法などが候補になります。
ここで、ブートストラップの中でCVを組む場合、同じ元対象者の複製が学習側と評価側へ分かれると、その対象者の情報を共有します。そのため、再標本化された行番号だけで分割せず、元の対象者の同一性をどう扱うかまで設計する必要があります。
また、OOF予測値を固定して、その観測値と予測値の組だけをブートストラップする方法も、元データから学習し直す手順とは異なります。さらに、OOF予測同士には学習データの共有に由来する依存があるため、単に組を再標本化すれば学習全体の不確実性が得られるとは言えません。
同じnested CVという名称
紛らわしい名称として、Batesらは、CVの誤差推定の不確実性を評価するための、入れ子型の手法を提案しています。名称はnested CVですが、ハイパーパラメータの選択と最終評価を分ける通常の二重CVを一度実行するだけで、その区間推定法を使ったことにはなりません。
そして、今回のR例ではモデル選択と評価を分離していますが、Batesらの区間推定法、CV全体を包むブートストラップ、独立標本による信頼区間の作成は行っていません。したがって、本文の3.115というRMSEに95%信頼区間は付けていません。
成績の幅をどう書くか
ここまでで扱った指標は、どれも一つの数字に要約されます。ただし、その数字がどれくらい動くかは、fold間の標準偏差をそのまま標準誤差として扱っても求まりません。
そもそも、五つのfoldは同じ300人を共有していて、独立な五つの研究ではないためです。幅を書きたい場合には、CVの外側を再標本化で包むなど、何の変動を評価したいのかに合わせた構成を先に決めます。
次の記事では、複数のモデルを比べるとき、最終モデルを全データで学習するとき、そして論文の方法と結果に何を書くかを扱います。
参考文献
Van Calster B, McLernon DJ, van Smeden M, Wynants L, Steyerberg EW. Calibration: the Achilles heel of predictive analytics. BMC Medicine. 2019;17(1):230.
Uno H, Cai T, Pencina MJ, D’Agostino RB, Wei LJ. On the C-statistics for evaluating overall adequacy of risk prediction procedures with censored survival data. Statistics in Medicine. 2011;30(10):1105-1117.
Graf E, Schmoor C, Sauerbrei W, Schumacher M. Assessment and comparison of prognostic classification schemes for survival data. Statistics in Medicine. 1999;18(17-18):2529-2545.
Bengio Y, Grandvalet Y. No unbiased estimator of the variance of K-fold cross-validation. Journal of Machine Learning Research. 2004;5:1089-1105.
Bates S, Hastie T, Tibshirani R. Cross-validation: what does it estimate and how well does it do it? Journal of the American Statistical Association. 2024;119(546):1434-1445.
