マルチステートモデル（Multi-state Model）とは ― illness-death モデルと状態遷移確率をRで実装する ―
この記事でわかること
この記事を読むと、次のことがわかります。
- マルチステートモデル（Multi-state Model）とは何か ― 患者が複数の「状態」を行き来する過程を、遷移（transition）としてモデル化する枠組みであること
- illness-death モデルの考え方 ― 「初期状態 → 発症 → 死亡」と「初期状態 → 死亡」という2経路で、病気の”途中経過”を表現できること
- 標準的な生存時間解析・競合リスクとの違い ― 状態が2つや「1回のイベント」に限られないため、発症後の予後まで一気通貫で扱えること
- 遷移確率と Aalen-Johansen 推定量 ― 各時点で「患者がどの状態にいる確率」をノンパラメトリックに推定する、Kaplan-Meier の多状態への一般化
- Rでの実装イメージ ―
survivalパッケージのmgus2データを使い、illness-death モデルを実際に解析する流れ
はじめに
がんや慢性疾患の臨床経過を考えると、患者は「生きている／死んでいる」の2状態だけを行き来するわけではありません。多くの場合、患者は「診断された初期状態」から「合併症や別の疾患を発症した状態」を経て「死亡」へと至る、というように複数の状態を順に経由していきます。
標準的なイベント時間解析（Kaplan-Meier 法や Cox 回帰）は、生存か死亡かという2状態と「最初に起きた1つのイベント」だけに注目します。競合リスク解析はイベントの種類を区別できますが、やはり「最初に起きた1回のイベント」で観察が終わってしまい、その後の予後は追いかけられません。つまりどちらも、発症してから死亡するまでといった疾患の”途中経過”を捉えきれないのです。
そこで役立つのがマルチステートモデルです。本記事では、製薬企業の生物統計担当者やRで生存時間解析を一歩進めたい実務家、大学院生の方を対象に、illness-death モデルを軸に基本概念と数理を解説し、後半で mgus2 データを用いたRでの実装までを示します。
マルチステートモデルとは ― 状態と遷移で「途中経過」を捉える
状態（state）と遷移（transition）
マルチステートモデルでは、患者が取りうる臨床的な立場を状態（state）として定義し、状態間の移り変わりを遷移（transition）として表現します。状態全体の集合を状態空間（state space）と呼びます。
状態は、大きく2種類に分けられます。
- 推移状態（transient state）：そこから別の状態へ出ていける状態。例えば「初期状態」や「発症した状態」がこれにあたります。
- 吸収状態（absorbing state）：一度入ると二度と出ていかない状態。典型例は「死亡」で、そこから他の状態へ移ることはありません。
代表的な例が illness-death モデル（疾病-死亡モデル） です。これは次の3状態・2経路からなります。
- 状態1「初期状態（診断時）」→ 状態2「PCM発症（形質細胞悪性腫瘍の発症など）」：病気の”途中経過”を表す遷移
- 状態2「PCM発症」→ 状態3「死亡」：発症した後に死亡へ至る遷移
- 状態1「初期状態」→ 状態3「死亡」：発症を経ずに死亡する遷移
このように、初期状態と発症は推移状態、死亡は吸収状態として、遷移の「向き」まで含めて臨床経過を描けるのがマルチステートモデルの特徴です。
標準生存解析・競合リスクとの違い
マルチステートモデルの位置づけを理解するには、標準的な生存時間解析や競合リスク解析と並べて比べるのが近道です。下の表に、3つの手法の違いを整理します。
| 手法 | 状態数 | 捉えられるもの | 代表的推定量 | 臨床での例 |
|---|---|---|---|---|
| 標準生存解析 | 2状態（生存・死亡） | 単一イベントまでの時間 | Kaplan-Meier 生存曲線 | 全死亡までの時間 |
| 競合リスク | 複数の吸収状態（1回のイベント） | 種類の異なる最初のイベント | 累積発生関数（CIF） | がん死 vs 他病死の区別 |
| マルチステート | 3状態以上（推移＋吸収） | 状態間の移り変わり（途中経過を含む） | 遷移確率・Aalen-Johansen 推定量 | 発症 → 死亡の経過全体 |
競合リスクは「最初に起きたイベントが何か」を区別できますが、そのイベントが起きた時点で観察は終わります。一方マルチステートモデルは、発症という推移状態を経てから死亡という吸収状態へ至る過程まで一続きで扱えます。競合リスクは、実はマルチステートモデルにおいて「1つの初期状態から複数の吸収状態へ枝分かれする」特別な場合とみなせます。
遷移強度・遷移確率と Aalen-Johansen 推定量
ここからは、マルチステートモデルを支える数理を見ていきます。
まず、状態 \(h\) から状態 \(k\) へ移る「瞬間的な移りやすさ」を表すのが遷移強度（transition intensity）です。
\[
\alpha_{hk}(t) = \lim_{\Delta t \to 0} \frac{P\big(X(t+\Delta t)=k \mid X(t)=h\big)}{\Delta t}
\]
ここで \(X(t)\) は時刻 \(t\) に患者がいる状態を表します。この式は、「今 状態 \(h\) にいる患者が、ごく短い時間 \(\Delta t\) のあいだに状態 \(k\) へ移る確率」を \(\Delta t\) あたりに直したもので、標準生存解析のハザード関数を状態のペアごとに拡張した概念です。
次に、時刻 \(s\) に状態 \(h\) にいた患者が、時刻 \(t\)（\(s \le t\)）に状態 \(k\) にいる確率を要素にもつのが遷移確率行列（transition probability matrix） \(P(s,t)\) です。
\[
P_{hk}(s,t) = P\big(X(t)=k \mid X(s)=h\big)
\]
この \(P_{hk}(s,t)\) は、「時刻 \(s\) の状態 \(h\) を出発点として、時刻 \(t\) までにどの状態へ落ち着くか」を確率で表します。各行の合計は必ず1になります（どこかの状態には必ずいるため）。
ここで重要になるのがマルコフ仮定（Markov assumption）です。これは「次にどの状態へ移るかは、今いる状態だけで決まり、過去にどの状態をどれだけの時間かけて通ってきたかには依存しない」という仮定です。この仮定のもとでは、遷移確率行列は遷移強度からプロダクト積分（product integral）によって次のように組み立てられます。
\[
P(s,t) = \prod_{(s,\,t]} \big\{ I + d A(u) \big\}
\]
ここで \(I\) は単位行列、\(A(u)\) は各遷移の累積遷移強度（Nelson-Aalen 型の累積ハザード）を並べた行列です。この式は「微小な時間区間ごとの遷移をすべて掛け合わせて（積み上げて）、区間 \((s,t]\) 全体の遷移確率を作る」という意味を持ちます。
これを実データから推定したものが Aalen-Johansen 推定量（Aalen-Johansen estimator） です。各時点でのNelson-Aalen累積強度をプロダクト積分することで、各状態の占有確率（state occupation probability、その時点で患者がその状態にいる確率）を推定します。状態が2つ（生存・死亡）だけの場合、Aalen-Johansen 推定量はちょうど Kaplan-Meier 推定量に一致します。つまり Aalen-Johansen 推定量は、Kaplan-Meier 法を多状態へ一般化したものと理解できます。
上記の遷移確率の組み立てはマルコフ仮定に依存しています。しかし現実には、「発症してから死亡へ移りやすいかどうかは、発症してからの滞在時間（経過時間）に依存する」ことが少なくありません。このように滞在時間が遷移強度に影響する場合、単純なマルコフ仮定は成り立たず、半マルコフモデル（semi-Markov model）など、より柔軟な枠組みの検討が必要になります。解析結果を解釈する際は、マルコフ仮定が妥当かどうかを常に意識することが極めて重要です。
Rによる実装 ― mgus2データでillness-deathモデルを解析する
ここからは実際にRでillness-deathモデルを動かします。題材は survival パッケージに同梱の mgus2（n=1384、単クローン性ガンマグロブリン血症＝MGUSの患者コホート）です。この集団では、一部が形質細胞悪性腫瘍（PCM）を発症し、その前後を問わず高齢のため多くが追跡中に死亡します。まさに「初期状態 → PCM発症 → 死亡」「初期状態 → PCMを経ない死亡」という2経路を持つillness-deathモデルの典型例です。
データの準備と状態遷移の確認
まず発症（PCM）と死亡を1つの時間軸に統合し、イベントを因子として定義します。
library(survival)
data(cancer, package = "survival") # mgus2 を含む
# 発症(PCM)と死亡を1つの時間軸に統合してイベント因子を作る
mgus2$etime <- with(mgus2, ifelse(pstat == 0, futime, ptime))
mgus2$event <- with(mgus2, ifelse(pstat == 0, 2 * death, 1)) # 0=打切, 1=PCM, 2=死亡
mgus2$event <- factor(mgus2$event, 0:2, labels = c("censor", "pcm", "death"))
table(mgus2$event)
> censor pcm death
> 409 115 860
続いて survcheck() で状態遷移の件数を確認します。id=id で1人1IDとして遷移を数えます。
survcheck(Surv(etime, event) ~ 1, id = id, data = mgus2)
> Transitions table:
> to
> from pcm death (censored)
> (s0) 115 860 409
> pcm 0 0 0
> death 0 0 0
1384人のうち、初期状態(s0)からPCMを発症したのが115人、PCMを経ずに死亡したのが860人、打切が409人です。ここで注目したいのは
pcm→death が 0 になっている点です。この単純なコーディングでは、PCMを発症した時点で追跡が打ち切られ、PCMが「吸収状態」として扱われます。つまり実質的には競合リスク（発症 vs 死亡）の構造で、「PCM発症後の死亡」という遷移を捉えられていません。真のillness-deathモデルにするには、発症後も死亡リスクにさらされ続けるようにデータを行分割する必要があります。そこで tmerge() を使い、1人を「発症前」「発症後」の複数行に分割してPCM→死亡の遷移を明示的に表現します。
# PCMと死亡が同月のタイは、PCMをわずかに前へずらして順序を確定
ptemp <- with(mgus2, ifelse(ptime == futime & pstat == 1, ptime - 0.1, ptime))
iddata <- tmerge(mgus2, mgus2, id = id,
death = event(futime, death),
pcm = event(ptemp, pstat))
iddata$event <- with(iddata, factor(pcm + 2 * death, 0:2,
labels = c("censor", "pcm", "death")))
survcheck(Surv(tstart, tstop, event) ~ 1, id = id, data = iddata)
> Transitions table:
> to
> from pcm death (censored)
> (s0) 115 860 409
> pcm 0 103 12
> death 0 0 0
行分割後は
pcm→death=103 という遷移が現れました。PCMを発症した115人のうち、103人がその後に死亡し、12人はPCM状態のまま打切になったという読み方です。(s0)→死亡の860は変わりません。これで「初期 → PCM → 死亡」と「初期 → 死亡」の2経路をもつ、正真正銘のillness-deathモデルのデータが完成しました。以降の解析はこの行分割済みデータ iddata を使います。Aalen-Johansen推定による状態占有確率
マルチステートモデルの主役は「ある時点で各状態にいる確率」＝状態占有確率です。survfit() にマルチステートのデータを渡すと、ノンパラメトリックなAalen-Johansen推定量が計算されます。性別で層別してみます。
fit_aj <- survfit(Surv(tstart, tstop, event) ~ sex, data = iddata, id = id)
# 120ヶ月(10年)・240ヶ月(20年)時点の各状態占有確率
summary(fit_aj, times = c(120, 240))
> sex=F
> time n.risk Pr((s0)) Pr(pcm) Pr(death)
> 120 223 0.446 0.0170 0.537
> 240 34 0.200 0.0064 0.794
>
> sex=M
> time n.risk Pr((s0)) Pr(pcm) Pr(death)
> 120 215 0.370 0.0079 0.623
> 240 27 0.156 0.0161 0.828
各時点で3状態の確率は必ず合計1になります（例：女性120ヶ月は \(0.446+0.017+0.537 \approx 1.00\)）。数式で書けば任意の時点 \(t\) で \(\sum_{k} P_k(t)=1\) です。ここで極めて重要なのは、Pr(pcm) が1〜2%と小さいこと。これは「累積で何割が発症したか」ではなく「その時点でPCM状態に留まっている確率（有病割合）」だからです。PCM発症後の死亡が速いため、PCM状態への滞在は短く、占有確率は一度上がってすぐ下がります。この「上がって下がる曲線」は、単調増加する競合リスクの累積発生関数とは本質的に異なる、マルチステートモデルならではの視点です。一方、死亡状態の占有確率は240ヶ月（20年）で女性0.79・男性0.83に達します。高齢コホートらしく大半が死亡へ吸収され、各時点で男性の方が死亡占有が高い（予後がやや不良）ことも読み取れます。
占有確率は図にすると直感的です。状態を1つ選んで性別で比較します。
plot(fit_aj[, "death"], col = c("#2E86C1", "#E67E22"), lwd = 2,
xlab = "診断からの月数", ylab = "死亡状態の占有確率")
legend("topleft", c("女性", "男性"), col = c("#2E86C1", "#E67E22"), lwd = 2)
mstateパッケージによる遷移別ハザード
「どの共変量が、どの遷移に、どれだけ効くか」を定量化するには、遷移別のCoxモデルが有効です。mstate パッケージの transMat() で遷移行列を定義し、msprep() でロング形式に変換してから、遷移を層別因子とした1本のCoxモデルを当てはめます。
library(mstate)
# 遷移行列：1=初期→PCM, 2=初期→死亡, 3=PCM→死亡
tmat <- transMat(x = list(c(2, 3), c(3), c()),
names = c("Diagnosis", "PCM", "Death"))
tmat
> to
> from Diagnosis PCM Death
> Diagnosis NA 1 2
> PCM NA NA 3
> Death NA NA NA
# 各遷移用の time / status を用意（タイは死亡を0.1ずらす）
mgus2$pcmtime <- with(mgus2, ifelse(pstat == 1, ptime, futime))
mgus2$deathtime <- mgus2$futime
tie <- mgus2$pstat == 1 & mgus2$pcmtime == mgus2$deathtime
mgus2$deathtime[tie] <- mgus2$deathtime[tie] + 0.1
ms <- msprep(time = c(NA, "pcmtime", "deathtime"),
status = c(NA, "pstat", "death"),
data = mgus2, trans = tmat, keep = c("age", "sex"))
events(ms)$Frequencies
> to
> from Diagnosis PCM Death no event total entering
> Diagnosis 0 115 860 409 1384
> PCM 0 0 103 12 115
> Death 0 0 0 963 963
遷移別イベント数は初期→PCM=115、初期→死亡=860、PCM→死亡=103で、survcheck の結果と一致します。続いて共変量（年齢・性別）を遷移別に展開し、遷移ごとにベースラインハザードを層別したCoxモデルを推定します。
# 共変量を遷移別に展開（age.1〜age.3, sexM.1〜sexM.3 が作られる）
ms <- expand.covs(ms, c("age", "sex"), append = TRUE)
cox_ms <- coxph(Surv(Tstart, Tstop, status) ~
age.1 + age.2 + age.3 +
sexM.1 + sexM.2 + sexM.3 + strata(trans),
data = ms)
summary(cox_ms)$conf.int
> exp(coef) exp(-coef) lower .95 upper .95
> age.1 1.0131 0.9870 0.9969 1.030
> age.2 1.0670 0.9372 1.0594 1.075
> age.3 1.0426 0.9591 1.0153 1.071
> sexM.1 0.9752 1.0255 0.6740 1.411
> sexM.2 1.4818 0.6749 1.2926 1.699
> sexM.3 1.0709 0.9338 0.7153 1.603
末尾の
.1 / .2 / .3 が遷移番号（1=初期→PCM、2=初期→死亡、3=PCM→死亡）です。年齢(age)は、初期→死亡（age.2、HR=1.067/年、95%CI 1.06–1.08、p<0.001）とPCM→死亡（age.3、HR=1.043/年、p=0.002）で有意に死亡ハザードを高める一方、初期→PCM発症（age.1、HR=1.013/年、p=0.11）への効果は弱く有意ではありません。診断1歳ごとに直接死亡のハザードが約6.7%上がる、という読み方です。性別は初期→死亡でのみ顕著で、男性は女性を基準にPCMを経ない死亡ハザードが約1.48倍（sexM.2、95%CI 1.29–1.70、p<0.001）。しかしPCM発症（sexM.1）やPCM後の死亡（sexM.3）には統計的な性差がほとんど見られません。「同じ共変量でも遷移ごとに効果の向き・大きさが違う」ことを1本のモデルで定量化できる点が、マルチステート回帰の最大の強みです。実務でのポイント
illness-deathモデルを現場で使いこなすうえで、時間軸の取り方・アウトカムの語り方・活用場面を整理しておきます。
① 時間軸（マルコフ vs セミマルコフ）を意識する。 遷移強度が「プロセス全体の時間 \(t\)」に依存するのがマルコフ（clock-forward：時計を進め続ける）、「状態に入ってからの経過時間」に依存するのがセミマルコフ（clock-reset：状態進入で時計をリセット）です。PCM後の予後が「発症からの月数」で決まると考えるならセミマルコフが自然で、Rでは
Surv(Tstop - Tstart, status) のように状態内時間を用います。② アウトカムを状態占有確率と状態内RMSTで語る。 ハザード比だけでなく、\(\mathrm{RMST}_k(\tau)=\int_0^{\tau} P_k(u)\,du\)（状態 \(k\) での制限付き平均滞在期間）を使うと、「治療がどの状態にどれだけ長く留まらせるか」で効果を表現できます。今回のデータでも
survfit の rmean から、女性は初期状態に平均約139ヶ月、男性は約123ヶ月留まる（女性の方が発症前の期間が長い）と読めます。③ 活躍場面。 腫瘍領域の「無増悪 → 増悪 → 死亡」（TTP・PFS・OSを一体で扱う）や、RWE（レジストリで疾患進行の途中経過を追う）で特に有用です。
④ サンプルサイズは遷移別イベント数で決まるため、設計段階で各遷移のイベント数を見積もることが欠かせません。
状態や経路を増やすほど、1つの遷移に割り当てられるイベント数は急激に少なくなります。今回のmgus2でもPCM→死亡は103イベントしかなく、この遷移に共変量を何個も投入すると推定が不安定になります。目安として「1変数あたり10イベント以上」を各遷移で確保できているかを必ず確認し、足りなければ遷移をまとめる・共変量を絞る・遷移間でパラメータを共有するといった対処を検討してください。
主要な時間軸の考え方は下表の通りです。
| 観点 | マルコフモデル（clock-forward） | セミマルコフモデル（clock-reset） |
|---|---|---|
| 遷移強度が依存する時間 | プロセス全体の時間 \(t\)（診断からの経過） | 現在の状態に入ってからの経過時間 |
| 基本の考え方 | 将来の遷移は「いつ発症したか」に関わらず現在時刻 \(t\) だけで決まる | PCM後の死亡リスクは「発症からの経過月数」で決まる |
| Rでの実装 | Surv(Tstart, Tstop, status) をそのまま使用 | Surv(Tstop - Tstart, status) で状態内時間を使用 |
| mgus2での相性 | 実装が簡単でAalen-Johansen推定と整合しやすい | 発症時期で予後が変わる場合に自然 |
この記事をより深く理解するための参考書籍
マルチステートモデルは、Cox回帰・競合リスク・臨床試験デザインの延長線上にあります。基礎から実務まで橋渡しになる書籍を3冊紹介します。
まとめ
本記事では、生存時間解析を「状態」と「遷移」に拡張するマルチステートモデルを、illness-deathモデルを軸に解説しました。Rの mgus2 データを使い、survcheck で状態遷移（初期→PCM 115、初期→死亡 860、PCM→死亡 103）を確認し、tmerge による行分割で単純な競合リスク構造を真のillness-deathモデルへ拡張しました。さらにAalen-Johansen推定で各時点の状態占有確率（合計は常に1）を求め、mstate の msprep＋Coxモデルで「年齢は死亡遷移に強く効き、男性はPCMを経ない死亡ハザードが約1.48倍」といった遷移別の効果を1本のモデルで定量化しました。
途中経過を丸ごと捉えられるマルチステートモデルは、腫瘍領域の複合エンドポイントやRWEでの疾患進行の記述において、単一イベントの解析では見えない構造を明らかにします。時間軸（マルコフ／セミマルコフ）と遷移別イベント数に注意しながら使えば、解析の表現力を大きく広げてくれます。
関連するトピックは、以下の記事もあわせてご覧ください。
- 競合リスク（Competing Risks）とは ― Fine-Grayモデルと累積発生関数をRで実装 ―
- Cox比例ハザード性のチェックとRMSTへの切り替え判断
- RMSTをRで実装する方法 ― survRM2による制限付き平均生存時間の解析 ―
- パラメトリック生存時間解析とは ― 指数・Weibull・AFTモデル ―
競合リスクやRMSTと組み合わせてマルチステートモデルを使いこなせれば、「イベントが起きたか」だけでなく「どの状態をどう経由したか」まで語れるようになり、臨床開発やRWE解析における確かな強みになります。ぜひご自身のデータでも、状態と遷移の視点からの解析を実現していただければと思います。
