カプラン＝マイヤー推定量

カプラン=マイヤー推定量（カプラン マイヤー すいていりょう、英: Kaplan-Meier estimator）は^[1][2]、積極限推定量（せききょくげんすいていりょう、英: product limit estimator）とも呼ばれ、生存データから生存関数を推定するために用いられるノンパラメトリック統計量である。医学研究では、治療後に一定期間生存している患者の割合を測定するためによく使われる。他の分野では、カプラン=マイヤー推定量を用いて、失業後に人々が失業している期間の長さや^[3]、機械部品の故障までの時間や、果実食動物（英語版）に食べられてしまうまでの肉果の残存期間を測定することができる。この推定量は、Edward L. KaplanとPaul Meierが米国統計学会誌（Journal of the American Statistical Association）に別々に原稿を提出したことにちなんで命名された^[4]。ジャーナル編集者のジョン・テューキーは、彼らの研究を1つの論文にまとめるよう説得した。この論文は1958年に発表されて以来、約61,000回も引用されている^[5][6]。

患者の生存に関連する2つの条件についてのカプラン=マイヤープロットの一例。

その生存関数 ${\displaystyle S(t)}$ （寿命が ${\displaystyle t}$ より長くなる確率）の推定量は次の式で与えられる。

${\displaystyle {\widehat {S}}(t)=\prod \limits _{i:\ t_{i}\leq t}\left(1-{\frac {d_{i}}{n_{i}}}\right),}$

ここに、 ${\displaystyle t_{i}}$ は少なくとも1つのイベントが発生した時刻、${\displaystyle d_{i}}$ は時刻 ${\displaystyle t_{i}}$ で発生したイベントの数（たとえば、死亡）、そして ${\displaystyle n_{i}}$ は時刻 ${\displaystyle t_{i}}$ まで生存していることが分かっている（まだイベントが発生していないか、打ち切られていない）個体の数である。

基本的な考え方

カプラン=マイヤー推定量のプロットは、一連の減少する水平ステップの系列であり、十分に大きな標本サイズの時に、その母集団の真の生存関数に近づく。連続する別個のサンプリングされた観測値（カチッと音がする）間の生存関数の値は一定であると仮定される。

カプラン=マイヤー曲線の重要な利点は、この手法がいくつかのタイプの打ち切りデータ、特に患者が研究から離脱した場合、またはフォローアップに失敗した場合、またはイベントなしで生存している場合に発生する「右側打ち切り」（right-censoring）を考慮に入れることができることである。プロット上では、小さな縦の目盛りが、生存時間が右側打ち切りされた個々の患者を示している。切り捨てや打ち切りが行われない場合、カプラン=マイヤー曲線は、経験分布関数（英語版）の補集合である。

医学統計学（英語版）では、一般的な応用例として、たとえば遺伝子Aプロファイルを持つ患者と遺伝子Bプロファイルを持つ患者のように、患者をカテゴリーに分類することがある。このグラフでは、遺伝子Bを持つ患者は、遺伝子Aを持つ患者よりも早く死亡する。2年後の生存率は、遺伝子Aの患者では約80%だが、遺伝子Bの患者では半分未満である。

カプラン=マイヤー推定量を作成するには、各患者（または各被験者）について、少なくとも2個のデータが必要である。それらは、最後の観察時の状態（イベント発生または右側打ち切り）と、イベント発生までの時間（または打ち切りまでの時間）の対からなる。2つ以上のグループ間の生存関数を比較する場合、3番目のデータが必要で、それは各被験者のグループ割り当てである^[7]。

問題の定義

確率変数 ${\displaystyle \tau }$（タウ）を、関心のあるイベントが起こるまでの時間と考えよう（ただし、${\displaystyle \tau \geq 0}$ である）。上に示したように、目的は、${\displaystyle \tau }$ の潜在的な生存関数 ${\displaystyle S}$ を推定することである。この関数は、

${\displaystyle S(t)=\operatorname {Prob} (\tau >t)}$

として定義され、${\displaystyle t=0,1,\dots }$ は時間であることを思い出そう。

${\displaystyle \tau _{1},\dots ,\tau _{n}\geq 0}$ を独立した同一分布の確率変数とし、その一般分布は ${\displaystyle \tau }$ で、${\displaystyle \tau _{j}}$ は、あるイベント ${\displaystyle j}$ が発生した確率的な時間としよう。${\displaystyle S}$ を推定するために利用できるデータは、${\displaystyle (\tau _{j})_{j=1,\dots ,n}}$ ではなく、ペアのリスト ${\displaystyle (\,({\tilde {\tau }}_{j},c_{j})\,)_{j=1,\dots ,n}}$ である。ここで、${\displaystyle j\in [n]:=\{1,2,\dots ,n\}}$ について、${\displaystyle c_{j}\geq 0}$ は固定の決定論的整数で、イベント ${\displaystyle j}$ の打ち切り時間（censoring time）であり、${\displaystyle {\tilde {\tau }}_{j}=\min(\tau _{j},c_{j})}$ である。特に、イベント ${\displaystyle j}$ のタイミングについて利用可能な情報は、イベントが固定時間 ${\displaystyle c_{j}}$ の前に起こったかどうかであり、もし起こった場合は、イベントの実際の時間も利用可能となる。課題は、このデータをもとに ${\displaystyle S(t)}$ を推定することである。

カプラン=マイヤー推定量の導出

ここでは、カプラン=マイヤー推定量の2つの導出方法を示す。どちらも生存関数を、ハザード（hazard）または死亡率（mortality rates）と呼ばれる観点で書き換えることに基づいている。ただし、これを行う前に、ナイーブ推定量を考える価値がある。

ナイーブ推定量

カプラン=マイヤー推定量の能力を理解するために、まず、生存関数のナイーブ推定量（naive estimator）を説明する価値がある。

${\displaystyle k\in [n]:=\{1,\dots ,n\}}$ とし、${\displaystyle t>0}$ とする。基本的な議論により、以下の命題が成立することがわかる。

命題1：イベント ${\displaystyle k}$ の打ち切り時間 ${\displaystyle c_{k}}$ が ${\displaystyle t}$ (${\displaystyle c_{k}\geq t}$) を超える場合、${\displaystyle \tau _{k}=t}$ である場合に限り、${\displaystyle {\tilde {\tau }}_{k}=t}$ になる。

${\displaystyle c_{k}\geq t}$ となるような ${\displaystyle k}$ があるとしよう。上記の命題から、

${\displaystyle \operatorname {Prob} (\tau _{k}\geq t)=\operatorname {Prob} ({\tilde {\tau }}_{k}\geq t).}$

が成り立つ。

${\displaystyle X_{k}=\mathbb {I} ({\tilde {\tau }}_{k}\geq t)}$ とし、${\displaystyle k\in C(t):=\{1\leq k\leq n\,:\,c_{k}\geq t\}}$ のものだけ、つまり時刻 ${\displaystyle t}$ 以前に結果が打ち切られなかった事象を考えよう。${\displaystyle m(t)=|C(t)|}$ を ${\displaystyle C(t)}$ の要素の数としよう。なお、集合 ${\displaystyle C(t)}$ は確率的ではないので、${\displaystyle m(t)}$ も確率的ではないことに注意を要する。さらに、${\displaystyle (X_{k})_{k\in C(t)}}$ は、共通パラメータ ${\displaystyle S(t-1)=\operatorname {Prob} (\tau \geq t)}$ を持つ独立同分布のベルヌーイ確率変数の列である。${\displaystyle m(t)>0}$ と仮定すると、

${\displaystyle {\hat {S}}_{\text{naive}}(t-1)={\frac {1}{m(t)}}\sum _{k:c_{k}\geq t}X_{k}={\frac {|\{1\leq k\leq n\,:\,{\tilde {\tau }}_{k}\geq t\}|}{|\{1\leq k\leq n\,:\,c_{k}\geq t\}|}}={\frac {|\{1\leq k\leq n\,:\,{\tilde {\tau }}_{k}\geq t\}|}{m(t)}},}$

を用いて ${\displaystyle S(t-1)}$ を推定することになる。ここで、${\displaystyle {\tilde {\tau }}_{k}\geq t}$ は ${\displaystyle c_{k}\geq t}$ を意味するため、2番目の等式が続く。最後の等式は単に表記法の変更である。

この推定量の質は、${\displaystyle m(t)}$ の大きさによって決まる。これは、${\displaystyle m(t)}$ が小さい場合に問題となる、これは定義上、多くのイベントが打ち切られた場合に起こる。この推定量の特に不快な特性は、おそらくそれが「最良」の推定量ではないことを示唆しており、それは打ち切り時間が ${\displaystyle t}$ より前のすべての観測を無視することである。直感的には、これらの観測はまだ ${\displaystyle S(t)}$ に関する情報を含んでいる。たとえば、${\displaystyle c_{k}<t}$ の多くのイベントで、${\displaystyle {\tilde {\tau }}_{k}<c_{k}}$ も成り立つ場合、イベントが早期に起こることが多いと推測できる。これは、${\displaystyle \operatorname {Prob} (\tau \leq t)}$ が大きいことを意味し、${\displaystyle S(t)=1-\operatorname {Prob} (\tau \leq t)}$ を介して、${\displaystyle S(t)}$ は小さくなければならないことを意味する。ただし、このナイーブ推定法では、この情報は無視される。そこで問題となるのは、すべてのデータをより有効に利用できる推定量が存在するかどうかである。これを実現したのが、カプラン=マイヤー推定量である。なお、打ち切りが行われていない場合には、ナイーブ推定量を改善することはできないので注意を要する。したがって、改善できるかどうかは打ち切りが行われているかどうかに決定的に依存する。

プラグインアプローチ

基本的な計算によって、

${\displaystyle {\begin{aligned}S(t)&=\operatorname {Prob} (\tau >t\mid \tau >t-1)\operatorname {Prob} (\tau >t-1)\\[4pt]&=(1-\operatorname {Prob} (\tau \leq t\mid \tau >t-1))\operatorname {Prob} (\tau >t-1)\\[4pt]&=(1-\operatorname {Prob} (\tau =t\mid \tau \geq t))\operatorname {Prob} (\tau >t-1)\\[4pt]&=q(t)S(t-1)\,,\end{aligned}}}$

となり、ここで最後の等式は ${\displaystyle \tau }$ が整数値であることを利用し、最終行で

${\displaystyle q(t)=1-\operatorname {Prob} (\tau =t\mid \tau \geq t).}$

を導いた。

等式 ${\displaystyle S(t)=q(t)S(t-1)}$ を再帰的に展開すると、

${\displaystyle S(t)=q(t)q(t-1)\cdots q(0).}$

となる。

ここで、${\displaystyle q(0)=1-\operatorname {Prob} (\tau =0\mid \tau >-1)=1-\operatorname {Prob} (\tau =0)}$ となることに注意すること。

カプラン=マイヤー推定量は、各 ${\displaystyle q(s)}$ がデータに基づいて推定され、${\displaystyle S(t)}$ の推定量はこれらの推定量の積として得られる「プラグイン推定量」（plug-in estimator）と見なすことができる。

あとは、${\displaystyle q(s)=1-\operatorname {Prob} (\tau =s\mid \tau \geq s)}$ をどのように推定するかを指定するだけである。命題1により、${\displaystyle c_{k}\geq s}$ となるような任意の ${\displaystyle k\in [n]}$ に対して、${\displaystyle \operatorname {Prob} (\tau =s)=\operatorname {Prob} ({\tilde {\tau }}_{k}=s)}$ と ${\displaystyle \operatorname {Prob} (\tau \geq s)=\operatorname {Prob} ({\tilde {\tau }}_{k}\geq s)}$ がともに成立する。したがって、${\displaystyle c_{k}\geq s}$ となるような任意の ${\displaystyle k\in [n]}$ に対して、

${\displaystyle \operatorname {Prob} (\tau =s|\tau \geq s)=\operatorname {Prob} ({\tilde {\tau }}_{k}=s)/\operatorname {Prob} ({\tilde {\tau }}_{k}\geq s).}$

となる。

上記のナイーブ推定量の構築と同様の推論により、

${\displaystyle {\hat {q}}(s)=1-{\frac {|\{1\leq k\leq n\,:\,c_{k}\geq s,{\tilde {\tau }}_{k}=s\}|}{|\{1\leq k\leq n\,:\,c_{k}\geq s,{\tilde {\tau }}_{k}\geq s\}|}}=1-{\frac {|\{1\leq k\leq n\,:\,{\tilde {\tau }}_{k}=s\}|}{|\{1\leq k\leq n\,:\,{\tilde {\tau }}_{k}\geq s\}|}}}$

という推定量が得られる（「ハザード率」${\displaystyle \operatorname {Prob} (\tau =s|\tau \geq s)}$ の定義において、分子と分母を別々に推定することを考えてみよ）。そして、カプラン=マイヤー推定量は、

${\displaystyle {\hat {S}}(t)=\prod _{s=0}^{t}{\hat {q}}(s).}$

で与えられる。

記事の冒頭で述べた推定量の形式は、さらにいくつかの代数を用いて得られる。そのためには、${\displaystyle {\hat {q}}(s)=1-d(s)/n(s)}$ と記述する。ここで、保険数理の用語を用いて、${\displaystyle d(s)=|\{1\leq k\leq n\,:\,{\tilde {\tau }}_{k}=s\}|}$ は時刻 ${\displaystyle s}$ における既知の死亡者数であり、${\displaystyle n(s)=|\{1\leq k\leq n\,:\,{\tilde {\tau }}_{k}\geq s\}|}$ は時刻 ${\displaystyle s-1}$ において生存している人の数とする。

なお、${\displaystyle d(s)=0}$ であれば、${\displaystyle {\hat {q}}(s)=1}$ であることに注意を要する。このことは、${\displaystyle {\hat {S}}(t)}$ を定義する積から、${\displaystyle d(s)=0}$ の項をすべて除外できることを意味する。そして、${\displaystyle d(s)>0}$、${\displaystyle d_{i}=d(t_{i})}$、${\displaystyle n_{i}=n(t_{i})}$ の時、${\displaystyle 0\leq t_{1}<t_{2}<\dots <t_{m}}$ を時間 ${\displaystyle s}$ とすると、冒頭に述べたカプラン=マイヤー推定量の形、

${\displaystyle {\hat {S}}(t)=\prod _{i:t_{i}\leq t}\left(1-{\frac {d_{i}}{n_{i}}}\right).}$

になる。

この推定量は、ナイーブ推定量とは対照的に、利用可能な情報をより効果的に利用していることがわかる。前述の特殊なケースでは、多くの初期イベントが記録されている場合、推定量は1未満の値を持つ多くの項を乗算するため、その結果、生存確率が大きくならないことを考慮に入れよ。

最尤推定量としての導出

カプラン=マイヤー推定量は、ハザード関数の最尤推定から導出できる^[8]。より具体的には、イベントの数を ${\displaystyle d_{i}}$、時刻 ${\displaystyle t_{i}}$ でのリスクのある個人の総数を ${\displaystyle n_{i}}$ とすると、離散ハザード率 ${\displaystyle h_{i}}$ は、時刻 ${\displaystyle t_{i}}$ でイベントが発生した個人の確率として定義できる。この場合、生存率は次のように定義でき、

${\displaystyle S(t)=\prod \limits _{i:\ t_{i}\leq t}(1-h_{i})}$

時刻 ${\displaystyle t_{i}}$ までのハザード関数に対する尤度関数は、

${\displaystyle {\mathcal {L}}(h_{j:j\leq i}\mid d_{j:j\leq i},n_{j:j\leq i})=\prod _{j=1}^{i}h_{j}^{d_{j}}(1-h_{j})^{n_{j}-d_{j}}}$

となり、したがって対数尤度は次のようになる。

${\displaystyle \log({\mathcal {L}})=\sum _{j=1}^{i}\left(d_{j}\log(h_{j})+(n_{j}-d_{j})\log(1-h_{j})\right)}$

${\displaystyle h_{i}}$ に対する対数尤度の最大値は、

${\displaystyle {\frac {\partial \log({\mathcal {L}})}{\partial h_{i}}}={\frac {d_{i}}{{\widehat {h}}_{i}}}-{\frac {n_{i}-d_{i}}{1-{\widehat {h}}_{i}}}=0\Rightarrow {\widehat {h}}_{i}={\frac {d_{i}}{n_{i}}}}$

と求められる。ここでハット記号（${\displaystyle {\hat {}}}$）は最尤推定を表すのに用いられている。この結果から、次のように書くことができる。

${\displaystyle {\widehat {S}}(t)=\prod \limits _{i:\ t_{i}\leq t}\left(1-{\widehat {h}}_{i}\right)=\prod \limits _{i:\ t_{i}\leq t}\left(1-{\frac {d_{i}}{n_{i}}}\right)}$

利点と限界

カプラン=マイヤー推定量は、生存分析で最も頻繁に使用される手法の1つである。この推定量は、回復率、死亡の確率、および治療の有効性を検討するのに有用である。その共変量（英語版）で調整された生存率を推定する能力には限界がある。共変量で調整された生存率を推定するには、パラメトリック生存モデルおよびCox比例ハザード検定が有用であろう。

統計学的考察

カプラン=マイヤー推定量は統計量であり、その分散を近似するためにいくつかの推定量が使用される。最も一般的な推定量の1つはGreenwoodの式である^[9]。

${\displaystyle {\widehat {\operatorname {Var} }}\left({\widehat {S}}(t)\right)={\widehat {S}}(t)^{2}\sum _{i:\ t_{i}\leq t}{\frac {d_{i}}{n_{i}(n_{i}-d_{i})}},}$

ここで、${\displaystyle d_{i}}$ は症例数、${\displaystyle n_{i}}$ は観測の総数で、${\displaystyle t_{i}<t}$ である。

上の式の数学的導出の「スケッチ」については、「表示」をクリックして表示される。

Greenwoodの式は、${\displaystyle n_{i}}$ 件のうち ${\displaystyle d_{i}}$ 件の故障が発生する確率が、故障確率 ${\displaystyle h_{i}}$ の二項分布に従うことに着目して導き出される^[10]。最尤ハザード率 ${\displaystyle {\widehat {h}}_{i}=d_{i}/n_{i}}$ に対する結果として、${\displaystyle E\left({\widehat {h}}_{i}\right)=h_{i}}$、および ${\displaystyle \operatorname {Var} \left({\widehat {h}}_{i}\right)=h_{i}(1-h_{i})/n_{i}}$となる。確率の乗算を回避するために、${\displaystyle {\widehat {S}}(t)}$ の対数の分散を計算し、デルタ法（英語版）を使用して元の分散に変換する。

${\displaystyle {\begin{aligned}\operatorname {Var} \left(\log {\widehat {S}}(t)\right)&\sim {\frac {1}{{{\widehat {S}}(t)}^{2}}}\operatorname {Var} \left({\widehat {S}}(t)\right)\Rightarrow \\\operatorname {Var} \left({\widehat {S}}(t)\right)&\sim {{{\widehat {S}}(t)}^{2}}\operatorname {Var} \left(\log {\widehat {S}}(t)\right)\end{aligned}}}$

マルチンゲール中心極限定理（英語版）を使用して、次の式の和の分散が、分散の和に等しいことが示される^[10]。:

${\displaystyle \log {\widehat {S}}(t)=\sum \limits _{i:\ t_{i}\leq t}\log \left(1-{\widehat {h}}_{i}\right)}$

結果として、次のように記述できる。

${\displaystyle {\begin{aligned}\operatorname {Var} ({\widehat {S}}(t))&\sim {{{\widehat {S}}(t)}^{2}}\operatorname {Var} \left(\sum _{i:\ t_{i}\leq t}\log \left(1-{\widehat {h}}_{i}\right)\right)\\&\sim {{{\widehat {S}}(t)}^{2}}\sum \limits _{i:\ t_{i}\leq t}\operatorname {Var} \left(\log \left(1-{\widehat {h}}_{i}\right)\right)\end{aligned}}}$

必要に応じて、デルタ法をもう一度使用して、

${\displaystyle {\begin{aligned}\operatorname {Var} ({\widehat {S}}(t))&\sim {{{\widehat {S}}(t)}^{2}}\sum _{i:\ t_{i}\leq t}\left({\frac {\partial \log \left(1-{\widehat {h}}_{i}\right)}{\partial {\widehat {h}}_{i}}}\right)^{2}\operatorname {Var} \left({\widehat {h}}_{i}\right)\\&={{{\widehat {S}}(t)}^{2}}\sum _{i:\ t_{i}\leq t}\left({\frac {1}{1-{\widehat {h}}_{i}}}\right)^{2}{\frac {{\widehat {h}}_{i}\left(1-{\widehat {h}}_{i}\right)}{n_{i}}}\\&={{{\widehat {S}}(t)}^{2}}\sum _{i:\ t_{i}\leq t}{\frac {{\widehat {h}}_{i}}{n_{i}\left(1-{\widehat {h}}_{i}\right)}}\\&={{{\widehat {S}}(t)}^{2}}\sum _{i:\ t_{i}\leq t}{\frac {d_{i}}{n_{i}(n_{i}-d_{i})}}\end{aligned}}}$

となる。

---

場合によっては、異なるカプラン=マイヤー曲線を比較したいことがある。これは、ログランク検定、およびCox比例ハザード検定によって行うことができる。

この推定量で使用できる他の統計量は、Hall-Wellnerバンド（Hall-Wellner band）^[11]および等精度バンド（equal-precision band）である^[12]。

ソフトウェア

• Mathematica: 組み込み関数 SurvivalModelFit で生存モデルを作成^[13]。
• SAS（英語版） - カプラン=マイヤー推定量は proc lifetest プロシージャで実装されている^[14]。
• R - カプラン=マイヤー推定量は、 survival パッケージの一部として利用可能である^[15][16][17]。
• Stata - コマンド sts は、カプラン=マイヤー推定量を返す^[18][19]。
• Python - lifelines パッケージにはカプラン=マイヤー推定量が含まれている^[20]。
• MATLAB - ecdf 関数に 'function','survivor' 引数を指定すると、カプラン=マイヤー推定量を計算したり、プロットしたりすることができる^[21]。
• StatsDirect - カプラン=マイヤー推定量は Survival Analysis メニューに実装されている^[22]。
• SPSS - カプラン=マイヤー推定量は、 Analyze > Survival > Kaplan-Meier... メニューに実装されている^[23]。
• Julia - Survival.jl パッケージには、カプラン=マイヤー推定量が含まれている^[24]。

参照項目

• 生存分析
• 超過頻度（英語版）
• 半数致死量
• ネルソン=アーラン推定量（英語版）