Propensity score analysis の sensitivity analysis

Propensity scoreをやってみたわけだが、観察研究を事後解析的になんやかんやするので、観察して得られたパラメータしか持っていない。
その限られたパラメータから、治療が行われる確率としてPSが計算されるので、もし、このモデルに組み込まれておらず、かつ、治療割付の確率を著しく変えてしまうような強力なパラメータが存在したとしたら、PSがほぼ同一の人をマッチングさせる、というこの解析の根本がガラガラと音を立てて崩れる。
というわけで、そのような強力な、未確認で進行形パラメータが存在した場合において、今回のPS マッチング解析がどれほどのものかを推定する方法がsensitivity analysisで、感度分析みたいな感じで訳される。
やり方としてはrboundsパッケージの解説文をなんとか読み解く感じで。
 
未確認パラメータ(unobserved confounder)がある、と仮定するわけだが、実際にこれが存在するとき…という状況を考えるのは面倒なので、unobserved confounderによって治療を行う確率が偏る場合を考える。つまり、PS調整を完全にうまくできたなら、同一のPSのペアを選んだときに、治療が行われる確率は1:1である、と言える。この状況であれば、「PSがRCTに近づいたぞー!!ヤター!!」となるわけだが、実際はunobserved confounderのせいで、PSマッチングをしたけれども、そのペアにおいて治療を行われる確率が微妙に異なるわけである。それをオッズとして考える。
 
治療が行われる確率の偏りを\Gammaとする。\Gamma=2ならば、あるPSで選ばれたけれども、治療を行われる確率が2倍高い、という意味になる。この\Gammaはよくわからない数値なので、この\Gammaを適当に増加させて分析するわけだが、social scienceの分野では6倍も偏る、ということは非常に考えにくいから、デフォルトでは6までとなっている。
PS マッチングしたあとで、ある人jさんに治療を行われる確率\pi_{j}は、jさんの共変量\mathbf{x}_{j}(これはいま手持ちの、すべてわかっているパラメータ)で決まる。で、もしunobserved confounderによるバイアスがあると、ある人jさんとある人kさんは、わかっているパラメータについては\mathbf{x}_{j}=\mathbf{x}_{k}だが、治療が行われる確率が異なるので\pi_{j}\not{=}\pi_{k}である。sensitivity analysisでは、\piがどれほど大きくなると、単純にPS マッチングした結果が統計学的有意でなくなってしまうかを考える。
\piは治療を行われる確率なので、治療を行われない確率とのオッズを考えると\frac{\pi}{1-\pi}であり、それをある人jさんとある人kさんで考えるとオッズ比になる。また、このオッズ比は、先の\Gammaを用いて、たかだか
\frac{1}{\Gamma}\leq\frac{\pi_{j}/(1-\pi_{j})}{\pi_{k}/(1-\pi_{k})}\leq \Gamma
の範囲に収まるらしい。
ここで、unobserved confounderをuとすると、PSの計算に用いたロジスティック回帰で、
log(\frac{\pi_{j}}{1-\pi_{j}})=\kappa(\mathbf{x})+\gamma u_{j}
となる。ただし\kappaは適当な関数、\gammaは適当な係数。
するとまあ、
\frac{\pi_{j}/(1-\pi_{j})}{\pi_{k}/(1-\pi_{k})}=\exp\{\gamma(u_{j}-u+k)\}
とかなって、\Gammaがunobserved confounderであるuの係数を対数処理した大きさで考えることができる。
 
要は偏りを考えましょう、ということで、こちらを見ながら、sensitivity analysisのexampleをやると、マッチした122ペアにおいて喫煙者と非喫煙者で肺癌死亡を調べたら、喫煙死亡(非喫煙者生存)は110ペア、非喫煙者死亡(喫煙者生存)は12ペアだった。
\begin{matrix}&&&death&\\&&smoke(+)&smoke(-)\\&smoke(+)&*&12\\alive&smoke(-)&110&*\end{matrix}
ここで*は、ただしイケメンに限るのではなく、喫煙者での生死と非喫煙者での生死を比較してもペアマッチングのうまみがないので、ここは計算に関係しない。
これをMcNemar検定すると

mcnemar.test(matrix(c(10, 110, 12, 10), 2))
	McNemar's Chi-squared test with continuity correction

data:  matrix(c(10, 110, 12, 10), 2)
McNemar's chi-squared = 77.123, df = 1, p-value < 2.2e-16

p値は1.6*10^{-18}となる。
ここで、\Gammaを用いて、喫煙者と非喫煙者を偏らせるようなunobserved confounderが存在するとしたら、としてsensitivity analysisを行う。
偏りによるp値のゆらぎを、上限T^{+}(upper bound)、下限T^{-}(lower bound)とすると、p^{+}=\frac{\Gamma}{1+\Gamma}を用いて
T^{+}=\sum_{a=110}^{122}\(\frac{122}{a}\)(p^{+})^{a}(1-p^{+})^{122-a}
となる。下限は興味がないので、\Gammaをいろいろ動かしたときに、upper boundであるT^{+}有意水準を満たせるかどうかが最終的に重要である。
で、この肺癌データでやると、\Gammaが5まではUpper boundが0.05以下なので、喫煙者に5倍偏るようなunobserved confounderがあっても統計学的有意は保てる、ということになる。

binarysens(12, 110)
 Rosenbaum Sensitivity Test 
 
Unconfounded estimate ....  0 

 Gamma Lower bound Upper bound
     1           0     0.00000
     2           0     0.00000
     3           0     0.00002
     4           0     0.00196
     5           0     0.02317
     6           0     0.09693

 Note: Gamma is Odds of Differential Assignment To
 Treatment Due to Unobserved Factors 

で、前の記事をやってみるわけだが、\Gamma=1.05でもうやばいっぽいからなんかこれunobserved confounderがありそうな気がしてならない。

binarysens(m.obj, Gamma=1.2, GammaInc=0.01) # 1.2 まで 0.01 刻みで
 Rosenbaum Sensitivity Test 
 
Unconfounded estimate ....  0.0139 

 Gamma Lower bound Upper bound
  1.00     0.01393     0.01393
  1.01     0.01007     0.01900
  1.02     0.00720     0.02549
  1.03     0.00510     0.03365
  1.04     0.00357     0.04373
  1.05     0.00248     0.05600
  1.06     0.00170     0.07067
  1.07     0.00116     0.08796
  1.08     0.00078     0.10800
  1.09     0.00052     0.13092
  1.10     0.00035     0.15673
  1.11     0.00023     0.18541
  1.12     0.00015     0.21683
  1.13     0.00010     0.25081
  1.14     0.00006     0.28707
  1.15     0.00004     0.32528
  1.16     0.00002     0.36507
  1.17     0.00002     0.40599
  1.18     0.00001     0.44760
  1.19     0.00001     0.48943
  1.20     0.00000     0.53101

 Note: Gamma is Odds of Differential Assignment To
 Treatment Due to Unobserved Factors