Top 5 This Week

関連記事

【医療統計版】決定曲線分析:ネットベネフィットで臨床的有用性を測る

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

予測モデルの論文では、C統計量が0.75前後、較正プロットもおおむね対角線に沿っている、という性能評価まで書かれていても、そのモデルに従って治療する人を選んでよいかという段になると、判断の根拠が急に薄くなることがあります。識別も較正も予測値と実測値の対応を測る指標であるため、モデルに従って対象者を選別した結果として患者集団全体でどれだけ得をするかには直接答えることはできません。

そこで、決定曲線分析(decision curve analysis)は、この差し引きを閾値確率の関数として一本の曲線に描く方法で、Vickers and Elkin(2006)が提案しました。この記事を読み終えると、ネットベネフィット(net benefit)が何を足して何を引いた量なのか、閾値確率をどの範囲に取るのか、決定曲線のどこを見てモデルの採否を判断するのかまでを、自分の研究データに当てはめて考えられるようになります。

識別と較正だけでは届かない判断

再発リスクの予測モデルを作り、リスクが高いと判定された患者にだけ追加の検査を行う、という方針を考えます。この方針の良し悪しは、見つかった再発例の数だけでは決まらず、再発しない患者に不要な検査を受けさせた数と合わせて初めて評価できます。

ところが、C統計量が答えるのは、再発した患者と再発しなかった患者を1人ずつ無作為に取り出したとき、前者により高いリスクを与える確率であり、どのリスク値を境に検査へ回すかという判断そのものを含みません。較正は予測確率の水準が観測割合と合っているかを測りますが、こちらも一致の度合いを表す量であって、検査を何人に行うかの損得へは換算されません。

この区別について、Steyerbergら(2010)は、識別、較正、臨床的有用性を予測モデルの性能の別々の軸として整理しています。C統計量そのものの読み方と検証の方法は予測モデルの検証:C統計量0.75は本当に0.75かで扱っており、ここではその先にある第三の軸だけを追います。

ネットベネフィットという一つの物差し

ネットベネフィットは、モデルに従って治療した人のうち本当にイベントを起こした人の割合から、治療してしまった非イベントの割合に重みを掛けて引いた量です。治療へ回す境目となるリスクを$t$と書くと、次の式になります。

$$\mathrm{NB}(t) = \frac{TP}{n} – \frac{FP}{n}\cdot\frac{t}{1-t} \tag{1}$$

式(1)は、対象者1人あたり何人分の真の陽性を、余計な治療を差し引いた後に手元へ残せるかを表しています。$n$は対象者数、$TP$は閾値以上と判定されてイベントを起こした人数、$FP$は閾値以上と判定されてイベントを起こさなかった人数です。

この引き算で効いているのは重み$t/(1-t)$で、これがあるために真陽性1人と偽陽性1人を同じ重さで足し引きせずに済みます。単位が真陽性の人数に揃っているので、閾値を固定すれば、異なるモデルの値も、モデルを使わない方針の値も、同じ目盛りの上で比べられます。

閾値確率が背負っている害と利益の比

閾値確率を、感度と特異度の釣り合いが最も良くなる点として選ぶ、と考えたくなりますが、決定曲線分析での閾値は最適化して決める量ではなく、治療を受け入れる境目として臨床的に決めたリスク値がそのまま入ります。閾値を$t$に置くことは、治療で避けられるイベント1件と、不要な治療$(1-t)/t$件とを等価とみなす判断を宣言したことに当たります。

たとえば閾値を0.10に置いた場合、重みは0.10を0.90で割った約0.111となり、イベント1件を防ぐために不要な治療を9件まで許容する、という判断を表します。閾値を0.30へ上げると重みは約0.429まで大きくなり、許容できる不要な治療は2件強まで減ります。

そのため、閾値の範囲は解析を始める前に決めておく必要があり、Vickersら(2019)も、臨床的に妥当な閾値の幅を先に定め、その範囲の中だけで曲線を比べる手順を勧めています。私なら、治療の侵襲性と見逃したときの重大さから上下2つの値を決め、そのあいだを報告範囲にします。

閾値を決める作業は、この記事で扱う計算のうち、統計ではなく臨床が決める唯一の部分です。

全員治療と誰も治療しないという二本の基準線

決定曲線には、モデルの曲線に加えて、全員を治療する方針と誰も治療しない方針を必ず重ねます。誰も治療しなければ真陽性も偽陽性も0になるため、そのネットベネフィットはどの閾値でも0の水平線になります。

一方、全員を治療する方針のネットベネフィットは、イベントを起こした人の割合$p$と閾値だけで書けます。

$$\mathrm{NB}_{\mathrm{all}}(t) = p – (1-p)\cdot\frac{t}{1-t} \tag{2}$$

式(2)は、全員治療の値が閾値確率とともに単調に下がり、閾値が$p$と一致したところでちょうど0を通ることを表しています。それより高い閾値では負になり、誰も治療しない方針に追い越されます。

この2本があるおかげで、モデルの曲線が高いか低いかを絶対値で悩まずに済み、判断は、対象とする閾値の範囲でモデルの曲線が両方の基準線より上にあるかどうかへ集約されます。どちらかの基準線を下回る閾値では、そのモデルに従うより、全員治療か経過観察のどちらかを選んだほうが集団全体の損得は良くなります。

2本の基準線を描かない決定曲線は、比べる相手のいないグラフになります。

Rで描く決定曲線

ここからは、ランダムに生成したデータセットで実際に曲線を描きます。年齢とバイオマーカーの2変数からイベント発生の確率が決まる構造を置き、開発用2,000例でロジスティック回帰を当てはめ、別に生成した検証用2,000例で予測確率を計算しました。

これに対して、同じ線形予測子に定数0.8を足しただけのモデルも比較対象として用意しました。順位を変えずに全員のリスクを高く見積もる変換なので、識別は完全に同じまま較正だけがズレる、という状況を作れます。

set.seed(42)
age_years <- rnorm(4000, mean = 65, sd = 10)
biomarker <- rnorm(4000, mean = 0, sd = 1)
event <- rbinom(4000, size = 1, prob = plogis(-2.0 + 0.05 * (age_years - 65) + 0.90 * biomarker))
dev_data <- data.frame(event, age_years, biomarker)[1:2000, ]
val_data <- data.frame(event, age_years, biomarker)[2001:4000, ]
fit_model <- glm(event ~ age_years + biomarker, family = binomial(), data = dev_data)
lp_val <- predict(fit_model, newdata = val_data, type = "link")
risk_calibrated <- plogis(lp_val)
risk_overestimated <- plogis(lp_val + 0.8)
outcome <- val_data[["event"]]

検証用データでは、イベントが307件、イベント割合は0.1535でした。この予測確率から式(1)をそのまま計算したのが次の部分で、全員治療の値は、全員のリスクを1とみなした場合として同じ関数から得ています。

# 閾値thresholdのもとでのネットベネフィット
net_benefit <- function(risk, outcome, threshold) {
  weight <- threshold / (1 - threshold)
  treated <- risk >= threshold
  (sum(treated & outcome == 1) - weight * sum(treated & outcome == 0)) / length(outcome)
}
show_row <- function(threshold) cat("threshold", formatC(threshold, format = "f", digits = 2),
  "| NB A:", formatC(net_benefit(risk_calibrated, outcome, threshold), format = "f", digits = 5),
  "| NB B:", formatC(net_benefit(risk_overestimated, outcome, threshold), format = "f", digits = 5),
  "| treat all:", formatC(net_benefit(rep(1, 2000), outcome, threshold), format = "f", digits = 5), "\n")
invisible(sapply(c(0.05, 0.15, 0.20), show_row))
#> threshold 0.05 | NB A: 0.10947 | NB B: 0.11018 | treat all: 0.10895
#> threshold 0.15 | NB A: 0.05526 | NB B: 0.03826 | treat all: 0.00412
#> threshold 0.20 | NB A: 0.03988 | NB B: 0.01275 | treat all: -0.05813

続いて、2つのモデルの識別と較正切片を並べました。較正切片は、線形予測子をオフセットに固定したときの切片として計算しています。

c_stat <- function(risk, event_flag) (sum(rank(risk)[event_flag == 1]) - sum(event_flag) * (sum(event_flag) + 1) / 2) / (sum(event_flag) * sum(1 - event_flag))
round(c(c_stat(risk_calibrated, outcome), c_stat(risk_overestimated, outcome)), 4)
round(unname(sapply(list(lp_val, lp_val + 0.8), function(lp) coef(glm(outcome ~ 1, family = binomial(), offset = lp)))), 4)
#> [1] 0.7456 0.7456
#> [1] -0.0034 -0.8034

出力の読み方の中心は、C統計量が両モデルとも0.7456で完全に一致している一方、較正切片はモデルAの-0.0034に対してモデルBが-0.8034と大きく離れている点です。切片が負へ振れているのは、予測リスクが観測割合より高い側へズレていることに対応します。

読み違えやすいのは閾値0.05での値で、ここではモデルBの0.11018がモデルAの0.10947をわずかに上回ります。ただし、この閾値では両モデルとも大半を治療側へ回すため、全員治療の0.10895との差もほとんど残りません。

thresholds <- seq(0.01, 0.50, by = 0.01)
nb_calibrated <- sapply(thresholds, function(t) net_benefit(risk_calibrated, outcome, t))
nb_overestimated <- sapply(thresholds, function(t) net_benefit(risk_overestimated, outcome, t))
nb_treat_all <- sapply(thresholds, function(t) net_benefit(rep(1, 2000), outcome, t))
plot(thresholds, nb_calibrated, type = "l", lwd = 2.5, col = "#1b5e9e",
     ylim = c(-0.05, max(nb_treat_all, nb_calibrated) + 0.01),
     xlab = "Threshold probability", ylab = "Net benefit",
     main = "Decision curves for two models with identical C statistic")
lines(thresholds, nb_overestimated, lwd = 2.5, col = "#c0562a", lty = 2)
lines(thresholds, nb_treat_all, lwd = 2, col = "#777777")
abline(h = 0, lwd = 2, col = "#333333", lty = 3)
legend("topright", bty = "n", cex = 0.85,
       legend = c("Model A (calibrated)", "Model B (risk overestimated)", "Treat all", "Treat none"),
       col = c("#1b5e9e", "#c0562a", "#777777", "#333333"),
       lty = c(1, 2, 1, 3), lwd = c(2.5, 2.5, 2, 2))
横軸が閾値確率、縦軸がネットベネフィットの折れ線グラフ。C統計量が同じ2つのモデルの曲線が、中程度の閾値で上下に離れている。
図1 ランダムに生成した検証用2,000例で描いた決定曲線。実線が較正のとれたモデルA、破線がリスクを一律に高く見積もるモデルB、灰色の線が全員治療、点線が誰も治療しない方針の値。閾値0.15から0.30の範囲で2本のモデル曲線がはっきり離れる。

C統計量が同じ二つのモデルの分かれ方

閾値0.15での値を並べると、モデルAが0.05526、モデルBが0.03826、全員治療が0.00412でした。識別が完全に同じ2つのモデルのあいだに、真陽性の人数へ直して100人あたり1.7人分の差が付いています。

この差の出どころは治療側に回る人数で、閾値0.20ではモデルAが25.3%、モデルBが58.65%を陽性と判定していました。イベント割合が0.1535の集団で6割近くを治療側へ回せば、増えた分の大半は偽陽性になり、重み0.25を掛けて引かれます。

そのため、閾値0.30ではモデルBのネットベネフィットが-0.00700まで下がり、誰も治療しない方針の0を下回りました。較正のズレが、順位付けの良さを一切損なわないまま、方針としての価値だけを消してしまった状態です。

ここで注意が必要なのは、この0.04前後という差の大きさが、今回置いた0.8という較正のズレと、この検証用データ1組から出た結果にすぎないという点です。ズレの向きが逆であれば、同じ閾値でも治療側に回る人数は減り、曲線の動き方も変わります。

識別の指標は1つの数字しか返しませんが、決定曲線は閾値ごとの損得の形を返します。

研究での位置づけと報告に求められる記載

決定曲線分析は、Vickersら(2016)が予測モデル、分子マーカー、診断検査を共通の枠組みで評価する方法として整理して以降、臨床予測研究で広く使われるようになりました。Collinsら(2024)のTRIPOD+AI声明でも、予測モデル研究の報告項目として臨床的有用性の評価が求められています。

報告のときに欠かせない記載として、Van Calsterら(2018)は、検討した閾値確率の範囲、比べた方針、そして曲線を描いたデータが開発データか検証データかを求めています。同じ論文は較正のズレた予測が決定曲線を歪めることも指摘しており、今回の計算で起きたのはまさにその状態です。

開発データで描いた曲線については、内部検証と同じ楽観性の問題が付きまといます。ブートストラップによる楽観性の補正は予測モデルの内部検証と報告:楽観性・欠測・生存時間解析で扱っており、決定曲線も同じ枠組みで補正の対象になります。

決定曲線分析を勧める場面と勧めない場面

私は、予測モデルで治療や追加検査の対象を絞る使い方を想定した研究では、C統計量と較正プロットに加えて決定曲線を載せることを勧めます。モデルの採否が最終的には治療する人数の損得で決まる以上、その損得を直接に計算した量を報告しない理由が見当たりません。

逆に、臨床的に妥当な閾値の範囲を誰も決められない状況で、0.01から0.99まで広く曲線を描いて見せる報告は適切ではありません。全域で他の方針を上回るモデルはまれで、範囲を先に決めずに描くと、自分に有利な閾値だけを後から選ぶ余地が残ります。

もう一点、較正が大きくズレたモデルについて、決定曲線が低いことを理由に予測因子そのものを捨てる、という結論も適切ではありません。今回のモデルBは線形予測子の切片を戻せばモデルAと一致する作りで、順位付けの情報は何も失われていないため、まず較正の再当てはめを試したほうが判断を誤りにくくなります。

自分の研究データで確かめる順序

自分の研究で決定曲線を使うかどうかは、そのモデルが誰かの治療方針を変える設計になっているかで決まります。変える設計であれば、臨床的に許容できる閾値の幅を先に定め、その範囲の中だけで全員治療と誰も治療しない方針の2本と比べます。

その比較で範囲の全域にわたってモデルが両方の基準線より上にあれば、ネットベネフィットの差を100人あたりの不要な治療の削減数へ直して報告すると、読者に届きやすくなります。今回のデータでは、閾値0.20でモデルAは全員治療と同じ真陽性数を保ったまま、100人あたり39.2人の不要な治療を減らしていました。

ただし、較正がズレたモデルでは、この人数はモデルの識別能力ではなく予測水準のズレを反映して動きます。曲線を描く前に較正切片と較正勾配を確かめておくと、決定曲線の解釈で迷う場面が減ります。

ここまで読んで、では自分の研究では閾値の範囲をどう定めればよいのか、と思われた方もいるかもしれません。実際の解析では、対象集団のイベント割合、治療の侵襲性、検証データの有無、そして事前に定めた解析計画まで含めて決めることになり、Dr.データサイエンスでは、こうした報告方針のご相談を承っています。

参考文献

  1. Vickers AJ, Elkin EB. Decision Curve Analysis: A Novel Method for Evaluating Prediction Models. Medical Decision Making. 2006;26(6):565-574.
  2. Vickers AJ, Van Calster B, Steyerberg EW. Net benefit approaches to the evaluation of prediction models, molecular markers, and diagnostic tests. BMJ. 2016;352:i6.
  3. Vickers AJ, van Calster B, Steyerberg EW. A simple, step-by-step guide to interpreting decision curve analysis. Diagnostic and Prognostic Research. 2019;3(1):18.
  4. Van Calster B, Wynants L, Verbeek JFM, Verbakel JY, Christodoulou E, Vickers AJ, Roobol MJ, Steyerberg EW. Reporting and Interpreting Decision Curve Analysis: A Guide for Investigators. European Urology. 2018;74(6):796-804.
  5. Steyerberg EW, Vickers AJ, Cook NR, Gerds T, Gonen M, Obuchowski N, Pencina MJ, Kattan MW. Assessing the Performance of Prediction Models: A Framework for Traditional and Novel Measures. Epidemiology. 2010;21(1):128-138.
  6. Collins GS, Moons KGM, Dhiman P, Riley RD, Beam AL, Van Calster B, Ghassemi M, Liu X, Reitsma JB, van Smeden M, et al. TRIPOD+AI statement: updated guidance for reporting clinical prediction models that use regression or machine learning methods. BMJ. 2024;385:e078378.

Popular Articles