生存分析

生存分析では、以下の用語が一般的に使用されます。^{[ 3 ]}

• イベント：死亡、疾病発生、疾病再発、回復、またはその他の注目すべき出来事
• 時間: 観察期間の開始 (手術や治療の開始など) から (i) イベント発生、(ii) 研究終了、(iii) 連絡の喪失または研究からの離脱までの時間。
• 検閲／検閲された観察：検閲とは、個々の生存時間に関する情報はあるものの、正確な生存時間がわからない場合に発生します。検閲された対象は、検閲時点以降、その対象について何も観察または把握されないという意味で検閲されます。検閲された対象は、観察終了時点以​​降に何らかの事象が発生する場合もあれば、発生しない場合もあります。
• 生存関数S(t): 対象者が時間tよりも長く生存する確率。

この例では、R の「survival」パッケージの急性骨髄性白血病生存データセット「aml」を使用します。データセットは Miller (1997) ^{[ 4 ]}によるもので、標準的な化学療法コースをさらにサイクル延長（「維持」）すべきかどうかが問題です。

生存時間順に並べられたAMLデータセットが枠内に表示されています。

• 時間は変数「time」で示され、これは生存時間または打ち切り時間である。
• イベント（AMLがんの再発）は変数「status」で示されます。0  ＝イベントなし（打ち切り）、1  ＝イベント（再発）
• 治療群：変数「x」は維持化学療法が実施されたかどうかを示します

最後の観察（11）は161週で打ち切りとなりました。打ち切りは、患者にイベント（AMLがんの再発なし）が発生しなかったことを示しています。別の被験者である観察3は13週で打ち切りとなりました（ステータス=0で示されています）。この被験者は13週間しか研究に参加しておらず、その13週間の間にAMLがんは再発しませんでした。この患者は研究の終盤近くに登録されたため、13週間しか観察できなかった可能性があります。また、患者は研究の初期に登録されたものの、追跡不能になったか、研究から脱落した可能性もあります。表には、他の被験者が16週、28週、45週で打ち切りになったことが示されています（ ステータス=0の観察17、6、9）。残りの被験者はすべて、研究中にイベント（AMLがんの再発）を経験しました。ここで注目すべき点は、維持療法を受けている患者では、維持療法を受けていない患者よりも再発が遅れて起こるかどうかである。

ログランク検定：AMLデータにおける生存率の差を検定する

ログランク検定は、 2 つ以上のグループの生存時間を比較します。^{[ 6 ]}この例では、AML データにおける維持療法群と非維持療法群の生存時間の差についてログランク検定を使用します。グラフは、データ内の変数「x」で示される治療群ごとに分けられた、AML データの KM プロットを示しています。

ログランク検定の帰無仮説は、グループの生存率が同じであるというものです。各時点における各グループの生存者の期待数は、各イベント時点におけるグループのリスクのある被験者数に合わせて調整されます。ログランク検定は、各グループで観察されたイベント数が期待数と有意に異なるかどうかを判定します。^{[ 7 ]}正式な検定はカイ二乗統計量に基づいています。ログランク統計量が大きい場合、グループ間の生存時間に差があるという証拠となります。ログランク統計量は、自由度1のカイ二乗分布に近似的に従い、p値はカイ二乗検定を使用して計算されます。

例のデータでは、生存率の差に関するログランク検定のp値はp=0.0653となり、α水準を0.05と仮定すると、治療群間で生存率に有意差はないことを示しています。23人の被験者というサンプルサイズは控えめなので、治療群間の差を検出する力はほとんどありません。カイ二乗検定は漸近近似に基づいているため、サンプルサイズが小さい場合はp値を慎重に扱う必要があります。^{[ 8 ]}

カプラン・マイヤー曲線とログランク検定は、予測変数がカテゴリ変数（例：薬剤 vs. プラセボ）である場合、または カテゴリ変数として扱える少数の値（例：薬剤投与量 0、20、50、100 mg/日）をとる場合に最も有用です。^{[ 9 ]}ログランク検定と KM 曲線は、遺伝子発現、白血球数、年齢などの定量的予測変数には容易に適用できません。定量的予測変数の場合、代替方法としてCox 比例ハザード回帰分析があります。Cox PH モデルは、{0,1} 指標またはダミー変数としてエンコードされたカテゴリ予測変数にも適用できます。ログランク検定は Cox PH 分析の特殊なケースであり、Cox PH ソフトウェアを使用して実行できます。

例：メラノーマに対するコックス比例ハザード回帰分析

この例では、ダルガード第14章のメラノーマデータセットを使用します。 ^{[ 10 ]}

データはRパッケージISwRに格納されています。Rを用いたCox比例ハザード回帰分析の 結果は、枠内に示したとおりです。

コックス回帰分析の結果は、以下のように解釈される。

• 性別は数値ベクトル（1：女性、2：男性）として符号化されます。Cox モデルのRサマリーは、第1グループ（男性対女性）に対する第2グループのハザード比（HR）を示します。
• coef = 0.662 は、男性と女性のハザード比の推定対数です。
• exp(coef) = 1.94 = exp(0.662) - ハザード比 (coef = 0.662) の対数を exp(coef) を使用してハザード比に変換します。Cox モデルの要約では、第 1 グループに対する第 2 グループ (男性対女性) のハザード比が示されています。推定ハザード比 1.94 は、このデータでは男性の方が女性よりも死亡リスクが高い (生存率が低い) ことを示しています。
• se(coef) = 0.265 は、対数ハザード比の標準誤差です。
• z = 2.5 = coef/se(coef) = 0.662/0.265。coefを標準誤差で割るとzスコアが得られます。
• p=0.013。性別に関するz=2.5に対応するp値はp=0.013であり、性別によって生存率に有意な差があることを示している。

要約出力には、ハザード比の95%信頼区間の上限と下限も表示されます。95%信頼区間の下限は1.15、上限は3.26です。

最後に、出力結果には、モデル全体の有意性に関する3つの代替検定のp値が表示されます。

• 尤度比検定 = 6.15 (自由度 1)、p = 0.0131
• ウォルド検定 = 6.24 (自由度 1)、p = 0.0125
• スコア（ログランク）検定 = 6.47、自由度1、p=0.0110

これら3つの検定は漸近的に同等です。Nが十分に大きい場合、同様の結果が得られます。Nが小さい場合は、多少異なる可能性があります。最後の行「スコア（ログランク）検定」はログランク検定の結果で、p=0.011です。ログランク検定はCox比例ハザード回帰の特殊なケースであるため、ログランク検定と同じ結果となります。尤度比検定はサンプルサイズが小さい場合に優れた挙動を示すため、一般的にこちらが好まれます。

メラノーマデータにおける共変量を用いたCoxモデル

Coxモデルは、追加の共変量を含めることを可能にすることで、ログランク検定を拡張します。^{[ 11 ]}この例では、予測変数に腫瘍の厚さ（変数名 = "thick")という連続共変量が含まれるメラノーマデータセットを使用します。

ヒストグラムでは、厚さの値は正に歪んでおり、ガウス分布のような対称的な確率分布を示しません。Coxモデルを含む回帰モデルは、一般的に正規分布する変数を使用した場合により信頼性の高い結果が得られます。この例では、対数変換を使用する場合があります。腫瘍の厚さの対数は正規分布に近いように見えるため、Coxモデルでは対数の厚さを使用します。Rを使用したCox PH分析の結果は以下のとおりです。

library ( ISwR ) # CRAN v2.0-11 library ( survival ) # CRAN v3.8-3# 共変量 log 腫瘍厚を使用してメラノーマデータセットの Cox PH 分析を実行しますf <- coxph ( Surv ( days , status == 1 ) ~ sex + log ( thick ), data = melanom )# Cox PH 出力の要約を表示します( f ) # 呼び出し: # coxph(formula = Surv(days, status == 1) ~ sex + log(thick), data = melanom) # # n= 205、イベント数= 57 # # coef exp(coef) se(coef) z Pr(>|z|) # sex 0.4580 1.5809 0.2687 1.705 0.0883 . # log(thick) 0.7809 2.1834 0.1573 4.963 6.94e-07 *** # --- # 有意水準コード: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 # # exp(coef) exp(-coef) lower .95 upper .95 # sex 1.581 0.6326 0.9337 2.677 # log(thick) 2.183 0.4580 1.6040 2.972 # # Concordance= 0.749 (se = 0.033 ) # Likelihood ratio test= 33.45 on 2 df, p=5e-08 # Wald test = 31 on 2 df, p=2e-07 # Score (logrank) test = 32.52 on 2 df, p=9e-08

3 つの全体検定 (尤度、ワルド検定、スコア検定) の p 値は有意であり、モデルが有意であることを示しています。${\displaystyle \log({\text{thick}})}$はハザード比は6.94 × 10 ⁻⁷${\displaystyle HR=\exp({\text{coef}})=2.183}$これは、腫瘍の厚さと死亡リスクの増加との間に強い相関関係があることを示している。

対照的に、p値は$\displaystyle {\text{セックス}}}$今は${\displaystyle p=0.088}$ハザード比${\displaystyle HR=\exp({\text{coef}})=1.581}$95%信頼区間は0.9337～2.677である。HRの信頼区間には1が含まれているため、これらの結果は、腫瘍の厚さを調整した後、性別がHRの差に与える影響が小さいことを示している。グラフの検討${\displaystyle log({\text{thick}})}$性別とt検定による${\displaystyle log({\text{thick}})}$性別による分析結果はいずれも、初診時の腫瘍の厚さにおいて、男性と女性の間に有意な差があることを示している。

Coxモデルは、ハザードが比例関係にあることを前提としています。比例ハザードの仮定は、survival RパッケージのR 関数を使用して検証できますcox.zph()。p値が0.05未満の場合、ハザードは比例関係にないことを示します。

cox.zph ( f ) # カイ二乗 自由度 p # 性別 1.33 1 0.248 # log(thick) 6.23 1 0.013 # グローバル 6.70 2 0.035

メラノーマのデータについては、${\displaystyle p=0.248}$性行為の場合。したがって、性行為の危険性が比例するという帰無仮説を棄却することはできません。

Cox比例ハザード回帰モデルは線形モデルです。線形回帰やロジスティック回帰と類似しています。具体的には、これらの手法は、単一の直線、曲線、平面、または曲面でグループ（生存、死亡）を分離したり、定量的な応答（生存時間）を推定したりするのに十分であると仮定しています。

場合によっては、代替の分割方法の方がより正確な分類や定量的推定値が得られることがあります。代替方法の1つとして、ツリー構造の生存モデル^{[ 12 ] [ 13 ] [ 14 ]}（生存ランダムフォレスト^{[ 15 ]}を含む）があります。ツリー構造の生存モデルは、Coxモデルよりも正確な予測値を与える可能性があります。特定のデータセットに対して両方のタイプのモデルを検討することは、妥当な戦略です。

生存ツリー分析の例

この生存ツリー分析の例では、R パッケージ「rpart」を使用しています。^{[ 16 ]}この例は、 rpart のデータセット stagec に含まれる 146 人のステージ C 前立腺がん患者に基づいています。rpart と stagec の例については、Atkinson と Therneau (1997) ^{[ 17 ]}で説明されており、rpart パッケージのビネットとしても配布されています。^{[ 16 ]}

各段階における変数は以下のとおりです。

• pgtime：病勢進行までの時間、または病勢進行が認められなかった最後の追跡調査時点
• pgstat：最終フォローアップ時のステータス（1＝進行、0＝打ち切り）
• 年齢：診断時の年齢
• eet：早期内分泌療法（1＝いいえ、0＝はい）
• 倍数性：二倍体／四倍体／異数体のDNAパターン
• g2： G2期にある細胞の割合
• グレード：腫瘍グレード（1～4）
• グリーソン：グリーソングレード（3～10）

解析によって得られた生存ツリーを図に示す。

ツリーの各枝は、変数の値による分割を示しています。たとえば、ツリーの根は、グレードが2.5未満の被験者とグレードが2.5以上の被験者を分割し​​ます。末端ノードは、ノード内の被験者数、イベントが発生した被験者数、および根と比較した相対イベント率を示します。一番左のノードでは、値1/33は、ノード内の33人の被験者のうち1人がイベントが発生したこと、および相対イベント率が0.122であることを示しています。一番右下のノードでは、値11/15は、ノード内の15人の被験者のうち11人がイベントが発生したこと、および相対イベント率が2.7であることを示しています。

深層表現学習の最近の進歩は、生存推定にも応用されています。DeepSurv ^{[ 19 ]}モデルは、CoxPH モデルの対数線形パラメータ化を多層パーセプトロンに置き換えることを提案しています。Deep Survival Machines ^{[ 20 ]}や Deep Cox Mixtures ^{[ 21 ]}などのさらなる拡張では、潜在変数混合モデルを使用して、イベント発生までの時間分布をパラメトリック分布またはセミパラメトリック分布の混合としてモデル化し、同時に入力共変量の表現を学習します。深層学習アプローチは、画像や臨床時系列などの複雑な入力データモダリティにおいて特に優れた性能を示しています。

主な関心対象は生存関数であり、慣例的にSと表記され、次のように定義される。$${\displaystyle S(t)=\Pr(T>t)}$$ ここで、 tはある時間、Tは死亡時間、またはより一般的には任意のイベントを表す確率変数であり、「Pr」は確率を表します。つまり、生存関数は、特定の時間tよりも大きい時間でイベントを観測する確率です。^{[ 7 ]} 生存関数は、生物学的生存の問題では生存関数または生存率関数とも呼ばれ、機械的生存の問題では信頼性関数とも呼ばれます。^{[ 22 ]}後者の場合、信頼性関数はR ( t ) で表されます。

通常はS (0)=1と仮定されるが、 即時死亡または失敗の可能性がある場合は1未満になる可能性がある。

生存関数は非増加でなければならない。u ≥ tの場合、S ( u ) ≤ S ( t )である。この性質は、 T > u がT > tを意味することから直接導かれる。これは、より若い年齢をすべて達成した場合にのみ、より後の年齢まで生存できるという考え方を反映している。この性質を考慮すると、寿命分布関数とイベント密度 (以下のFとf ) は適切に定義される。^{[ 2 ]}

生存関数は通常、年齢が無限に増加するにつれてゼロに近づくと仮定されます（つまり、t → ∞ のときS ( t ) → 0 ）。ただし、永遠の命が可能であれば、その極限値はゼロより大きくなる可能性があります。例えば、安定炭素同位体と不安定炭素同位体の混合物に生存分析を適用することができます。不安定同位体は遅かれ早かれ崩壊しますが、安定同位体は永久に存在し続けます。

関連する量は、生存関数を用いて定義される。

寿命分布関数は、慣例的にFと表記され、生存関数の補関数として定義される。

$${\displaystyle F(t)=\Pr(T\leq t)=1-S(t).}$$Fが微分 可能であれば、寿命分布の密度関数である導関数は慣例的にfと表記される。

$${\displaystyle f(t)=F'(t)={\frac {d}{dt}}F(t).}$$ 関数fは、イベント密度と呼ばれることもあり、単位時間あたりの死亡または故障イベントの発生率を表します。

生存関数は、確率分布と確率密度関数を用いて表現することができる。

$${\displaystyle S(t)=\Pr(T>t)=\int _{t}^{\infty }f(u)\,du=1-F(t).}$$ 同様に、生存イベント密度関数は次のように定義できます。

$${\displaystyle s(t)=S'(t)={\frac {d}{dt}}S(t)={\frac {d}{dt}}\int _{t}^{\infty }f(u)\,du={\frac {d}{dt}}[1-F(t)]=-f(t).}$$ 統計物理学などの他の分野では、生存事象の密度関数は、初回通過時間密度として知られています。

ハザード関数${\displaystyle h}$は、時刻におけるイベント発生率として定義される。${\displaystyle t,}$時点での生存を条件とする${\displaystyle t.}$

ハザード関数の同義語には、さまざまな分野において、ハザード率、強度関数、^{[ 23 ]}死亡力（人口統計学およびアクチュアリー科学、以下で示される）などがあります。${\displaystyle \mu }$破壊力、または破壊率（工学用語では、${\displaystyle \lambda }$例えば、アクチュアリー科学では、${\displaystyle \mu (x)}$は、18歳以上の人の死亡率を表します。${\displaystyle x}$信頼性工学では${\displaystyle \lambda (t)}$一定時間稼働後の部品の故障率を示す${\displaystyle t}$。

ハザード関数は、対象者が時間tまで生存しているという条件の下で、次の短い時間間隔で事象を経験する確率を、この時間間隔の長さで割った値を表します。${\displaystyle t}$[ ^{23 ]正式には、これは} 次のように記述できます。

$${\displaystyle h(t)=\lim _{dt\rightarrow 0}{\frac {\Pr(t<T<t+dt|T>t)}{dt}}=\lim _{dt\rightarrow 0}{\frac {\Pr(t<T<t+dt)}{dt\cdot S(t)}}={\frac {f(t)}{S(t)}}=-{\frac {S'(t)}{S(t)}},}$$

ベイズの定理、そして識別${\displaystyle \Pr(T>t)}$最初の等式では生存関数として、2番目の等式では寿命分布の密度関数の定義が使用されています。

任意の機能${\displaystyle h}$ハザード関数であるのは、以下の性質を満たす場合に限る。

1. ${\displaystyle \forall x\geq 0\left(h(x)\geq 0\right)}$、
2. $\displaystyle \int _{0}^{\infty }h(x)dx=\infty }$。

実際、ハザード率は、寿命分布の他の表現方法よりも、故障の根本的なメカニズムについてより多くの情報を提供することが多い。

ハザード関数は非負でなければならない。$\displaystyle \lambda (t)\geq 0$、そしてその積分範囲は${\displaystyle [0,\infty ]}$は無限でなければならないが、それ以外に制約はない。増加関数、減少関数、非単調関数、不連続関数の場合がある。例として、バスタブ曲線ハザード関数がある。これは、の値が小さい場合に大きくなる。${\displaystyle t}$減少して最小値に達した後、再び増加する。これは、機械システムが動作直後に故障するか、システムの経年劣化に伴ってずっと後になって故障するかという特性をモデル化することができる。

ハザード関数は、慣例的に累積ハザード関数で表されることもあります。${\displaystyle \Lambda }$または${\displaystyle H}$:

$${\displaystyle \,\Lambda (t)=-\log S(t)}$$ 符号を入れ替えて指数化すると

$${\displaystyle \,S(t)=\exp(-\Lambda (t))}$$ または（連鎖律を用いて）微分する

$$\displaystyle {\frac {d}{dt}}\Lambda (t)=-{\frac {S'(t)}{S(t)}}=\lambda (t).}$$ 「累積ハザード関数」という名称は、

$${\displaystyle \Lambda (t)=\int _{0}^{t}\lambda (u)\,du}$$ これは、時間の経過とともに危険が「蓄積」していくことを意味します。

定義から$\displaystyle \Lambda (t)$t が無限大に近づくにつれて、無限に増加することがわかります(${\displaystyle S(t)}$（ゼロに近づく）。これは、${\displaystyle \lambda (t)}$定義上、累積ハザードは発散しなければならないため、あまり急速に減少してはならない。例えば、${\displaystyle \exp(-t)}$これは、その積分が1に収束するため、いかなる生存分布のハザード関数でもありません。

生存関数${\displaystyle S(t)}$累積ハザード関数$\displaystyle \Lambda (t)$密度${\displaystyle f(t)}$ハザード関数${\displaystyle \lambda (t)}$、そして寿命分布関数${\displaystyle F(t)}$は $${\displaystyle S(t)=\exp[-\Lambda (t)]={\frac {f(t)}{\lambda (t)}}=1-F(t),\quad t>0.}$$

ある時点における将来の寿命$t0$は、年齢まで生存した場合の、死亡までの残り時間です。$t0$したがって、${\displaystyle T-t_{0}}$現在の表記では、期待される将来寿命は将来寿命の期待値です。年齢で、またはそれ以前に死亡する確率${\displaystyle t_{0}+t}$年齢まで生存した場合$t0$、それはただ

$${\displaystyle P(T\leq t_{0}+t\mid T>t_{0})={\frac {P(t_{0}<T\leq t_{0}+t)}{P(T>t_{0})}}={\frac {F(t_{0}+t)-F(t_{0})}{S(t_{0})}}.}$$ したがって、将来の寿命の確率密度は

$${\displaystyle {\frac {d}{dt}}{\frac {F(t_{0}+t)-F(t_{0})}{S(t_{0})}}={\frac {f(t_{0}+t)}{S(t_{0})}}}$$ そして、予想される将来の寿命は

$${\displaystyle {\frac {1}{S(t_{0})}}\int _{0}^{\infty }t\,f(t_{0}+t)\,dt={\frac {1}{S(t_{0})}}\int _{t_{0}}^{\infty }S(t)\,dt,}$$ ここで、2番目の式は部分積分を用いて得られる。

のために${\displaystyle t_{0}=0}$つまり、出生時においては、これは平均寿命に相当します。

信頼性問題において、期待寿命は平均故障時間と呼ばれ、期待される将来の寿命は平均残存寿命と呼ばれます。

個体が年齢t以降まで生存する確率はS ( t ) であるため、定義により、すべての個体に対して同じ生存関数を仮定すると、初期集団nの新生児のうち年齢tで生存する個体の期待値はn × S ( t ) となります。したがって、生存者の期待割合はS ( t ) です。異なる個体の生存が独立している場合、年齢tでの生存者数はパラメータnとS ( t ) を持つ二項分布に従い、生存者の割合の分散はS ( t ) × (1- S ( t ))/ nとなります。

生存者のうち特定の割合が残っている年齢は、方程式S ( t ) = q をtについて解くことで求められます。ここで、q は対象となる分位数です。通常、関心があるのは中央値寿命で、この場合q = 1/2 となります。また、 q = 0.90 やq = 0.99などの他の分位数も対象となります。