区間打ち切り多重イベントデータの
周辺 Cox 比例ハザードモデル

　糖尿病や高血圧の発症など、複数種類のイベントについて発生までの時間を分析したいものの、正確な発生時点がわからない場合があります。新しい stmgintcox コマンドを使用すると、このような区間打ち切り多重イベントデータを分析し、異なるイベントの発生時間どうしに存在しうる相関を考慮できます。比例ハザード仮定の評価、全イベントに共通する共変量効果の検定、共変量別の生存関数・ハザード関数などのグラフ化も可能です。

　stmgintcox は、新しい周辺 Cox 比例ハザードモデル（Xu, Zeng, and Lin 2023）を実装しています。

区間打ち切り多重イベントデータとは

　区間打ち切り多重イベントデータ（より正確には、複数イベントの発生時間が区間打ち切りされたデータ）は、縦断研究でよく現れます。各対象者が複数種類のイベントを経験し、それぞれの発生時点を直接観測できず、ある時間区間内に発生したことだけがわかるためです。

　この形式のデータは、医学、疫学、生物学、社会学など多くの分野で生じます。たとえば疫学研究では、複数回の診察を通じて心疾患や代謝疾患を追跡します。社会学では、転職や結婚などのライフイベントを定期調査で記録します。生態学では、営巣や出産などの繁殖サイクルを定期観測します。研究者は要因がイベント発生時間に及ぼす影響に関心を持ちますが、正確な発生時間が観測されず、異なる発生時間間の依存構造も不明なことが多いため、分析は容易ではありません。

　周辺比例ハザードモデルは、このようなデータの分析に利用できます。イベント間の依存構造をモデル化する必要がないため、より頑健な推論が可能であり、パラメータは母集団平均効果として解釈できます。また、ランダム効果モデルより計算が速い場合も少なくありません。

　Stata 17では、単変量の区間打ち切りイベント時間データに対して真正のセミパラメトリック Cox モデルを当てはめる stintcox を導入しました。Stata 18では、時間依存共変量（TVC）に対応しました。新しい stmgintcox は、区間打ち切り多重イベントデータに周辺比例ハザードモデルを当てはめます。イベント当たり1レコードと複数レコードの両形式に対応し、すべてまたは一部のイベントにTVCを指定できます。イベント別共変量の柔軟な指定、全イベントに共通する共変量効果の推定・検定、イベント別予測、生存関数などのグラフ、適合度グラフも利用できます。

操作例

• イベント当たり1レコード形式のモデル
• イベント別共変量を含むモデル
• 全イベントに共通する共変量効果の推定と検定
• 生存関数のグラフ化
• 適合度グラフによるモデル全体の評価
• イベント当たり複数レコード形式のモデル

イベント当たり1レコード形式の区間打ち切りデータにモデルを当てはめる

　Xu, Zeng, and Lin（2023）で説明されている Atherosclerosis Risk in Communities（ARIC）研究をもとにシミュレーションした仮想データを使用します。データには米国4地域の200人が含まれ、追跡期間中の複数回の診察で糖尿病と高血圧の有無が評価されています。診察は定期的にしか行われないため、疾患の正確な発症時間は不明ですが、診察間のどの区間に発症したかはわかります。変数 ltime は発症前最後の診察時間、rtime は発症後最初の診察時間を記録しています。

　糖尿病と高血圧の発症時間に影響する要因を調べます。対象とするのは、人口統計学的変数 （race）、性別（male）、居住地域（community） と、5つのベースラインのリスク因子、年齢（ age）、BMI（bmi）、血糖値（glucose）、収縮期血圧（sysbp）、拡張期血圧（diabp）です。各対象者・各イベントにつき1レコードがあり、イベント時間は区間データとして記録されています。対象者IDが91と92の一部を確認します。

　上記の要因に応じて糖尿病と高血圧の発症時間が変化する周辺 Cox 比例ハザードモデルを当てはめます。イベント当たり1レコード形式では、id()、event()、interval() オプションを指定します。このコマンドは計算負荷が高いため、実行に時間がかかる場合があります。

　係数表の上には打ち切り情報の要約が表示されます。糖尿病については、Whiteの人で発症リスクが低く、BMIと血糖値が高いほど発症リスクが高いことがわかります。高血圧については、収縮期血圧と拡張期血圧が高いほど発症リスクが高くなっています。

イベント別共変量を含むモデルを当てはめる

　上のモデルから、BMIと血糖値は糖尿病の主要なリスク因子ですが、高血圧には重要でないことがわかります。一方、収縮期血圧と拡張期血圧は高血圧には重要ですが、糖尿病には重要でないようです。そこで、2つのイベントに異なる共変量集合を使用します。また、オプションについて反復計算ログを省略する nolog と、計算を高速化する favorspeed も指定します。

　age、male、community、race は両イベントに共通するため、共通共変量とイベント別共変量を組み合わせて次のようにも指定できます。

　上の2通りの指定からは同じ推定結果が得られます。

全イベントに共通する共変量効果を推定・検定する

　上のモデルを当てはめた後、age の効果がすべてのイベントで0であるという仮説を検定します。estat common は、全イベントにわたる age の最適な加重平均効果を推定し、帰無仮説のもとで平均効果が0かどうかを z 検定します。

　全イベントにわたる age の加重平均効果が0であるという帰無仮説は棄却できません（p 値 = 0.518）。イベント間で効果が似ている場合、この平均効果は共通効果を推定し、test コマンドによる従来の多変量 Wald 検定より高い検出力を持ちます。

生存関数をグラフ化する

　stcurve で推定生存関数を描画できます。survival オプションを指定すると、デフォルトでは各イベントの共変量を全体平均に設定して生存関数を評価し、両イベントをサブグラフとして表示します。

　sepevents を追加すると、各イベントの推定生存関数を別々のグラフに表示できます。地域別に、平均的な人の糖尿病生存曲線を比較するには、at() に community の複数の値を指定し、events() にイベント値ラベル "Diabetes" を指定します。

　生存曲線から、Forsyth（青線）と Minneapolis（緑線）の平均的な人は糖尿病発症リスクが同程度に高く、Washington（黄線）と Jackson（赤線）ではリスクが低いことがわかります。

適合度グラフでモデル全体の適合を評価する

　estat gofplot はイベント別の適合度グラフを作成し、モデル全体の適合を視覚的に評価します。イベントごとに Cox–Snell 型残差に基づく累積ハザード関数の経験的推定値を計算し、累積ハザード率を残差に対してプロットします。モデルがデータに適合していれば、プロットは基準線の近くに位置します。

　デフォルトでは、すべてのイベント別適合度グラフが1つのグラフ内にサブグラフとして表示されます。左のプロットから、糖尿病については周辺 Cox 比例ハザードモデルがデータによく適合していることがわかります。右のプロットでは、高血圧についても概ね良好ですが、外れ値の影響で裾の部分が基準線から外れています。

イベント当たり複数レコード形式の区間打ち切りデータにモデルを当てはめる

　区間打ち切りデータは、イベント当たり複数レコード形式でも記録できます。この形式では、対象者の各イベントについて、複数の診察時間と各時点のTVCを複数レコードで保持できます。ここでは先ほどのデータを拡張し、ベースライン共変量に加えて、bmi、glucose、sysbp、diabp の4つの時間依存共変量を含むデータを使用します。診察時間は time、前回診察以降にイベントが発生したかどうかは status に記録されています。

　糖尿病と高血圧の発症時間が、ベースライン共変量 age、male、community、race と、時間依存共変量 bmi、glucose、sysbp、diabp に依存する周辺 Cox 比例ハザードモデルを当てはめます。イベント当たり複数レコード形式では、id()、event()、time()、status() を指定します。さらに detail を指定し、各イベントの詳細な打ち切り情報を表示します。

　時間とともに変化するリスク因子を考慮した結果、血糖値、BMI、収縮期血圧が高いほど糖尿病の発症リスクが高いことがわかります。Whiteは糖尿病発症リスクの低下と関連しています。また、収縮期血圧と拡張期血圧が高いほど高血圧発症リスクが高くなっています。

参考文献

Xu, Y., D. Zeng, and D. Lin. 2023. Marginal proportional hazards models for multivariate interval-censored data. Biometrika 110: 815–830.

The ARIC investigators. 1989. The Atherosclerosis Risk in Communities (ARIC) study: Design and objectives. American Journal of Epidemiology 129: 687–702. https://doi.org/10.1093/oxfordjournals.aje.a115184.

さらに詳しく

　多変量区間打ち切りデータに周辺 Cox 比例ハザードモデルを当てはめる方法については、[ST] stmgintcox および Stata Survival Analysis Reference Manual をご覧ください。