SASユーザー総会2022：Time-varying treatmentsに対するIPTW法による因果効果の推定

Time-varying treatmentsに対するIPTW法による 因果効果の推定 ○鈴木 徳太1、折原隼一郎2 1横浜市立大学 2東京医科大学 データサイエンス学部 医療データサイエンス分野 Estimation of causal effect on time-varying treatment using IPTW Norihiro Suzuki1, Shunichiro Orihara2 1School 2Department of Data Science, Yokohama City University of Health Data Science, Tokyo Medical University

要旨： 介入が複数の時点で発生し、かつその値が時間の経過とともに変わりうる場合 の因果効果の推定手法について概説する。特に周辺構造モデルにおける IPTW法に注目し、その統計的理論をSASでの実装方法と共に簡単に紹介する。 キーワード： causal inference, time-varying treatment, IPTW, MSMs 2

Agenda • Introduction • Sequential treatmentの導入 • IPTW of Marginal Structural Models • SASによる実装 3

統計的因果推論を取り巻く状況の変化 • 臨床試験で収集されたデータに加え、RWD（リアルワールドデータ）の利用が 活発化しており、近年ではRWDを承認申請へ活用することも期待されている • 状況に応じた適切な手法の選択が必要 – 1時点の治療に対する因果推論のprocedureは存在するものの、経時的な治療 には対応はしていない • • • – CAUSALMED (SAS/STAT 14.3, 2017~) CAUSALTRT (SAS/STAT 14.2, 2016~) PSMATCH (SAS/STAT 14.2, 2016~) Time-varying treatmentに対する因果推論を紹介 出典：日本製薬工業協会. (2021). リアルワールドデータ・疾患レジストリを活用した医薬品開発. https://www.jpma.or.jp/information/evaluation/results/allotment/bd_rwd_202105.html 4

• https://www.jpma.or.jp/information/evaluation/results/allotment/bd_rwd_202105.html

治療がベースラインでのみ行われる例 • 新型コロナワクチンの開発（第Ⅲ相臨床試験） – – 被験薬 vs プラセボ 多くの臨床試験ではランダム化や盲検化が実施 治療群 試験終了 被験薬 1 試験期間 対照群 プラセボ 注） 投与の不遵守や脱落等は考慮しない 5

治療が経時的に行われる例 • 新型コロナワクチンの実際の接種（ファイザー社製 or 武田/モデルナ社製） – – – 1, 2回目とそれ以降では、接種するワクチンが変わりうる (time-varying treatment) 複数のワクチンを接種した場合（交互接種）の効果にも注目 異なる治療計画 (treatment regime) の比較 個人A ファイザー ファイザー モデルナ モデルナ 1 2 3 4 モデルナ モデルナ ファイザー モデルナ 個人B 6

臨床試験で不遵守が存在するとき • 想定する群比較 (time-fixed treatment) 治療群 対照群 • 被験薬 対照薬 1 2 対照薬 対照薬 被験薬 3 被験薬 被験薬 4 被験薬 不遵守 一部の被験者が実際にとりうる治療計画 被験者A 被験者B – 被験薬 対照薬 1 2 対照薬 対照薬 対照薬 3 被験薬 Time-varying treatmentとみなすことが可能 対照薬 4 被験薬 7

Agenda • Introduction • Sequential treatmentの導入 • IPTW of Marginal Structural Models • SASによる実装 8

記法 • 𝐙! , … , 𝐙" : サイズがnのi.i.d.サンプル , $ , 𝐴̅ $ , Y ) – Z = ( 𝐗 # , 𝐴# , 𝐗! , 𝐴! , … , 𝐗 $ , 𝐴 $ ,Y ) = ( 𝐗 Ø𝐴% : 時点t における二値治療 (t = 0, … , T ) – 0:untreated , 1:treated – 𝐴̅! : 治療歴 ( 𝐴̅ " = 𝐴̅ ) Ø𝐗 % : 時点t における共変量 $ ! : 共変量歴 ( 𝐗 $" = 𝐗 $) –𝐗 Ø𝑌: 研究の最終時点において観測されるアウトカム – 𝑌 $# : 治療計画 𝐴̅ = 𝑎' を受ける場合の潜在アウトカム , % , 𝐴̅%&! ) : 時点t における治療までの共変量歴および治療歴 Ø𝐇% = (𝐗 9

想定するDAGと記号の対応 𝐗0 • • • • • 𝐴0 𝐗 # ：ベースライン共変量 𝐗! ：t = 1での共変量 𝐴# ：ベースラインでの治療 𝐴! ：t = 1での治療（２回目の治療） 𝑌：観測されるアウトカム 𝐗1 𝐴1 𝑌 𝐇! = (𝐗 # , 𝐗! , 𝐴# )：t = 1での治療以前の情報 10

治療が二時点で行われる場合の因果効果 𝐴0 𝐴1 𝑌 • 想定される治療計画 1. 2. 3. 4. • (𝐴%, 𝐴&) = (0, 0) ：ともに治療を受けない (𝐴%, 𝐴&) = (0, 1) ：t = 0では治療を受けず、t = 1では治療を受ける (𝐴%, 𝐴&) = (1, 0) ：t = 0では治療を受け、t = 1では治療を受けない (𝐴%, 𝐴&) = (1, 1) ：ともに治療を受ける それぞれの場合の潜在アウトカムの期待値の比較により平均因果効果を定義 Ex) 𝔼 𝑌 $!'&, $"'& − 𝔼 𝑌 $!'%, $"'% 11

因果効果の同定のために必要な条件 • 因果効果の同定のために必要となる仮定 1. 一致性 • If 𝐴̅ = 𝑎$ for a given subject, 𝑌 $# = 𝑌 for that subject Ø 観測されるアウトカムが潜在アウトカムと一致 2. 条件付き交換可能性 • 𝑌 $# ∐ 𝐴% | 𝐡% for all regime 𝑎$ and all 𝒙% Ø 未測定交絡が存在しないこと (No unmeasured confounders) と同義 3. 正値性 • ℙ 0 < 𝜋% 𝒉% < 1 = 1 for all 𝑎$ Ø 各時点において、治療を受ける確率が0 or 1になるような被験者は存在しない Ø 𝜋% 𝒉% = ℙ 𝐴% = 1 𝐇% = 𝐡% ) » 各時点において治療を受ける確率（傾向スコア） 12

時間依存性交絡の問題 𝐴0 • 𝐗1 𝐴1 𝑌 時間依存性交絡因子 (Time-dependent confounders) 𝐗! の存在 𝐴%と 𝑌の中間因子 (𝐴% → 𝐗& → 𝐴& → 𝑌) ※ 1. • 条件付けると因果効果の推定にバイアス 𝐴&と 𝑌の交絡因子 (𝐴& ← 𝐗& → 𝑌) 2. • 条件付けないと因果効果の推定にバイアス • 標準的な手法（回帰、層別、マッチング）では対応不可 – No feed back (𝐗& → 𝐴&がない)という仮定があれば可能だが、非現実的 ※または𝐴!と共通の原因を持つ 13

時間依存性交絡下での因果効果の推定手法 • Robins’ g-methods 1. 2. 3. g-computation algorithm formula (“g-formula”) IPTW of marginal structural models (MSMs) g-estimation of structural nested models (SNMs) 各手法の特徴（一部） メリット デメリット G-formula パラメトリックモデルが正しく特定 されれば推定精度が高い パラメトリックモデルを用いた際には Null paradoxが発生する可能性 IPTW of MSMs 直感的な理解、実装が容易 推定が不安定になりやすい G-estimation of SNMs 治療と時間依存性共変量の 交互作用を直接モデル化可能 （effect modificationの評価が容易） 時間依存性共変量による効果修飾が 存在する場合、non-dynamic regimeの 比較が難解 14

Agenda • Introduction • Sequential treatmentの導入 • IPTW of Marginal Structural Models • SASによる実装 15

周辺構造モデル① • 平均因果効果：𝜓 = 𝔼 𝑌 '", '# − 𝔼 𝑌 '")#, '#)# – average causal effect / marginal causal effect – 個人レベルの因果効果を、集団レベルに周辺化 (marginalize) しているため • 周辺構造モデル (MSMs: Marginal Structural Models) Ø 潜在アウトカムの周辺期待値と、治療計画𝑎の関係を示したもの ' • 𝔼 𝑌 $! , $" = 𝛽& + 𝜓& 𝑎& + 𝜓' 𝑎' – 最も単純なモデル • 𝔼 𝑌 $! , $" = 𝛽& + 𝜓& 𝑎& + 𝜓' 𝑎' + 𝜓( 𝑎& 𝑎' – それぞれの治療の交互作用項を含む 16

周辺構造モデル② 今回は周辺構造モデルを 𝔼 𝑌 (), (* = 𝛽* + 𝜓* 𝑎* + 𝜓+ 𝑎+ と考える • t = 0, 1でともに治療を受けた場合の因果効果 𝜓 – 𝜓 = 𝔼 𝑌 $!'&, $"'& − 𝔼 𝑌 $!'%, $"'% = ( 𝛽%+ 𝜓% + 𝜓& ) − 𝛽% = 𝜓% + 𝜓& – 𝜓%：t = 0（ベースライン）での治療による因果効果 – 𝜓&：t = 1での治療による因果効果 • 傾向スコアを用いた逆確率重み付け (IPTW) によるパラメータ推定 17

IPTW法の概要 • 逆確率重み付け (IPTW; Inverse Probability of Treatment Weighting) – 各個人のアウトカムをそれぞれが受けた実際の治療 𝐴̅ を受ける確率の逆数で重み付け – 各個人iへの重み𝑊+ 𝑊* = • 1 𝑓 𝐴*,# 𝒉*,# ) × 1 𝑓 𝐴*,! 𝐡*,! ) ! =; +)# 1 𝑓 𝐴*,+ 𝐡*,+ ) 重み付けした後の集団 = 疑似母集団 (pseudo population) – – 疑似母集団に対してモデルの当てはめ 𝔼 𝑌 $!, $" = 𝔼,- 𝑌|𝐴% = 𝑎%, 𝐴& = 𝑎& – 𝔼 𝑌 $!,'& $"'& − 𝔼 𝑌 $!'%, $"'% = 𝔼,- 𝑌|𝐴% = 1, 𝐴& = 1 −𝔼,- 𝑌|𝐴% = 0, 𝐴& = 0 18

IPTW法の直感的なイメージ • 一般には各治療計画をとった集団の共変量の分布は不揃い – 比較可能な集団ではない 𝐴" = 0 𝐴! = 0 𝐴" = 1 𝐴! = 1 𝐴" = 0 𝐴" = 1 (𝐴! , 𝐴" ) = (0, 0)の群 (𝐴! , 𝐴" ) = (0, 1)の群 (𝐴! , 𝐴" ) = (1, 0)の群 (𝐴! , 𝐴" ) = (1, 1)の群 19

IPTW法の直感的なイメージ • 一般には各治療計画をとった集団の共変量の分布は不揃い – 比較可能な集団ではない 人数 𝐴" = 0 (𝐴! , 𝐴" ) = (0, 0)の群 𝐴! = 0 被験者の状態 𝐴! = 1 𝐴" = 1 (𝐴! , 𝐴" ) = (1, 1)の群 人数 20

IPTW法の直感的なイメージ • 集団全員がある治療計画を受けた場合の仮想的な結果を生成 – 擬似母集団 (pseudo population)の作成 𝐴" = 0 𝐴! = 0 𝐴" = 1 𝐴! = 1 𝐴" = 0 𝐴" = 1 (𝐴! , 𝐴" ) = (0, 0)の群 (𝐴! , 𝐴" ) = (0, 1)の群 (𝐴! , 𝐴" ) = (1, 0)の群 𝐴! = 0 𝐴" = 0 𝐴! = 0 𝐴" = 1 𝐴! = 1 𝐴" = 0 𝐴! = 1 𝐴" = 1 (𝐴! , 𝐴" ) = (1, 1)の群 21

IPTW法の直感的なイメージ • 例えば次のような擬似母集団の比較では… 𝐴! = 0 𝐴" = 0 𝐴! = 1 𝐴" = 0 VS → 𝜓* （ベースラインでの治療による因果効果）がわかる！ • IPW推定量に帰着 – 具体的な数式展開、証明については参考文献を参照 22

Stabilized weight • 個人iへの重み𝑊* , 𝑊* = ; +)# 1 𝑓 𝐴*,+ 𝐡*,+ ) – 傾向スコアが0か1に近い場合には𝑊は非常に極端な値となり、推定が不安定になる – Stabilized weightの使用が推奨 • Stabilized weight (SW) , 𝑆𝑊* = ; +)# 𝑓 𝐴*,+ 𝐴̅+&! ) 𝑓 𝐴*,+ 𝐡*,+ ) – 分母と分子が相関するため極端な値を取りにくい • 想定するMSMにベースライン共変量𝑿# を組み込む場合（≒ ベースライン共変量による効果修飾に興味） にはSWの分子の条件付けに入れる必要 23

Agenda • Introduction • Sequential treatmentの導入 • IPTW of Marginal Structural Models • SASによる実装 24

想定するDAG 𝐗0 • • • • • 𝐴0 𝐗1 𝐴1 𝑌 𝐗 # ：ベースライン共変量 𝐗! ：t = 1での共変量（時間依存性交絡因子） 𝐴# ：ベースラインでの治療 𝐴! ：t = 1での治療（２回目の治療） 𝑌：観測されるアウトカム 25

想定するデータ ID Trt0 Trt1 Outcome Gender Age Press TG 1 1 1 518 1 49 100 62 2 1 0 410 1 50 153 107 3 1 0 472 0 54 98 58 ・ ・ 𝐴* ※ 1被験者1オブザベーションを想定 𝐴+ Y 𝐗* 𝐗+ 26

シミュレーションの設定 • 変数の生成部分（抜粋） Gender = rand('BINORM',0.5,1); Age = round(20+45*rand('UNIFORM'),1); if Gender=0 then Press = round(rand('NORMAL',100,20)+0.6*Age,1); if Gender=1 then Press = round(rand('NORMAL',110,20)+0.6*Age,1); TG = round(rand('NORMAL',90,15) +0.6*Age + 10*Gender,1); PS0_True=logistic(0.5+0.8*Gender-0.03*Age); Trt0 = rand('BINOMINAL',PS0_True,1); PS1_True=logistic(0.6-0.4*Gender-0.01*Age+0.5*Trt0+0.01*Press-0.015*TG); Trt1 = rand('BINOMINAL',PS1_True,1); epsilon = rand('NORMAL',0,10); Outcome = round(130 - 10*Trt0 - 4*Trt1 +Gender +Age +Press+ TG+ epsilon,1); • 推定したい真の因果効果（パラメータ） – 𝜓% = −10 – 𝜓& = −4 27

因果効果の推定の3ステップ ① 各時点の傾向スコアの算出 – 今回推定しなければならない傾向スコアは2つ • • 𝜋& 𝒉& = ℙ 𝐴& = 1 𝒙& ) 4' , 𝑎& ) 𝜋' 𝒉' = ℙ 𝐴' = 1 𝒙 – ロジスティック回帰モデルによる推定 & Ex) ℙ 𝐴% = 1 𝒙% ) = &./ ((*!+*",-.+*/0.12.3) ② 各個人への重みの計算 – 今回はSWを用いるため、各時点におけるℙ 𝐴+,0 = 𝑎 𝐴̅01&)も計算 ③ パラメータの推定 ①、②の結果を用いてMSMのパラメータ推定 28

因果効果の推定の3ステップ (1/3) ① 各時点の傾向スコアの算出 /*各時点の傾向スコアの算出*/ /*t=0*/ proc logistic data=df descending; class Gender(ref="0")/param=ref ref=first; model Trt0=Gender Age; output out=df predicted=PS0; run; /*t=1*/ proc logistic data=df descending; class Gender(ref="0")/param=ref ref=first; model Trt1=Gender Age Trt0 Press TG; output out=df predicted=PS1; run; 29

因果効果の推定の3ステップ (2/3) ② 各個人への重みの計算 /*SW用*/ proc logistic data=df descending; model Trt0=; output out=df predicted=Pr0; run; proc logistic data=df descending; class Trt0/param=ref ref=first; model Trt1=Trt0; output out=df predicted=Pr1; run; /*weightの計算*/ data df_wt; set df; if Trt0=1 then wt0=Pr0/PS0; else wt0=(1-Pr0)/(1-PS0); if Trt1=1 then wt1=Pr1/PS1; else wt1=(1-Pr1)/(1-PS1); wt=wt0*wt1; run; 30

因果効果の推定の3ステップ (3/3) ③ パラメータの推定 /*因果パラメータの推定*/ proc genmod data=df_wt; class ID Trt0 Trt1 /param=ref ref=first; model Outcome = Trt0 Trt1/ error=normal link=id; weight wt; repeated subject=ID/ type=ind; run; 31

推定結果 • 推定された因果効果（パラメータ） – 𝜓> # = −10.8851(−13.2644, −8.5057) – 𝜓>! = −4.0081(−6.3476, −1.6687) 32

IPTW of MSMsの限界（一部） • モデルの誤特定の可能性 – 周辺構造モデル、傾向スコアモデルが共に正しく特定される必要 • 推定精度、およびpositivityの問題 – – 正値性が守られない場合には因果効果の同定不可 (positivity violation) 時点が増えるほど極端な重みを取りやすい – Structural nested mean modelやpositivityを必要としない方法など • 感度解析を行う際に必要となる仮定が他の方法よりも多く、有用性に欠ける 33

おわりに • • 在学中は統計的因果推論の諸手法の実装方法の紹介や、CDISC標準に準拠 したデータセットに対応するmacroの開発・改良を行う予定です。実務上の観点 から、ご意見いただきたいと思っています。 上記以外にも研究課題等ありましたら、お気軽にご連絡ください！ • • 連絡先：[email protected] GitHub：https://github.com/Norihiro-Suzuki/SAS-Causal-inference – 今回紹介した実装例のコードの保管場所です – 他の手法の実装方法等についてもアップロード予定です 34

• https://github.com/Norihiro-Suzuki/SAS-Causal-inference

Reference • • • • • 日本製薬工業協会. (2021). リアルワールドデータ・疾患レジストリを活用した医薬品開発. Robins, J. M., & Hernán, M. A. (2009). Estimation of the causal effects of time-varying exposures. Longitudinal data analysis, 553-599. Robins, J. M., Hernan, M. A., & Brumback, B. (2000). Marginal structural models and causal inference in epidemiology. Epidemiology, 11(5), 550-560. Naimi, A. I., Cole, S. R., & Kennedy, E. H. (2017). An introduction to g methods. International journal of epidemiology, 46(2), 756-762. Kennedy, Edward H. (2019). "Nonparametric causal effects based on incremental propensity score interventions." Journal of the American Statistical Association, 114(526), 645-656. 35

Appendix 36

暗に仮定する条件 • Fixed study population – – • フォローアップ以後の治療（曝露）効果のみを考慮 – • Pre-baseline treatmentの因果効果を考える場合には、追加の仮定を置く必要 打ち切り、脱落、測定誤差は存在しない – • 全ての被験者に対しフォローアップが完全に行われる 十分に定義された閉鎖コホート (well defined closed cohort) 打ち切りが発生する場合には、IPTWだけではなく、IPCWも用いる必要 共変量、治療は一定間隔で測定される 37

治療計画の分類 • 治療計画のタイプは以下の2つの観点から分類することが可能 – Dynamic / Static • Dynamic regime – • Static regime – – 任意の時点tでの治療𝑎% が時間のみに依存 Random / Deterministic • Random regime – • 任意の時点において、いずれかの介入を受ける確率が確率的に変動する Deterministic regime – • 任意の時点tでの治療𝑎% が時間と他の変数𝒉% に依存 任意の時点において、いずれかの介入を受ける確率は100%（決定論的に定まる） なお今回の発表はstatic regimeを想定した 38