時間依存共変量を伴う区間打ち切りCoxモデル
Stata 17 では、stintcoxコマンドが導入され、区間打ち切りのあるイベント時間データに対してセミパラメトリックCoxモデルを適合させる機能が追加されました。
Stata 18 では、フォローアップ期間中に変化する特性を記録するために一般的に使用される時間依存共変量（TVCs）のサポートが追加されました。
したがって、生物マーカー変数や政策変数など、時間とともに変化する特性をモデルに組み込むことができます。例えば、がんの再発やフォローアップ中に記録される雇用状況の変化など、
ある時間区間内にのみ発生する興味のあるイベントをモデル化することができます。いずれも、イベント時間が正確に観察されていない区間打ち切りのあるイベント時間データの例です。
stintcoxを使用して、既存の共変量を時間の指定された確定関数と対話させてTVCを自動的に作成するか、新しいサブジェクトごとの複数レコードのデータ形式で独自のTVCを指定します。
TVCを使用して比例ハザード仮定をWald検定や尤度比検定でテストします。そして、TVCを予測や生存者関数などのプロットに組み込みます。
TVC（時変共変量）を使用して比例ハザードの仮定を検証しよう
ここでは、TVC（時変共変量）を使用して、Coxモデルの比例ハザード仮定を検証する方法を示します。Zengら（2016）に従い、注射薬物使用者のHIV感染までの時間をモデル化するために、区間打ち切りCox回帰を使用します。
1,124人の対象者は最初にHIV陰性であることが確認されました。約4か月ごとに血液検査を行い、HIV-1陽性かどうかを評価しました。対象者は定期的に検査されたため、HIV感染の正確な時間は観察されませんでしたが、
感染した時間は血液検査の特定の間隔に含まれることが分かっています。対応する下限時間と上限時間は、変数ltimeとrtimeに記録されています。
HIV-1陽性をモデル化するために使用する基準となる要因は、リクルート時の中心化された年齢（age_mean）、性別（male）、針の共有履歴（needle）、薬物注射の記録（inject）、およびリクルート時に被験者が
刑務所にいたかどうか（jail）です。記録は被験者ごとに1行にまとめられており、このデータセットの一部である271から274行目の観測値を表示します。
. use https://www.stata-press.com/data/r18/idu. format ltime rtime age_mean %6.2f. list ltime rtime age_mean male needle inject jail if _n >= 271 & _n <= 274
上記の基準要因に依存するHIV感染までの時間について、Cox比例ハザードモデルを適合させます。各対象者に1つのレコードがある形式で、stintcoxのinterval()オプションで
下限および上限のイベント時間区間を指定します。
. stintcox age_mean i.male i.needle i.inject i.jail, interval(ltime rtime)
年齢が上がるとHIV感染のリスクが低くなり、また、刑務所に滞在しているとHIV感染のリスクが高まるといった傾向が観察されます。
共変量の比例ハザード仮定を検証する方法の一つは、その共変量に関連する係数が時間に依存しないかどうかを検定することです。これは、モデル内でこの共変量と時間の関数との相互作用を含め、
対応する係数がゼロと等しいかどうかを検定することで実現できます。
この例で比例ハザード仮定を検証するために、tvc()オプションに全ての共変量を含め、それらと解析時間 _t との相互作用をモデルに追加します（デフォルトのオプション）。
この分析では、nohrオプションを使用してハザード比のデフォルトの報告を抑制します。
. stintcox age_mean i.male i.needle i.inject i.jail, interval(ltime rtime) ///
> tvc(age_mean i.male i.needle i.inject i.jail) nohr
推定結果の表の上段（main）は、時間と変動しない共変量の係数を報告しています。推定結果の表の下段（tvc）は、時間と相互作用する共変量の結果を報告しています。 そして、tvcで報告されたp値に基づいて、共変量maleおよびjailに関しては比例ハザード仮定が破られているように見えます。 推定結果の表の下部にあるWald検定も示しているように、比例ハザード仮定が全体的に成り立たないことが示されています。
複数の記録データに対してTVCsを使用して区間打ち切りCoxモデルを適合させる
前述の例では、比例ハザード仮定を検証するために、stintcoxのtvc()オプションを使用しました。
いくつかの応用では、データセットにすでにTVCsが存在していることがあります。注射薬使用者のデータを使用して、刑務所にいるかどうかはTVCです。
TVCsは、個体ごとに複数のレコードがあるデータ形式でのみ記録できます。
この形式では、各個体には複数の検査回数、各検査回数でのイベントステータス、および各検査回数での潜在的な時間変動共変量が含まれる可能性があります。
前回のidu.dtaデータセットの拡張バージョンであるidu2.dtaを使用しましょう。このデータにはベースラインの共変量としてのすべてが含まれており、
時間変動する刑務所への入所を示す変数（jail_vary）も含まれています。jail_varyは、前回のクリニック訪問以来、対象者が刑務所に入所したかどうかを示します。
データセットにはさらに、対象者の識別子（id）、血液検査の検査時間（time）、および検査時間での血液検査が陽性であるかどうか（is_seropos）も記録されています。
以下は、対象者271から274までのサブセットの例です。
. use https://www.stata-press.com/data/r18/idu2
. format time age_mean %6.2f
. list id time is_seropos age_mean male needle inject jail_vary ///
> if id >= 271 & id <= 274, sepby(id) noobs abbreviate(10)
前回のCoxモデルを再適合させますが、今回はベースライン変数の代わりに時間変動変数jail_varyを使用します。
サブジェクトごとの複数のレコード形式では、id()オプションで対象者の識別子、time()オプションで検査時間、
status()オプションでイベント指標のステータスを指定します。
. stintcox age_mean i.male i.needle i.inject i.jail_vary, id(id) time(time) ///
> status(is_seropos)
前の例と比較して、時間変動する刑務所への入所の影響を考慮すると、注射のハザード比は1.25から1.37に増加しますが、 刑務所への入所の影響は、ベースラインのjailの場合の1.57から、時間変動するjail_varyの場合は1.44に減少します。
生存関数のグラフ化
モデルを適合させた後、stcurveを使用して推定された生存、故障、ハザード、または累積ハザード関数をプロットできます。
デフォルトでは、stcurveは共変量の全体の平均で関数を評価します。
. stcurve, survival
複数のデータレコードが時間変動する共変量を持っている場合、関数をプロットする際には、それらの共変量が時間とともに変化する性質を反映させたい場合があるかもしれません。
その際、オプションattmeansを使用して、関数を時間ごとの平均で評価することができます。
. stcurve, survival attmeans
TVC（時間変動共変量）の値を指定して、生存関数を評価するためには、オプションatframe()を使用します。
例えば、データセット内のidが1と同じ共変量パターンを持つ個体の生存曲線を描きたい場合、新しいフレーム（id1と呼ぶ）を作成し、
frame putを使ってsubject 1の関連情報をこの新しいフレームにコピーします。フレームid1に保存したデータを一覧表示します。
. frame put time age_mean male needle inject jail_vary if id==1, into(id1)
. frame id1: list
そして、下記のコマンドを実行して生存曲線のグラフを表示します。
. stcurve, survival atframe(id1)
