比例ハザード モデルは、統計学における生存モデルの一種です。生存モデルは、あるイベントが発生するまでの経過時間を、その時間量に関連付けられる可能性がある1 つ以上の共変量に関連付けます。比例ハザード モデルでは、共変量の単位増加の固有の効果は、ハザード率に対して乗法的です。時刻におけるハザード率は、時刻までにイベントがまだ発生していない場合に、イベントが から の間に発生する短い時間 d tあたりの確率です。たとえば、薬を服用すると、脳卒中が発生するハザード率が半分になる場合があり、製造された部品の材料を変更すると、故障のハザード率が 2 倍になる場合があります。加速故障時間モデルなどの他の種類の生存モデルでは、比例ハザードは示されません。加速故障時間モデルは、イベントの生物学的または機械的なライフ ヒストリーが加速 (または減速) される状況を説明します。
背景
生存モデルは、2 つの部分から構成されると考えられます。1 つは、多くの場合 と表記される基礎となるベースラインハザード関数で、単位時間あたりのイベントのリスクが共変量のベースラインレベルで時間の経過とともにどのように変化するかを説明します。もう 1 つは、説明的な共変量に応じてハザードがどのように変化するかを説明する効果パラメーターです。典型的な医学的な例では、変動性を低減したり交絡を制御したりするために、治療割り当てなどの共変量や、研究開始時の年齢、性別、研究開始時の他の疾患の存在などの患者特性が含まれます。
比例ハザード条件[1]は、共変量がハザードと乗法的に関連していることを示しています。たとえば、最も単純な定常係数のケースでは、薬物による治療により、ある時点で被験者のハザードが半分になる場合がありますが、ベースラインのハザードは変化する可能性があります。ただし、これによって被験者の寿命が 2 倍になるわけではないことに注意してください。共変量が寿命に及ぼす正確な影響は、 のタイプによって異なります。共変量はバイナリ予測子に限定されません。連続共変量 の場合、通常、ハザードは指数関数的に応答すると想定されます。 の各単位が増加すると、ハザードが比例的にスケーリングされます。
コックスモデル
導入
デイビッド・コックス卿は、比例ハザード仮定が成り立つ(または、成り立つと仮定される)場合、完全なハザード関数を考慮せずに、以下に示す効果パラメータを推定できることに気づきました。生存データへのこのアプローチは、コックス比例ハザードモデルの適用と呼ばれ、[2]コックスモデルまたは比例ハザードモデルと略されることもあります。[3]しかし、コックスは、比例ハザード仮定の生物学的解釈は非常に難しい場合があることにも言及しました。[4] [5]
X i = ( X i 1 , … , X ip )を被験者iのp個の共変量の実現値とします。Cox 比例ハザード モデルのハザード関数は、次の形式になります。 この式は、共変量ベクトル (説明変数) X iを持つ被験者iの 時刻tにおけるハザード関数を示します。被験者間では、ベースライン ハザードは同一である ( iに依存しない) ことに注意してください。被験者のハザード間の唯一の違いは、ベースライン スケーリング係数 から生じます。
なぜ「比例」と呼ばれるのか
まず、共変量が 1 つだけ、したがって係数が 1 つだけあると仮定します。モデルは次のようになります。
1増加した場合の効果を考えてみましょう。
共変量を 1 増やすと、元のハザードが定数倍されることがわかります。少し整理すると、次のようになります。
右辺は時間の経過に伴って一定です(を含む項はありません)。この関係 は比例関係と呼ばれます。
より一般的には、それぞれ共変量とを持つ2 つの被験者iとj を考えます。それらのハザード比を考えます。
唯一の時間依存因子である が相殺されたため、右辺は時間に依存しません。したがって、2 つの被験者の危険の比率は定数、つまり危険は比例します。
切片項の不在
多くの場合、回帰モデルでは切片項(定数項またはバイアス項とも呼ばれる)が使用されます。Cox モデルには切片項がありません。これは、ベースライン ハザード が切片項の代わりとなるためです。では、切片項 を含めた場合に何が起こるかを見てみましょう。切片項 は と表され、新しいベースライン ハザード として 再定義されています。したがって、ベースライン ハザードには、被験者の共変量に依存しないハザードのすべての部分が組み込まれており、切片項(定義によりすべての被験者で一定)も含まれます。
ユニークな時間の可能性
以下に示すCox の部分尤度は、ベースライン ハザード関数の Breslow の推定値を使用して、それを全尤度に代入し、結果が 2 つの要因の積であることを確認することによって得られます。最初の要因は、以下に示す部分尤度であり、ベースライン ハザードが「相殺」されています。これは、発生時刻のセットと被験者の共変量が与えられた場合に、被験者が実際に発生した順序でイベントを経験する確率です。2 番目の要因は回帰係数を持たず、打ち切りパターンを通じてのみデータに依存します。したがって、比例ハザード モデルによって推定される共変量の影響は、ハザード比として報告できます。
部分尤度、つまりイベントの順序の確率を計算するには、イベントがすでに発生しているM 個のサンプルを発生時刻の昇順でインデックス付けします ( Y 1 < Y 2 < ... < Y M )。イベントが発生していない他のすべての被験者の共変量には、インデックスM +1、...、Nが付けられます。部分尤度は、発生したイベントごとに 1 つの因子に因数分解できます。i番目の因子は、時間Y iより前にイベントが発生していないすべての被験者 ( i、i +1、...、N ) のうち、時間Y iに実際に発生したイベントが被験者iのイベントである確率です。 ここで、θ j = exp( X j ⋅ β ) であり、合計は時間Y iより前にイベントが発生していない被験者jの集合(被験者i自身を含む) に対して行われます。明らかに、0 < L i (β) ≤ 1 です。
被験者を統計的に互いに独立しているものとして扱うと、 イベントの順序の 部分尤度[6]は、 イベントが発生した被験者がC i = 1 で示され、その他すべての被験者がC i = 0 で示される。対応する対数部分尤度は、 上で紹介したインデックスをより一般的な方法で使用して と書き表した である。重要なのは、時間の経過に伴うハザード関数を指定しなくても、共変量の影響を推定できることである。部分尤度をβにわたって最大化することで、モデルパラメータの最大部分尤度推定値を生成することができる。
部分スコア関数は
そして部分対数尤度の ヘッセ行列は
このスコア関数とヘッセ行列を使用すると、ニュートン・ラプソンアルゴリズムを使用して部分尤度を最大化できます。 βの推定値で評価されたヘッセ行列の逆行列は、推定値の近似分散共分散行列として使用でき、回帰係数の 近似標準誤差を生成するために使用できます。
同点の場合の可能性
時間データに同点がある場合の対処法として、いくつかのアプローチが提案されています。Breslow法は、同点がある場合でも上記の手順をそのまま使用するアプローチです。より良い結果が得られると考えられる別のアプローチは、Efron 法です。[7] t j を一意の時間、H j をY i = t jかつC i = 1となるインデックスiのセット 、m j = | H j | とします。Efron 法は、次の部分尤度を最大化します。
対応する対数部分尤度は スコア関数はで あり 、ヘッセ行列は
H jが空の場合(時刻t jのすべての観測が打ち切られている場合)、これらの式の加数はゼロとして扱われることに注意してください。
例
以下は Cox モデルの実践例です。
単一のバイナリ共変量
関心のあるエンドポイントが、手術後 5 年間の観察期間中の患者の生存であるとします。患者は 5 年の期間内に死亡する可能性があり、その場合は死亡した日時を記録します。または、患者が 5 年以上生存する可能性があり、その場合は 5 年以上生存したことのみを記録します。手術は 2 つの病院AまたはBのいずれかで実施され、病院の場所が 5 年生存率と関連しているかどうかを知りたいと考えています。具体的には、病院 A で実施された手術と病院 B で実施された手術のハザードの相対的な増加 (または減少) を知りたいと考えています。提供される (偽の) データでは、各行が患者を表します。T は、患者が死亡する前に観察された期間、つまり 5 年間 (月単位で測定)、C は患者が 5 年間に死亡したかどうかを示します。病院は、病院Aの場合は 1 、病院Bの場合は 0というバイナリ変数Xとしてエンコードされています。
私たちの単一共変量 Cox 比例モデルは次のようになります。病院の影響を表し、各患者にインデックスを付けます。
統計ソフトウェアを使用すると、2.12 と推定できます。ハザード比はこの値の指数です。理由を理解するには、ハザード比を具体的に考えてみましょう。
したがって、病院 A と病院 B のハザード比は です。統計的有意性はさておき、病院 A の患者は病院 B の患者に比べて、短期間で死亡するリスクが 8.3 倍高いと言えます。
解釈に関して言及すべき重要な注意点があります。
- 死亡リスクが 8.3 倍高いということは、病院で死亡する患者が 8.3 倍多いことを意味するわけではありません。A: 生存分析では、イベントが発生するかどうかではなく、イベントが発生する速さを調べます。
- より具体的に言うと、「死亡リスク」は速度の尺度です。速度には、メートル/秒などの単位があります。ただし、相対速度には単位がありません。つまり、自転車は別の自転車 (基準自転車) の 2 倍の速さで走ることができますが、単位は指定されません。同様に、病院Aでの死亡リスク (自転車の速度に相当) は、病院B (基準グループ)での死亡リスクより 8.3 倍 (速い) 高くなります。
- 逆数は病院Aに対する病院Bのハザード比です。
- 病院間の生存確率については、まだ推測していません。これは、推定値だけでなく、ベースライン ハザード率の推定値も必要になるからです。ただし、Cox 比例ハザード モデルの標準的な推定では、ベースライン ハザード率を直接推定することはできません。
- モデルの唯一の時間変動要素であるベースライン ハザード率を無視したため、推定値は時間スケールに依存しません。たとえば、時間を月ではなく年で測定した場合、同じ推定値が得られます。
- 病院が2 つのグループ間の危険性の差を引き起こしたと言いたくなりますが、私たちの研究は因果関係を問うものではないため (つまり、データがどのように生成されたかがわからないため)、「関連」などの用語を使用します。
単一の連続共変量
生存分析のあまり伝統的ではない使用例を説明するために、次の例は経済学の質問です。企業の IPO 初周年における株価収益率 (P/E) と将来の生存率にはどのような関係があるのでしょうか。より具体的には、企業の「誕生イベント」を IPO 初周年と見なし、破産、売却、非公開化などを企業の「死」イベントと見なすと、企業の「誕生」(IPO 初周年) における P/E 比率が企業の生存率にどのような影響を与えるかを知りたいと考えます。
12 社の生存データを含む (偽の) データセットが提供されています。Tは、 IPO 1 周年から消滅までの日数 (消滅しなかった場合は終了日 2022-01-01) を表します。Cは、会社が 2022-01-01 より前に消滅したかどうかを表します。P/E は、IPO 1 周年における会社の株価収益率を表します。
バイナリ変数があった前の例とは異なり、このデータセットには連続変数 P/E があります。ただし、モデルは似ています。 ここで、 は会社の P/E 比率を表します。このデータセットを Cox モデルに通すと、未知の値の推定値が生成され、-0.34 になります。したがって、ハザード全体の推定値は次のとおりです。
ベースラインハザードは推定されていないため、ハザード全体を計算することはできません。ただし、企業iとj のハザードの比率を考えてみましょう。
右側のすべての項は既知であるため、企業間のハザード比を計算することができます。右側には時間に依存する項がないため(すべての項は定数)、ハザードは互いに比例します。たとえば、会社 5 と会社 2 のハザード比は です。これは、調査期間中、会社 5 の「死亡」リスクが 0.33 ≈ 会社 2 の死亡リスクの 1/3 であることを意味します。
解釈に関して言及すべき重要な注意点があります。
- ハザード比は、上記の例にある 量 です。上記の最後の計算から、これは、変数が 1 単位異なる 2 つの「対象」間のハザード比として解釈できます。つまり、 の場合、 です。「1 単位異なる」という選択は、 の値を正確に伝えるので便利です。
- ベースライン ハザードは、スケーリング係数が 1、つまり のときに表すことができます。ベースライン ハザードを、P/E がたまたま 0 である「ベースライン」企業のハザードとして解釈できますか? ベースライン ハザードを「ベースライン サブジェクトのハザード」として解釈することは不完全です。共変量が 0 になることはこのアプリケーションでは不可能だからです。P/E が 0 であることは意味がありません (つまり、企業の株価が 0 であること、つまり、企業が「死んでいる」ことを意味します)。より適切な解釈は、「すべての変数がゼロである場合のハザード」です。
- のような値を理解して解釈し、企業のリスクを表すようにしたくなるのは当然です。しかし、これが実際に何を表しているか考えてみましょう。ここでは、企業i のリスクを P/E が 0 の架空のベースライン企業と比較する、リスクの比率が暗黙的に存在します。しかし、上で説明したように、このアプリケーションでは P/E が 0 になることはあり得ないため、この例では意味がありません。ただし、考えられるリスク間の比率には意味があります。
時間変動予測子と係数
時間依存変数、時間依存層、および被験者ごとの複数のイベントへの拡張は、アンダーセンとギルの計数プロセスの定式化によって組み込むことができます。[8]時間変動回帰変数を使用したハザードモデルの使用例の1つは、失業保険が失業期間に与える影響を推定することです。[9] [10]
Cox モデルは、時間とともに変化する共変量(予測変数)を許容するだけでなく、時間とともに変化する係数にも一般化できます。つまり、治療の比例効果は時間とともに変化する可能性があります。たとえば、薬は発病から 1 か月以内に投与された場合は非常に効果的ですが、時間が経つにつれて効果が低下します。その後、係数が時間とともに変化しない (定常性) という仮説をテストできます。詳細とソフトウェア ( R パッケージ) は Martinussen と Scheike (2006) で入手できます。[11] [12]
この文脈では、理論的には、加法ハザード[13]を使用して共変量の影響を特定することも可能であることも言及できます。つまり 、 このような加法ハザードモデルを(対数)尤度最大化が目的である状況で使用する場合、非負の値に制限するように注意する必要があります。おそらくこの複雑さの結果として、このようなモデルはほとんど見られません。目的が最小二乗である場合、非負の制限は厳密には必要ありません。
ベースラインハザード関数の指定
ベースライン ハザードが特定の形式に従うと仮定する理由がある場合、Cox モデルを特殊化することができます。この場合、ベースライン ハザードは特定の関数に置き換えられます。たとえば、ハザード関数をワイブルハザード関数と仮定すると、ワイブル比例ハザード モデルが得られます。
ちなみに、ワイブル ベースライン ハザードを使用することは、モデルが比例ハザードと加速故障時間モデルの両方を満たす唯一の状況です。
パラメトリック比例ハザード モデルという一般的な用語は、ハザード関数が指定されている比例ハザード モデルを表すために使用できます。対照的に、Cox 比例ハザード モデルはセミパラメトリック モデルと呼ばれることもあります 。
一部の著者は、基礎となるハザード関数を特定する場合でも、Cox比例ハザードモデルという用語を使用し、 [14]この分野全体がDavid Coxに負っていることを認めています。
Cox 回帰モデル(比例ハザードを省略)という用語は、時間依存因子を含むように Cox モデルを拡張したものを表すために使用されることがあります。ただし、Cox 比例ハザード モデル自体が回帰モデルとして説明できるため、この用法は潜在的に曖昧です。
ポアソンモデルとの関係
比例ハザード モデルとポアソン回帰モデルの間には関係があり、ポアソン回帰のソフトウェアで近似比例ハザード モデルを当てはめるために使用されることがあります。通常、これを行う理由は、計算がはるかに高速になるためです。これは、コンピュータが低速だった時代にはより重要でしたが、特に大規模なデータ セットや複雑な問題には今でも役立ちます。Laird と Olivier (1981) [15] は数学的な詳細を示しています。彼らは、「[ポアソン モデル] が正しいと仮定するのではなく、単に尤度を導き出す手段としてそれを使用する」と述べています。McCullagh と Nelder [16]の一般化線型モデルに関する本には、比例ハザード モデルを一般化線型モデルに変換する章があります。
高次元設定下
高次元では、共変量数pがサンプルサイズnに比べて大きい場合、LASSO法は古典的なモデル選択戦略の1つです。Tibshirani(1997)は、比例ハザード回帰パラメータのLasso手順を提案しました。[17]回帰パラメータβのLasso推定量は、 L 1ノルム型制約の下でCox部分対数尤度の反対を最小化するものとして定義されます。
最近、このテーマに関して理論的な進歩がありました。[18] [19] [20] [21]
ソフトウェア実装
- Mathematica :
CoxModelFit関数。 [22] - R : survival
coxph()パッケージにある関数。 - SAS :
phreg手順 - Stata :
stcoxコマンド - Python : lifelines
CoxPHFitterライブラリにあります。statsmodelsライブラリにあります。phreg - SPSS : Cox 回帰で利用できます。
- MATLAB :
fitcoxまたはcoxphfit関数 - Julia : Survival.jlライブラリで利用可能です。
- JMP :比例ハザードの適合プラットフォームで利用できます。
- Prism : 生存分析と多変量分析で利用可能
参照
注記
- ^ Breslow, NE (1975). 「比例ハザードモデルによる生存データの分析」.国際統計評論 / Revue Internationale de Statistique . 43 (1): 45–57. doi :10.2307/1402659. JSTOR 1402659.
- ^ Cox, David R (1972). 「回帰モデルと生命表」.英国王立統計学会誌、シリーズ B. 34 ( 2): 187–220. JSTOR 2985181. MR 0341758.
- ^ Kalbfleisch, John D.; Schaubel, Douglas E. (2023年3月10日). 「Coxモデルの50年」.統計とその応用の年次レビュー. 10 (1): 1–23. Bibcode :2023AnRSA..10....1K. doi : 10.1146/annurev-statistics-033021-014043 . ISSN 2326-8298.
- ^ Reid, N. (1994). 「サー・デイヴィッド・コックスとの対話」統計科学9 ( 3): 439–455. doi : 10.1214/ss/1177010394 .
- ^ Cox, DR (1997).生存データの分析に関するいくつかのコメント。第1回シアトル生物統計シンポジウム:生存分析。
- ^ 「各失敗は尤度関数に寄与する」、Cox (1972)、191 ページ。
- ^ エフロン、ブラッドリー (1974)。「検閲データに対するコックス尤度関数の効率」アメリカ統計学会誌。72 ( 359 ) : 557–565。doi : 10.1080/01621459.1977.10480613。JSTOR 2286217。
- ^ Andersen, P.; Gill, R. (1982). 「Coxの回帰モデルによる計数プロセス、大規模サンプル研究」Annals of Statistics . 10 (4): 1100–1120. doi : 10.1214/aos/1176345976 . JSTOR 2240714.
- ^ マイヤー、BD (1990)。「失業保険と失業期間」( PDF)。エコノメトリカ。58 (4): 757–782。doi :10.2307/2938349。JSTOR 2938349 。
- ^ Bover, O.; Arellano, M .; Bentolila, S. (2002). 「失業期間、給付期間、およびビジネスサイクル」(PDF) . The Economic Journal . 112 (479): 223–265. doi :10.1111/1468-0297.00034. S2CID 15575103.
- ^ マルティヌッセン;シャイケ (2006)。生存データの動的回帰モデル。スプリンガー。土井:10.1007/0-387-33960-4。ISBN 978-0-387-20274-7。
- ^ 「timereg: 生存データのための柔軟な回帰モデル」CRAN。
- ^ Cox, DR (1997).生存データの分析に関するいくつかのコメント。第1回シアトル生物統計シンポジウム:生存分析。
- ^ Bender, R.; Augustin, T.; Blettner, M. (2006). 「Cox比例ハザードモデルをシミュレートするための生存時間の生成」。Statistics in Medicine . 24 (11): 1713–1723. doi : 10.1002/sim.2369 . PMID 16680804. S2CID 43875995.
- ^ Nan Laird および Donald Olivier ( 1981)。「対数線形分析手法を用いた打ち切り生存データの共分散分析」。アメリカ統計学会誌。76 (374): 231–240。doi :10.2307/2287816。JSTOR 2287816。
- ^ P. McCullagh および JA Nelder (2000)。「第 13 章: 生存データのためのモデル」。一般化線形モデル(第 2 版)。フロリダ州ボカラトン: Chapman & Hall/ CRC。ISBN 978-0-412-31760-6。(第 2 版 1989 年、CRC による最初の再版 1999 年)
- ^ Tibshirani, R. (1997). 「Coxモデルにおける変数選択のためのLasso法」.医学統計. 16 (4): 385–395. CiteSeerX 10.1.1.411.8024 . doi :10.1002/(SICI)1097-0258(19970228)16:4<385::AID-SIM380>3.0.CO;2-3. PMID 9044528.
- ^ Bradić, J.; Fan, J.; Jiang, J. (2011). 「NP次元によるCoxの比例ハザードモデルの正規化」Annals of Statistics . 39 (6): 3092–3120. arXiv : 1010.5233 . doi :10.1214/11-AOS911. PMC 3468162. PMID 23066171 .
- ^ Bradić, J.; Song, R. (2015). 「ノンパラメトリックCoxモデルにおける構造化推定」.電子統計ジャーナル. 9 (1): 492–534. arXiv : 1207.4510 . doi :10.1214/15-EJS1004. S2CID 88519017.
- ^ Kong, S.; Nan, B. (2014). 「Lasso による高次元 Cox 回帰の非漸近的オラクル不等式」. Statistica Sinica . 24 (1): 25–42. arXiv : 1204.1992 . doi :10.5705/ss.2012.240. PMC 3916829. PMID 24516328 .
- ^ Huang, J.; Sun, T.; Ying, Z.; Yu, Y.; Zhang, CH (2011). 「Cox モデルにおける Lasso の Oracle 不等式」. The Annals of Statistics . 41 (3): 1142–1165. arXiv : 1306.4847 . doi :10.1214/13-AOS1098. PMC 3786146. PMID 24086091 .
- ^ "CoxModelFit". Wolfram言語およびシステムドキュメントセンター。
参考文献
- Bagdonavicius, V.; Levuliene, R.; Nikulin, M. (2010). 「左切断データおよび右打ち切りデータからの Cox モデルの適合度基準」。Journal of Mathematical Sciences . 167 (4): 436–443. doi :10.1007/s10958-010-9929-6. S2CID 121788950.
- Cox, DR; Oakes, D. (1984)。生存データの分析。ニューヨーク:チャップマン&ホール。ISBN 978-0412244902。
- コレット、D. (2003)。『医学研究における生存データのモデリング(第2版)』ボカラトン:CRC。ISBN 978-1584883258。
- Gouriéroux, Christian (2000)。「持続モデル」。質的従属変数の計量経済学。ニューヨーク: Cambridge University Press。pp. 284–362。ISBN 978-0-521-58985-7。
- Singer, Judith D.; Willett, John B. (2003)。「Cox 回帰モデルの適合」。応用縦断的データ分析: 変化とイベント発生のモデル化。ニューヨーク: Oxford University Press。pp. 503–542。ISBN 978-0-19-515296-8。
- Therneau, TM; Grambsch, PM (2000).生存データのモデリング: Cox モデルの拡張. ニューヨーク: Springer. ISBN 978-0387987842。
