周辺構造モデル
Marginal Structural Models (MSM)とは疫学で主に用いられる統計モデルの一つで、縦断データの分析に有用である。ここでは原典とも言える以下のEpidemiologyの論文をじっくり読みながらMSMについての理解を深めていく。
要約
経時的に変化する曝露や治療を伴う観察研究では、時間依存性の交絡因子が存在し、それが以前の治療にも影響される場合、交絡の調整のための標準的なアプローチはバイアスがかかる。本論文では、このような状況において交絡の調整を改善することを可能にする因果モデルの新しいクラスである周辺構造モデルを紹介する。周辺構造モデルのパラメータは、新しい推定量クラスであるinverse-probability-of-treatment重み付け推定量を用いて一貫して推定することができる。
キーワード：因果推論、半事実、疫学手法、縦断データ、構造モデル、交絡、中間変数
はじめに
周辺構造モデル（MSM）は、観察データから、同時に交絡因子や中間変数となりうる時間依存性の共変量が存在する場合の、時間依存性の曝露の因果効果を推定するための新しいクラスの因果モデルである。MSMのパラメータは、IPTW（inverse-probability-of-treatment weighted）推定量という新しい推定量を用いて一貫して推定することができる。MSMは、構造的入れ子モデル（SNM）に代わるものであり、そのパラメータはg推定の手法によって推定される。
時間的に変化する曝露や治療の効果を推定するための通常のアプローチは、層別解析やそのパラメトリック類似体（例えば、ロジスティック回帰や比例ハザード回帰）のような解析手法を用いて、病気の確率を過去の曝露や過去の交絡因子の履歴の関数としてモデル化することである。セクション4と7.1で、(1)対象イベントの危険因子または予測因子であり、その後の曝露も予測する時間依存性の共変量が存在し、(2)過去の曝露履歴がその共変量のその後のレベルを予測する場合、分析において過去の交絡因子の履歴をさらに調整するかどうかにかかわらず、これらの標準的なアプローチが偏る可能性があることを示す。条件1を満たす共変量を時間依存交絡因子と呼ぶ。条件1と2は、多くの観察研究、特に徴候による交絡がある研究で当てはまる。例えば、ジドブジン（AZT）治療がヒト免疫不全ウイルス（HIV）感染者の死亡率に及ぼす影響に関する研究では、時間依存性の共変量であるCD4リンパ球数は、生存率とAZTによる治療開始の独立した予測因子であると同時に、それ自体がAZT治療歴の影響を受ける。死亡率に対する肥満の影響に関する研究では、臨床的心疾患または呼吸器疾患の発症は、死亡率とその後の体重減少の独立した予測因子であり、事前の体重増加に影響される。交絡因子と中間変数である時間依存性の共変量が同時に存在する場合、条件1と2は常に成立する。標準的な方法のバイアスについてのより詳細な説明と、時間依存性交絡の疫学的な例がいくつか追加されている。
1. 時間依存交絡
HIV感染者の追跡調査を考えてみよう。\(A_k\)を、追跡調査開始後\(k\)日目の、目的の治療または曝露（例えばAZT）の投与量とする。\(K+1\)日目の追跡調査終了時に測定された、関心のある二値の結果を\(Y\)とする（例えば、血液中のHIV RNAが検出されなければ\(Y = 1\)、そうでなければ\(Y = 0\)）。我々の目的は、アウトカム\(Y\)に対する時間依存の治療\(A_k\)の因果効果を推定することである。
図1は、\(k = 1\)での我々の研究を表すDAGである。DAGは、グラフの頂点（ノード）が変数を表し、有向辺（矢印）が直接的な因果効果を表す有向無サイクルグラフである。図1では、\(L_k\)は、年齢、CD4リンパ球数、白血球数（WBC）、ヘマトクリット、後天性免疫不全症候群（AIDS）の診断、さまざまな症状や口腔カンジダ症などの日和見感染の有無など、転帰に関するすべての測定された危険因子のベクトルの\(k\)日目の値を表す。同様に，\(U_k\)は，\(Y\)に関するすべての未測定の原因リスク因子の\(k\)日目の値を表す．図1のbは，未測定の原因リスク因子から治療変数への矢印が削除されている点でのみ、図1のaと異なる．図1のbのように、測定されていない原因リスク因子から治療変数への矢印がない場合、測定された交絡因子\(L_k\)のデータがあれば、測定されていない交絡因子はないと言える。図1のcは、\(Y\)のリスク因子（測定されたもの、測定されていないもの）のどれもが、どの治療変数への矢印も持っていないという点で、図1のa,bと異なる。しかし，初期の治療\(A_0\)は，後の治療\(A_1\)に因果的に影響しうることに注意する。図1のcのように、どの（非治療）危険因子からどの治療変数への矢印もない場合、測定された因子または測定されていない因子のいずれによる交絡もなく、この場合治療は非交絡であると言う。
上記の区別は、治療が時間に依存しない、より馴染みのある点治療研究にも同様に適用される。図2に示されるように、点治療研究は、\(k = 0\)の一般的なセットアップの特別なケースである。図2のa-cは、点治療研究の図1のa-cのアナログを含んでいる。
どの観察研究でもそうであるように、\(L_k, A_k, Y\)に関する観察データから、測定されていない危険因子による交絡があるかどうかを判断することはできない。\(U_k\)による交絡が残っていても、それが小さいことを願うしかない。しかし、\(L_k\)があれば他に測定不能の交絡はないという検証不可能な仮定の下では、治療が無交絡であるかどうかをデータから経験的に検証することができる。具体的には、治療が無交絡であるための十分条件は、各時刻\(k\)において、同じ過去の治療履歴\(A_0, \cdots, A_{k-1}\)を持つ被験者の間で、治療\(A_k\)が測定された共変量\(L_0, \cdots, L_k\)の過去の履歴と関連していないことである。例えば、図2のポイント治療研究で、\(A_0\)が\(L_0\)と関連していなければ、治療は無交絡である。
2. ポイント治療研究における反事実性
我々は、図2の点治療研究において、\(Y\)に対する\(A_0\)の効果をどのように推定するかを検討することから始める。さらに、図2のcが真の因果グラフであり、測定された共変量も測定されていない共変量も治療と結果の関係を交絡させないとする。そして、粗リスク差、リスク比、オッズ比はそれぞれ異なる尺度ではあるが、結果\(Y\)に対する治療\(A_0\)の因果効果を測定する。粗リスク差は\(cRD = pr[Y = 1| A_0 = 1]-pr[Y = 1|A_0 = 0]\)であり、粗リスク比は\(cRR = pr[Y = 1|A_0 = 1]/pr[Y = 1|A_0 = 0]\)である。粗オッズ比は、\(\displaystyle cOR = \frac{pr[Y = 1|A_0 = 1]pr[Y = 0|A_0 = 0]}{pr[Y = 1|A_0 = 0]pr[Y = 0|A_0 = 1]}\)であり、例えば、\(pr[Y = 1|A_0 = 1]\)は、治療された被験者\(A_0 = 1\)の中で\(Y = 1\)である確率である。調査対象者は、大規模な、おそらく仮想の、元集団からの無作為標本であると仮定する。確率は元集団における比率を意味する。
これらの関連パラメータに対応する因果対比には、反事実変数が含まれる。具体的には、変数\(Y_{a_0 = 1}\)は、治療した場合の被験者の結果を表し、\(Y_{a_0 = 0}\)は、治療しなかった場合の被験者の結果を表す。ある被験者について、差尺度で測定した治療の因果効果は、\(Y_{a_0 = 1}- Y_{a_0 = 0}\)である。対象者が処置された場合（\(A_0 = 1\)）、対象者の観察結果\(Y\)は\(Y_{a_0 = 1}\)に等しく、\(Y_{a_0 = 0}\)は観察されない。\(A_0 = 0\)の場合、\(Y\)は\(Y_{a_0 = 0}\)に等しく、\(Y_{a_0 = 1}\)は観察されない。それぞれ\(pr(Y_{a_0 = 1} = 1)\)と\(pr(Y_{a_0 = 0} = 1)\)を、\(Y_{a_0 = 1}\)が\(1\)に等しい確率と\(Y_{a_0 = 0}\)が\(1\)に等しい確率とすると、治療\(A_0\)が無交絡であれば、粗RDは、元集団における因果リスク差\(pr(Y_{a_0 = 1} = 1)-pr(Y_{a_0 = 1} = 1)\)に等しい。同様に、粗RRは因果RR、\(pr(Y_{a_0 = 1} = 1)/pr(Y_{a_0 = 1} = 1)\)に等しく、粗ORは因果OR、\(\displaystyle \frac{pr(Y_{a_0=1} = 1)pr(Y_{a_0=0} = 0)}{{pr(Y_{a_0=1} = 0)pr(Y_{a_0=0} = 1)}}\)に等しい。効果修正の可能性があるので、母集団の因果パラメータは、治療が無交絡であっても、測定された危険因子\(L_0\)の層内の因果パラメータと等しくなる必要はない。効果修正については第9章で検討する。
3. 点治療研究のモデル
因果RD、RR、ORは、2つの反事実確率\(pr(Y_{a_0 = 1} = 1)\)と\(pr(Y_{a_0 = 0} = 1)\)に対して、以下の線形、対数線形、線形ロジスティックモデルのパラメータで表すこともできる。$$\begin{eqnarray}pr[Y_{a_0} = 1] & = & \psi_0 + \psi_1a_0 \ \tag{1}\label{1}\\ \log{pr[Y_{a_0} = 1]} & = & \theta_0 + \theta_1a_0\ \tag{2}\label{2}\\ logit\ {pr[Y_{a_0} = 1]} & = & \beta_0 + \beta_1a_0\ \tag{3}\label{3}\end{eqnarray}$$
ここで、もし\(a_0 = 1\)ならば\(Y_{a_0} = Y_{a_0 = 1}\)であり、もし\(a_0 = 0\)ならば\(Y_{a_0} = Y_{a_0 = 0}\)である。特に、因果\(\text{RD} = \psi_1\)であり、因果\(\text{RR} = e^{\theta_1}\)であり、あるいは因果\(\text{OR} = e^{\beta_1}\)である。モデル\eqref{1}, \eqref{2}, \eqref{3}は飽和MSMである。それらは周辺モデルであり、共同分布ではなく、反事実確率変数\(Y_{a_0 = 1}\)と\(Y_{a_0 = 0}\)の周辺分布をモデルするからである（つまり、モデル\eqref{1}, \eqref{2}, \eqref{3}は、\(Y_{a_0 = 1}\)と\(Y_{a_0 = 0}\)の相関をモデルしない）。これらは構造モデルであるが、それは反実仮想変数の確率をモデル化するからであり、計量経済学や社会科学の文献では反実仮想変数のモデルはしばしば構造モデルと呼ばれる。最後に、それぞれが2つの未知のパラメータを持ち、したがって各モデルは2つの未知の確率\(pr(Y_{a_0 = 1} = 1)\)と\(pr(Y_{a_0 = 0} = 1)\)の取りうる可能な値に制限を設けないので、これらは飽和モデルである。これらのモデルは共変量を含まないことに注意する。なぜなら、これらのモデルは定義上、元集団全体に対する因果効果のモデルであり、観察された関連性のモデルではないからである。
粗RD, RR, ORは，観察された結果\(Y\)に関する次の飽和線形モデル、対数線形モデル、線形ロジスティックモデルのパラメータで表すこともできる。$$\begin{eqnarray}pr[Y = 1 | A_0 = a_0] & = & {\psi_0}^{\prime} + {\psi_1}^{\prime}a_0 \ \tag{4}\label{4}\\ \log{pr[Y = 1 | A_0 = a_0]} & = & {\theta_0}^{\prime} + {\theta_1}^{\prime}a_0 \ \tag{5}\label{5}\\ logit\ pr[Y = 1|A_0 = a_0] & = & \beta_{0}^{\prime} + {\beta_1}^{\prime}a_0\ \tag{6}\label{6}\end{eqnarray}$$
これらは、元集団の部分集団（治療のレベルによって定義される）を比較するときに観察される関連に関するモデルである。粗RDは\({\psi_1}^{\prime}\)に等しく、粗RRは\(e^{{\theta_1}^{\prime}}\)に等しく、粗ORは\(e^{{\beta_1}^{\prime}}\)に等しい。関連モデル\eqref{4}, \eqref{5}, \eqref{6}のパラメータは、治療が無交絡の場合を除き、MSM\eqref{1}, \eqref{2}, \eqref{3}のパラメータとは異なる。モデル\eqref{4}, \eqref{5}, \eqref{6}は観察されたデータに対するモデルなので、モデル・パラメータの（漸近的に）不偏推定値は標準的な統計ソフトウエアを用いて得ることができる（選択バイアスや測定誤差がないと仮定して）。治療が無交絡である場合、これらの同じ推定値は、モデル\eqref{1}, \eqref{2}, \eqref{3}の対応する因果パラメータについても不偏となる。たとえば、モデル\eqref{4}, \eqref{5}, \eqref{6}を適合させるために、結果\(Y\)を二値変数として指定したモデル文\(Y = A_0\)を用いて、汎用のSASプログラムProc Genmodを使うことができる。\({\psi_1}^{\prime}\)を推定するには、identityリンクを指定し、\({\theta_1}^{\prime}\)を推定するには、logリンクを指定し、\({\beta_1}^{\prime}\)については、logitリンクを指定する。Proc Genmodに類似したプログラムは、S-Plus、Gauss、Stataなどの他のパッケージにもある。
4. 未測定の交絡なし
ここで治療が交絡しているとする。すると、粗関連パラメータは対応する因果パラメータと等しくならない。同様に、MSMのパラメータは、対応する観測データモデルのパラメータと等しくならない（例えば、\({\beta_0}\ne {\beta_0}^{\prime}, {\beta_1}\ne {\beta_1}^{\prime}\)）。測定された交絡因子\(L_0\)に関するデータが与えられ、測定されていない交絡因子がないと仮定すると、Proc Genmodを使用して重み付け分析を実行することにより、原因パラメータ\(\psi_1, \theta_1, \beta_1\)の不偏推定を得ることができる。具体的には、Proc Genmodのweight文（つまり、オプションSCWGT）を用いて、各被験者\(i\)に、自身の治療を受ける条件付き確率の逆数に等しい重み\(w_i\)を割り当てる。すなわち\(\displaystyle w_i = \frac{1}{pr[A_0 = a_{0i} | L_0 = l_{0i}]}\)であり、たとえば\(l_{0i}\)は被験者\(i\)の変数\(L_0\)の観察値である。真の重み\(w_i\)は未知であるが\(L_0\)上の\(A_0\)の予備ロジスティック回帰のデータから推定できる。たとえば、\(A_0\)がAZT治療で、\(L_0\)が、年齢、CD4カウント、WBCカウント、ヘマトクリット、および症状の有無を成分とする共変量の列ベクトルで、\(\alpha_1\)が未知のパラメータの行ベクトルであるロジスティック回帰モデルを指定することができる。$$logit\ pr[A_0 = 1|L_0 = l_0] = {\alpha}_0 + {\alpha}_1l_0\ \tag{7}\label{7}$$
標準的なロジスティック回帰ソフトウェアを用いて\(({\alpha}_0, {\alpha}_1)\)の推定値\((\hat{\alpha_0}, \hat{\alpha_1})\)を得ることができる。そして、\(A_0 = 0\)で\(L_0 = l_{0i}\)の被験者\(i\)について、\(\displaystyle w_i = \frac{1}{pr[A_0 = 0|L_0 = l_{0i}]}\)を\(\displaystyle 1/\frac{1}{1+\exp{(\hat{\alpha_0}+\hat{\alpha_1}l_{0i})}} = 1 + \exp{(\hat{\alpha_0}+\hat{\alpha_1}l_{0i})}\)によって推定する。\(A_0 = 1\)で\(L_0 = l_{0i}\)の被験者\(i\)については、\(\displaystyle w_i = \frac{1}{pr[A_0 = 1|L_0 = l_{0i}]}\)を\(\displaystyle \frac{1+\exp{(\hat{\alpha_0}+\hat{\alpha_1}l_{0i})}}{\exp{(\hat{\alpha_0}+\hat{\alpha_1}l_{0i}})} = 1+\exp{(-\hat{\alpha_0}-\hat{\alpha_1}l_{0i})}\)によって推定する。
要約すると、\(L_0\)に関するデータが与えられ、測定されていない交絡因子がない場合、各被験者\(i\)を\(w_i\)で重み付けすることによって粗分析を修正することで、（\(L_0\)による）交絡をコントロールすることができる。\(w_i\)の分母は、非公式には、被験者が自身の観察された治療を受けた確率である。したがって、これらの重み付き推定量をIPTW推定量と呼ぶ。
なぜこの方法が有効なのか？Proc Genmodにおける重み付けの効果は、各被験者\(i\)のコピー\(w_i\)からなる擬似集団を作ることである。つまり、ある被験者について\(w_i = 4\)とすると、その被験者は擬似集団に自分自身のコピーを4つ寄与する。この新しい擬似集団には、次の2つの重要な性質がある。第一に、擬似集団では、実際の集団とは異なり、\(A_0\)は測定された共分散\(L_0\)に束縛されない。第2に、擬似母集団における\(pr(Y_{a_0 = 1} = 1)\)と\(pr(Y_{a_0 = 0} = 1)\)は、真の調査母集団と同じなので、因果RD、RR、ORは両母集団で同じである。したがって、疑似集団における標準的な粗分析によって、原因RD、RR、ORを不偏に推定できることになる。しかし、重み\(w_i\)は必要に応じて各被験者のコピー\(w_i\)を作成するのに役立つので、これこそが我々のIPTW推定器が行うことなのである。Appendixでは、IPTWの手法をより明確にするために、詳細な数値例を示し、RosenbaumとRubinの密度スコア手法と比較する。
5. 未測定の交絡因子
測定されていない交絡因子\(U_0\)が存在する場合でも、Proc Genmodで解析を実施する際に、重み\(\displaystyle w_i = \frac{1}{pr[A_0 = a_{0i}| L_0 = l_{0i}, U_0 = u_{0i}]}\)を用いれば、上記のように原因リスク差、リスク比、オッズ比を不偏に推定することができる。しかし、\(U_0\)というデータは観測されていないので、これらの重みを不偏に推定することはできない。実際、測定されていない交絡因子が存在する場合、どのような方法でも不偏に推定することは、強力な付加的仮定なしでは不可能である。
6. マルチレベル処理と不飽和MSM
治療が図2のcで表されているように無交絡であり、\(A_0\)が100mgのAZTを単位とする被験者の1日の投与量を表す順序変数であるとする。\(A_0\)の取りうる値は、\(0, 1, \cdots, 14, 15\)である。この場合、各被験者に関連する潜在的な結果の数は\(16\)となる。具体的には、被験者が観察された投与量ではなく、投与量\(a_0\)を受けた場合に観察されたであろう\(Y\)の値を\(Y_{a_0}\)とする。したがって原理的には、被験者は\(16\)の可能なAZT投与量\(a_0\)ごとに別々の反事実変数を持つことになる。被験者の観察結果\(Y\)は、被験者の観察用量と等しい用量\(a_0\)に対応する結果\(Y_{a_0}\)である。説明の便宜上、観察された用量に等しい\(a_0\)では\(Y_{a_0}\)が事実変数\(Y\)であるにもかかわらず、\(Y_{a_0}\)のすべてを反事実と呼び続ける。むしろ、我々は通常、$$logit\ pr[Y_{a_0} = 1] = \beta_0 + \beta_1a_0\ \tag{8}\label{8}$$のような線形ロジスティックMSMを指定することによって、簡略化された用量反応関係を仮定する。 このモデルは、すべての被験者が用量\(a_0\)で治療されていた場合の成功確率は、傾きパラメータ\(\beta_1\)と切片\(\beta_0\)を持つ用量の線形ロジスティック関数であり、\(e^{\beta_1}\)はAZT用量100mgの増加に関連する因果ORであると言う。
観察されたデータについて、MSMモデル\eqref{8}を次の線形ロジスティック関連モデルと対比する。$$logit\ pr[Y = 1|A_0 = a_0] = {\beta_0}^{\prime} + {\beta_1}^{\prime}a_0\ \tag{9}\label{9}$$選択バイアスや測定誤差がないと仮定すると、SAS Proc Logistic や Proc Genmodのような標準的なロジスティック回帰ソフトウェア・パッケージを用いて線形ロジスティックモデル\eqref{9}をフィッティングすることにより、関連パラメータ\({\beta_0}^{\prime}\)と\({\beta_1}^{\prime}\)を不偏的に推定できる。もし治療が無交絡であれば、モデル\eqref{8}と\eqref{9}のパラメータは等しくなる。結果として、\({\beta_1}^{\prime}\)のロジスティック回帰推定は、我々の因果パラメータ\(\beta_1\)の不偏推定でもある。
もし治療が測定変数\(L_0\)によって交絡されるなら、\(\beta_1 \ne {\beta}^{\prime}\)で、\({\beta}^{\prime}\)の標準ロジスティック回帰推定は\(L_0\)による交絡のために因果パラメータ\(\beta_1\)の偏った推定となる。しかし治療が交絡される場合でも、\(L_0\)が与えられて測定されない交絡因子がない場合、被験者固有の重み\(\displaystyle w_i = \frac{1}{pr[A_0 = a_{0i} | L_0 = l_{0i}]}\)を用いれば，Proc Genmodでロジスティックモデル\eqref{9}をフィッティングすることによって、モデル\eqref{8}の因果パラメータ\(\beta_1\)の不偏推定を得ることができる。繰り返すが、実際には\(w_i\)は未知であり、モデルを指定してデータから推定しなければならない。
例えば、パラメータ\({\alpha}_{01}, {\alpha}_{02}, \cdots, {\alpha}_{015}\)および\(\alpha_1\)の推定値を得るために、Proc Logistic または Genmodを使用してSASで適合できる次のような多値ロジスティックモデルを指定することができる。$$\begin{eqnarray}pr[A_0 = a_0|L_0=l_0] & = & \frac{\exp{(\alpha_{0a_0}+\alpha_1l_0)}}{1+\sum_{j=1}^{15}{\exp{(\alpha_{0j}+\alpha_1l_0)}}}, \alpha_0 = 1, \cdots, 15;\\ pr[A_0=0|L_0=l_0] & = & \frac{1}{1+\sum_{j=1}^{15}{\exp{(\alpha_{0j}+\alpha_1l_0)}}} \ \tag{10}\label{10}\end{eqnarray}$$
6.1 安定したウエイト
\(L_0\)の成分が\(A_0\)と強く関連している場合、確率\(pr[A_0=a_{0i}|L_0 = l_{0i}]\)は被験者によって大きく異なる可能性がある。このばらつきにより、少数の被験者では重み\(w_i\)の値が非常に大きくなることがある。これらの少数の被験者は、自分自身の非常に多くのコピーを擬似集団に寄与するため、重み付け分析を支配することになり、その結果、IPTW推定量は大きな分散を持ち、ほぼ正規分布にならない。MSMモデルが飽和している場合（例えば、モデル\eqref{1}, \eqref{2}, \eqref{3}）、この変動は避けられないが、それは交絡因子\(L_0\)が治療\(A_0\)と高い相関を持つ結果、データの情報不足を反映しているからである。しかし、モデル\eqref{8}のような不飽和MSMの場合、このばらつきは、重み\(w_i\)を「安定化重み」\(\displaystyle sw_i = \frac{pr[A_0 = a_{0i}]}{pr[A_0 = a_{0i}|L_0 = L_{0i}]}\)で置き換えることによって、かなりの程度まで緩和することができる。安定化重みを理解するために、\(A_0\)が未交絡で、\(A_0\)と\(L_0\)が関連せず、\(pr[A_0 = a_{0i}] = pr[A_0 = a_{0i}| L_0 = l_{0i}]\)であったと仮定する。そして、\(sw_i = 1\)となり、各被験者は同じ重みを寄与する。\(A_0\)が交絡している場合、\(sw_i\)は一定ではなく、被験者の\(L_0\)の値によって、数\(1\)の周りで変化する。さらに、Robinsは、Proc Genmodで重み\(w_i\)ではなく重み\(sw_i\)を使用する場合、MSMのパラメータ\(\beta\)の推定値は不偏性を保ち、不飽和MSMの場合、一般に変動が少なくなることを示している。飽和MSMの場合、\(\beta\)の推定値のばらつきは、安定化された重みを使っても安定化されていない重みを使っても同じである。
もちろん、\(pr[A_0 = a_{0i}]\)および\(pr[A_0 = a_{0i}| L_0 = l_{0i}]\)は未知であり、推定されなければならない。\(pr[A_0 = a_{0i}| L_0 = l_{0i}]\)は、上記のように推定することができ、\(pr[A_0 = a_{0i}]\)は、\(A_0\)が\(a_{0i}\)に等しい調査標本中の被験者の割合として推定することができる。この推定は、多値モデルをフィッティングすることによって得られる推定と等価であり、パラメータ\(\alpha^{*}_{0a_{0}}\)にアスタリスクを置き、\(A_0\)が交絡している場合、このパラメータがモデル\eqref{10}のパラメータ\(\alpha_{0a_{0}}\)と異なることを示す。$$\begin{eqnarray}pr(A_0 = a_0) & = & \frac{\exp{(\alpha^{*}_{0a_0})}}{1+\sum_{j=1}^{15}{\exp{(\alpha^{*}_{0j})}}}, \alpha_0 = 1, \cdots, 15;\\ pr(A_0=0) & = & \frac{1}{1+\sum_{j=1}^{15}{\exp{(\alpha^{*}_{0j})}}} \ \tag{11}\label{11}\end{eqnarray}$$
文献2では、Robinsが拡張IPTW推定量を紹介している。これらの推定量は、安定化重みを使用するIPTW推定量よりもさらに効率的であるが、計算がより困難である。
6.2 連続治療
被験者の1日のAZT摂取量を10分の1ミリグラム単位で測定することができたとすると、AZTは本質的に連続的な治療となる。さらに、説明の便宜上、AZT投与量が\(0\)に近い被験者がおらず、\(A_0\)の分布を正規分布として効果的にモデル化できると仮定する。今、各個人は、非常に多数の反事実結果\(Y_{a0}\)を持つ。もし安定化重み\(\displaystyle sw_i = \frac{f(a_{0i})}{f(a_{0i}|l_{0i})}\)を使用するなら、モデル\eqref{8}の因果パラメータ\(\displaystyle \beta_1\)の不偏推定を、Proc Genmodでロジスティック・モデル\eqref{9}をフィッティングすることによって得ることができる。\(f(a_0|l_0)\)を推定するには、\(L_0\)が与えられたとき、\(A_0\)は平均\(\alpha_0+\alpha_1L_0\)、分散\({\zeta}^2\)を持つ正規分布であると指定することができる。そして、\((\hat{\alpha}_{0}, \hat{\alpha}_{1}, \hat{\sigma}^2)\)の不偏推定\((\alpha_0, \alpha_1, {\sigma}^2)\)は、たとえば、SASのProc REGを用いた\(L_0\)での\(A_0\)の通常の最小2乗回帰によって得ることができる。そして\(f(a_{0i}|l_{0i})\)は、正規密度\(\exp{\{-[a_{0i}-(\hat{\alpha}_0+\hat{\alpha}_1l_{0i})]^2/2\hat{\sigma}^2\}}\)によって推定される。安定化重み\(sw_i\)(stabilized weight)の分子\(f(a_{0ui})\)を推定するために、\(A_0\)が平均\(\alpha^{*}_{0}\)と分散\({\sigma}^{*2}\)を持つ正規分布であることを指定することができる。\(f(a_{0i})\)は正規密度\(\displaystyle \frac{\exp{[-(\alpha_{0i}-\hat{\alpha}^{*}_{0})^2/2\hat{\sigma}^{*2}]}}{\sqrt{2\pi\hat{\sigma}^{*2}}}\)（\(\hat{\alpha}^{*}_{0}\)は観測された\(A_0\)の平均、\(\hat{\sigma}^{*2}\)はそれらの経験的分散）によって推定することができる。\(A_0\)が連続的であるとき、安定化されていない重み\(\displaystyle w_i = \frac{1}{f_{a_{0i}|l_{0i}}}\)に基づく推定は無限分散を持つので、使用できない。
6-3. 信頼区間
上述したように、安定化重み\(sw_i\)の推定値を用いて Proc Genmod で関連モデル をフィッティングすることにより、MSM\eqref{8}のパラメータを推定する。オプション “repeated “を選択し、独立作業相関行列を指定すると、Proc Genmod プログラムは、\(\displaystyle \hat{\beta}_1\pm1.96\sqrt{var(\hat{\beta}_1)}\)で与えられる\(\beta_1\)の\(95\%\)”ロバスト “Wald 信頼区間も出力する。ここで\(var(\hat{\beta}_1)\)は、\(\hat{\beta}_1\)の分散のいわゆる “ロバスト “または “サンドイッチ “推定量である。 Robinsは、”ロバスト”Wald 区間は、少なくとも95%のカバレッジ確率を持つことを示しているが、いくつかの追加プログラミングでより狭い有効区間を得ることができる。ほとんどの重み付きロジスティック回帰プログラムによって出力される通常の非ロバスト・モデルベースのWald 信頼区間は、少なくとも95%のカバレッジを提供することは保証されないので、避けるべきである。S-Plus、Gauss、STATAなどの他のソフトウェア・パッケージも “ロバスト “分散推定量を提供し、Proc Genmodの代わりに使用できる。
7. 時間依存性の治療
ここで第1節の設定に戻り、\(A_k\)は追跡開始から\(k\)日目のAZT投与量、\(Y\)は\(K+1\)日目の追跡終了時に測定された\(0/1\)の転帰変数とする。同様に、\(L_k\)は、転帰に関するすべての測定された危険因子のベクトルの\(k\)日目の値を表す。\(\bar{A_k} = (A_0, A_1, \cdots, A_k)\)を\(k\)日目までの治療または曝露履歴とし、\(\bar{A} = \bar{A_k}\)とする。\(\bar{L_k}\)と\(\bar{L}\)も同様に定義する。\(Y_{\bar{a}}\)を、すべての被験者が観察された投与歴\(\bar{A}\)ではなく、投与歴\(\bar{a} = (a_0, a_1, \cdots, a_K)\)を受けていたら観察されたであろう\(Y\)の値とする。各日\(k\)で\(a_k\)が二値であっても（つまり、各日で被験者は治療を受けているか受けていないかである）、\(2^{K}\)個の投与歴\(\bar{a}\)が存在し、したがって\(2^{K}\)個の可能な反事実が存在するが、そのうちの1つだけが観察される（つまり、事実である）ことに注意する。したがって、飽和MSMを推定することはできないかもしれない。一般的には、\(\displaystyle \text{cum}(\bar{a}) = \sum_{k = 0}^{K}{a_k}\)が投与歴\(\bar{a}\)に関連する追跡終了までの累積投与量である以下のような線形ロジスティックMSMを指定することにより、ある種の簡略化された用量反応関係を仮定する。$$logit\ {pr[Y_{\bar{a} = 1}]} = \beta_{0} + \beta_1\text{cum}(\bar{a}) \ \tag{12}\label{12}$$
$$logit\ pr[Y= 1 | \bar{A} = \bar{a}] = \beta^{\prime}_0 + \beta^{\prime}_1cum(\bar{a})\ \tag{13}\label{13}$$追跡不能選択バイアスや測定誤差がないと仮定すると、標準的なロジスティック回帰ソフトウェア・パッケージを用いて線形ロジスティックモデル\eqref{13}をフィッティングすることにより、パラメータ\(\beta^{\prime}_{1}\)を不偏的に推定することができる。治療が無交絡であれば、モデル\eqref{12}と\eqref{13}のパラメータは等しい。結果として、\(\beta^{\prime}_{1}\)のロジスティック回帰推定は、我々の因果パラメータ\(\beta_1\)の不偏推定でもある。もし治療が交絡しているなら\(\beta_1\ne \beta^{\prime}_1\)で、\(\beta^{\prime}_1\)の標準ロジスティック回帰推定は\(\bar{L}_k\)による交絡の結果として、因果パラメータ\(\beta_1\)の偏った推定となる。しかし、治療が交絡されるとき\(L_k\)が与えられて測定されない交絡因子がない場合、 \(\displaystyle \prod_{k = 0}^{K}{b_k} = b_0\times b_1\times b_2\times \cdots \times b_K\)および\(\bar{A}_{-1}\)が\(0\)と定義される安定化重みでロジスティックモデル\eqref{13}を適合させることによって、モデル\eqref{12}の因果パラメータ\(\beta_1\)の不偏推定を得ることができる。\(K = 0\)の特別なケース（つまり、点治療研究）では、モデル\eqref{12}と\eqref{13}は、我々の以前のモデル\eqref{8}と\eqref{9}に、\eqref{14}は、我々の以前の\(sw_i\)に帰着する。$$\begin{eqnarray}sw_i & = & \prod_{k = 0}^{K}{pr(A_k = a_{ki}| \bar{A}_{k-1} = \bar{a}_{(k-1)i})}/\left\{\prod_{k = 0}^{K}{pr(A_k = a_{ki}| \bar{A}_{k-1} = \bar{a}_{(k-1)i}, \bar{L}_k = \bar{l}_{ki}})\right\}\ \tag{14}\label{14}\end{eqnarray}$$\(sw_i\)の分母は、非公式にはある被験者が時間\(K\)を通して自分自身の観察された治療履歴を持つ条件付き確率である。時間に依存する治療では、安定化されていない重みの変動はしばしば莫大になり、その結果\(\beta\)の推定量は著しく非正規なサンプリング分布で大きく変動する可能性がある。したがって、安定化重みの使用を強く推奨する。
治療が交絡する場合、政策的に重要なのは、関連モデル\eqref{13}のパラメータ\(\beta^{\prime}_1\)ではなく、MSM\eqref{12}のパラメータ\(\beta_1\)であることを強調する。その理由を知るために、\(N\)人の研究対象者と交換可能なソース集団またはターゲット集団からの新しい対象者を考えよう。この被験者には、追跡調査終了時に血清からHIVが検出される確率、すなわち\(pr(Y_{\bar{a}} = 1)\)を最小化する治療\(\bar{a}\)を行いたい。したがって、例えばMSM\eqref{12}のパラメータ\(\beta_1\)が正の場合（すなわち、血液中にHIVが存在する確率は、AZT治療の期間が長くなるにつれて増加する）、私たちは、被験者からAZT治療を差し控えるであろう。対照的に、モデル\eqref{13}のパラメータ\(\beta^{\prime}_1\)は、共変量と治療との関連によって交絡する可能性がある。例えば、予後因子の履歴（例えば、CD4数）によって示されるように、予後が悪い被験者に対して医師が優先的にAZTを開始したとする。さらに、AZTが\(Y\)に因果的効果を及ぼさなかったとする（つまり、\(\beta_1 = 0\)）。そして、パラメータ\(\beta^{\prime}_1\)（したがって、重みなしロジスティック回帰からの我々の推定値）は、正になるが、\(Y\)でのAZTの効果としての因果的な解釈を持たない。
7.1 治療の影響を受ける変数でコントロールすることによって生じるバイアス
測定された共変量による交絡をコントロールする代替的なアプローチは、以下の式のような交絡履歴\(\bar{L}\equiv \bar{L}_k\)を調整する加重なしロジスティックモデルであると考えられる。$$logit\ pr[Y = 1 | \bar{ A} = \bar{a}, \bar{L} = \bar{l}] = \beta^{\prime\prime}_{0} + \beta^{\prime\prime}_{1}cum(\bar{a}) + {\beta}^{\prime\prime}_2cum(\bar{l}) + {\beta}^{\prime\prime}_3l_{k} + {\beta}^{\prime\prime}_{4}l_{k-1} + {\beta}^{\prime\prime}_{5}l_0$$ここでは簡単のため、\(L_k\)は時刻\(k\)における単一の共変量CD4カウントから構成されていると仮定する。さらに悪いことに、たとえ\(pr[Y = 1| \bar{A} = \bar{a}, \bar{L} = \bar{l}]\)のモデルが正しく規定されていたとしても、パラメータ\(\beta^{\prime\prime}_1\)は一般的に因果関係の解釈を持たない。これは、\(cum(\bar{A})\)が\(A_0\)を含む被験者の治療履歴全体に依存し、\(A_0\)が時間に依存する共変量\(L_k\)と\(L_{k-1}\)に影響する可能性があるためである。治療によって影響を受け、結果のリスク因子でもある共変量を調整するロジスティックモデルをあてはめると、関連パラメータ\(\beta^{\prime\prime}_1\)の不偏推定値は得られるが、因果パラメータ\(\beta_1\)の推定値は偏る。これは、図1のように\(L_k\)の成分（例えば、赤血球数）とアウトカム\(Y\)に測定されていない共通の原因\(U_0\)（例えば、ベースラインの骨髄幹細胞数）がある場合、直接的・間接的、または全体的な治療効果がないという帰無仮説（モデル\eqref{12}の\(\beta_1\)が\(0\)に等しくなる）の下でも当てはまる。
要約すると、標準的な回帰手法は共変量を回帰変数としてモデルに含めることによって調整する。なぜなら(1) \(L_k\)は後の治療の交絡因子かもしれないので調整されなければならないが、(2) 以前の治療にも影響されるかもしれないので、標準的な手法では調整されるべきではないからである。この難問の解決策は、\(L_k\)を回帰変数として回帰モデルに追加するのではなく、重み\(sw_i\)を計算するためにそれらを使用することによって、時間に依存する共変量\(L_k\)を調整することである。
8. ウエイトの推定
ここで、重み\(sw_i\)の推定方法を説明する。簡単のために、各時刻\(k\)における治療\(A_k\)は二値であると仮定する。まずモデル\eqref{14}の分母を考えよう。各人日を観測値として扱うプールロジスティックモデルを用いて、未知の確率\(pr[A_k = 1 | \bar{A_k} = 1 | \bar{A_{k-1}} = \bar{a_{k-1}}, \bar{L_k} = \bar{l_k}] \)を推定することから始める。たとえば、\(l_k\)がCD4カウント、WBC、ヘマトクリット、および時間\(k\)での症状のインジケータのベクトルであり、\(\alpha_4, \alpha_5, \alpha_6\)および\(\alpha_7\)が行ベクトルである以下のモデルをあてはめることができるかも知れない。$$logit\ pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{k-1}}, \bar{L_k} = \bar{l_k}] =\alpha_0 +\alpha_1k + \alpha_2a_{k-1} + \alpha_3a_{k-2}+\alpha_4l_k + \alpha_5l_{k-1}+\alpha_6a_{k-1}l_k + \alpha_7l_0 \ \tag{15}\label{15}$$このモデルは、\(k\)日目に治療される確率が、\(k\)日目、前の\(2\)日間の治療、現在と前の日の共変量、昨日の治療と今日の共変量の間の交互作用、およびベースラインの共変量に線形ロジスティックに依存することを言う。
標準的なロジスティック回帰プログラムを用いて、モデル\eqref{15}を適合させることができる。分子確率も同様に推定できるが、モデル\eqref{15}を適合させる際、共変量の関数である最後の4項を除去することを除いては、以下のように適合させる。すなわち，我々は次のようにモデルを適合させる。$$logit\ pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{k-1}}] = \alpha^{*}_{0} + \alpha^{*}_{1}k + \alpha^{*}_{2}a_{k-1}+\alpha^{*}_{3}a_{k-2}\ \tag{16}\label{16}$$
そして、各被験者\(i\)について、我々はロジスティックプログラムに、モデル\eqref{15}の適合からの推定予測値\(\hat{p_{0i}}, \cdots, \hat{p_{Ki}}\)を出力させ、これは\(pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{(k-1)i}}, \bar{L_k} = \bar{l_{ki}}]\)の最尤推定値である。同様に、モデル\eqref{16}から予測値\(\hat{p^{*}_{1i}}, \cdots, \hat{p^{*}_{Ki}}\)を出力し、これは量\(pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{k-1}i}]\)の推定値である。例えば、\(1-\hat{p}^{*}_{ki}\)は、\(a_{ki} = 0\)のときの確率\(pr[A_k = 1|\bar{A_{k-1}} = \bar{a_{(k-1)i}}]\)の推定値である。データ分析者は、\eqref{15}や\eqref{16}のようなモデルの適合から出力された予測値から、各被験者の\(sw_i\)を計算する小さなプログラムを書く必要がある。$$sw_i = \frac{\prod_{k = 0}^{K}{(\hat{p_{ki}^{*}})^{a_{ki}}(1-\hat{p_{ki}^{*}})^{1-a_{ki}}}}{\left\{\prod_{k=0}^{K}{(\hat{p_{ki}})^{a_{ki}}(1-\hat{p_{ki}})^{1-a_{ki}}}\right\}}\ \tag{17}\label{17}$$
未測定の交絡因子がないという我々の仮定の下では、\(pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{k-1}}, \bar{L_k} = \bar{l_k}]\)のモデル\eqref{15}が正しく規定されていれば、結果として得られる因果パラメータ\(\beta_1\)の推定値は不偏である。さらに、これらの同じ条件下で、95%ロバストWald 信頼区間は、少なくとも95%の確率で\(\beta_1\)をカバーすることが保証される。\(pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{k-1}}]\)のモデル\eqref{16}が誤って仕様化されていても、\(\beta_1\)の推定値は不偏のままである。実際、モデル\eqref{15}が正しく、治療が交絡している場合、ロジスティックモデルの非可換性のため、モデル\eqref{16}はいくらか誤仕様化されることが保証される。
9. 治療前の共変量による効果の修正
MSMは、治療前共変量を含めるように一般化できる。たとえば、モデル\eqref{12}は、\(V\)が測定された治療前共変量\(L_0\)のベクトルの成分であり、\(\beta_3\)が治療-共変量相互作用を表す以下の式に一般化できる。$$logit\ pr[Y_{\bar{a}} = 1 | V = v] = \beta_0 + \beta_1 cum(\bar{a}) + \beta_2v + \beta_3 cum(\bar{a})v \ \tag{18}\label{18}$$モデル\eqref{18}では、\(\beta_1 + \beta_3 v\)は、ベースライン変数\(V\)のレベル\(v\)内の線形ロジスティック尺度での累積治療の効果を表すことに注意する。我々のIPTW推定量は、すでに\(V\)による交絡を自動的に調整するので、研究者がモデル\eqref{18}に含めることを選択する\(L_0\)の特定のサブセット\(V\)は、研究者の実質的な関心を反映するだけでよい。例えば、ある変数\(V\)は、研究者が\(V\)が効果修飾因子であると考え、かつその共変量\(V\)のレベル内における治療の因果効果について、原集団全体よりも実質的な関心が高い場合にのみ、モデル\eqref{18}に含めるべきである。
未測定の交絡因子がないという仮定の下で、\eqref{17}の推定重み\(sw_i\)を用いてProc Genmodを用いて$$logit\ pr[Y = 1 | \bar{A} = \bar{a}, V = v] = \beta_0 + \beta_1cum(\bar{a}) + \beta_2v + \beta_3cum(\bar{a})v\ \tag{19}\label{19}$$のような関連モデルをフィットさせ、\(p^{*}_k\)が\eqref{19}のようなモデルのフィットから推定された予測値であるという点だけを修正して、モデル\eqref{18}のパラメータの不偏推定値を得る。$$logit pr[A_k = 1 | \bar{A_{k-1}} = \bar{a_{k-1}}, V = v] = \alpha_{0}^{*} + \alpha_{1}^{*}k + \alpha_{2}^{*}a_{k-1} + \alpha_{3}^{*}a_{k-2} + \alpha_{4}^{*}v$$ 初歩的な疫学の教科書は、効果修飾が交絡とは論理的に異なることを強調している。それにもかかわらず、多くの学生はこの区別を理解するのが難しい。なぜなら、交絡因子のコントロールと効果修飾の検出の両方に同じ統計的手法（層別化と回帰調整）が使われているからである。したがって、初歩的な疫学的手法を限界構造モデルを用いて教えることには、何らかの利点があると思われる。なぜなら、その場合、交絡因子制御のための手法（治療の逆確率重み付け）は、効果修飾の検出のための手法（MSMに治療の共変量相互作用項を加える）とは区別されるからである。
最後に重要な注意点として、MSMは治療と時変共変量との交互作用のモデル化には使用できない。そのためには、構造的入れ子モデルを用いるべきである。したがって、モデル\eqref{18}の共変量\(V\)に、任意の時間\(k > 0\)で測定された時間依存性共変量\(L_k\)の成分を含めることは無効である。
10. 追跡調査不能による打ち切り
これまでは、各被験者が\(K+1\)時点の追跡終了まで観察されると仮定してきた。本節では、追跡不能による打ち切りを許容する。具体的には、被験者が\(k\)日目までに追跡不能になった場合は\(C_k = 1\)とし、そうでない場合は\(C_k = 0\)とする。いったん追跡不能になった対象者は、その後追跡調査に再参加しないと仮定する。打ち切りを単なる時変治療として概念化すれば、打ち切りを説明するために新しいアイデアは必要ない。この観点からは、打ち切りを調整したいということは、我々の関心が、観察された打ち切り履歴に従うのではなく、事実に反して、すべての被験者が打ち切られずに残っていた場合の処置\(\bar{A}\)の因果効果にあるというだけのことである。我々の目的は、ロジスティックMSM\eqref{12}のパラメータ\(\beta_1\)を推定することに変わりはないが、今、\(Y_{\bar{a}}\)は、もし被験者が事実に反して、治療履歴\(\bar{a}\)に従い、一度も打ち切られなかった場合の被験者の結果を意味する。ここでも、治療と打ち切りの両方に測定されない交絡因子がなければ、そうすることができる。この考えを定式化するために、各時刻\(k\)で、\(L_k\)の直前と\(A_{k-1}\)の直後の図1のグラフに変数\(C_k\)を追加する。その場合、測定された共変量\(L_k\)は、交絡と追跡不能による選択バイアスの両方を調整するのに十分である。
ここでも、我々は、適切な重みを含む線形ロジスティック関連モデル\eqref{13}を適合させることによって、因果パラメータ\(\beta_1\)の不偏推定値を得ることができる。結\(\)Y果は被験者が脱落しない限り観察されない。つまり、\(\bar{C} = (C_0, \cdots, C_{K+1})\)なので、モデル\eqref{13}の重み付きロジスティック回帰のあてはめは、打ち切られていない被験者に制限される。必要な被験者固有の重みは\(sw_i = sw_i\times sw_{i}^{\dagger}\)で、ここで $$sw_{i}^{\dagger} = \prod_{k = 0}^{K+1}{pr(C_k = 0|\bar{C_{k-1}} = 0)}$$と、さらに\(sw_i\)の定義と推定では，モデル\eqref{14}-\eqref{16}の各条件づけイベントの右側にイベント\(\bar{C_k} = 0\)を追加する。\(sw_i^{\dagger}\)の未知の確率は、各人の日を1つの観測値として扱うプールされたロジスティックモデルを使って推定することができる。具体的には、logit pr[Ck = 0|¯Ck-1 = 0, Āk-1 = āk-1, ¯Lk = ¯lk]とlogit pr[Ck = 0|¯Ck-1 = 0, Āk-1 = āk-1]について、モデル\eqref{15}と\eqref{16}のアナログを適合させる。\(sw_i\times sw_{i}^{\dagger}\)の積の分母は、非公式には打ち切られていない被験者が、時間\(K+1\)まで観察された治療と打ち切りの履歴を持つ条件付き確率であることに注意する。したがって，我々は重み付きロジスティック推定量を逆-probability-of-treatment-and-censoring重み付き推定量と呼ぶ。もし\((A_k, C_k)\)を時間\(k\)における “共同治療 “と見なせば、この分母を単に被験者が彼または彼女自身の治療履歴をたどる確率として非公式に解釈することができる。
11. 周辺構造モデルの限界
文献2やAppendix2で示されているように、IPTW推定量にはバイアスがかかるため、各時刻\(k\)で共変量レベル\(l_k\)が存在し、そのレベルの共変量を持つ被験者全員が同一の治療\(a_k\)を受けることが確実であるような研究では、MSMを用いるべきではない。例えば、このような状況が想定される職業コホート研究ではMSMを用いるべきではないことを意味する。その理由を説明するために、時点\(k\)における工業化学物質への暴露レベルを\(A_k\)とし、時点\(k\)で被験者が休職していれば\(L_k = 1\)で、そうでなければ\(L_k = 0\)とする職業コホート研究を考えてみよう。その場合、\(L_k = 1\)の対象者はすべて\(A_k = 0\)となる。同様に、子宮頸がんによる死亡率に対するスクリーニングの効果に関する研究では、時刻\(k\)までに子宮頸部を手術で摘出した女性（これを\(L_k = 0\)とする）は、その時点で曝露（つまりスクリーニング）を受けることができないので、MSMを使用すべきではない。それにもかかわらず、構造的入れ子モデルのg推定は、MSMを使用できない研究であっても、暴露効果を推定するために常に使用できる。付随論文に記載されているMulticenter AIDS Cohort Studyのデータ解析のような多くの研究では、実質的な検討に基づいて、上記の困難は発生せず、MSMは実用的な方法であると私たちは信じている。
12. 結論
我々は、二値的な結果に対する時間的に変化する曝露や治療の因果効果を推定するために、MSMを使用する方法を説明した。我々のグループの論文では、結果を生存時間アウトカムに拡張し、MSMに基づく方法を、構造的入れ子モデルのg-推定やg-計算アルゴリズム式の推定に基づく、以前に提案された別の方法と比較対照している。
コメント