データサイエンティスト上がりのDX参謀・起業家

データサイエンティスト上がりのDX参謀・起業家のブログ。データ分析や事業について。自身はアーティスト、経営者、事業家。

操作変数(instrumental variable, IV)での交絡調整

傾向スコアでの調整は、未測定の交絡要因は調整できないという問題がありました。

そこで操作変数(IV)という方法が提案されており、この方法では未測定の交絡要因も調整できると言われています。

それによって、よりランダム化試験に近い結果が出ると期待されます。


操作変数の定義は次のようになります。

  • A variable that is related to treatment but neither directly nor indirectly related to outcome, except through the effect of the treatment itself
  • 治療には関連しているが、結果には治療を通してでしか直接的にも間接的にも関連していない変数
    • 治療群Xを予測できる変数Zは、Xを通してでしか結果Yに関連しない(未知の交絡要因Uも介さない)


Rで実行するには、パッケージSEMを使います。

IVは計量経済(econometrics)の分野でも使われており、構造方程式モデルを利用してパラメータ推定を行います。

具体的には次のような2段階モデルを組みます。

  • X = α0 + α1 * Z + α2 * C + ε1
  • Y = β0 + β1 * X + β2 * C + ε2
    • X:治療、Y:結果、C: (複数の)測定済み交絡要因、Z:操作変数
    • αi, βi:係数
    • ε1, ε2:誤差(2変量正規分布を仮定することが多い)


プログラムの実行例はこのようになります。

library(sem)

data(Klein)
Klein$P.lag <- c(NA, Klein$P[-22])
Klein$X.lag <- c(NA, Klein$X[-22])
summary(tsls(C ~ P + P.lag + I(Wp + Wg), instruments=~1 + G + T + Wg + I(Year - 1931) + K.lag + P.lag + X.lag, data=Klein))

> 2SLS Estimates
>
>Model Formula: C ~ P + P.lag + I(Wp + Wg)
>
>Instruments: ~1 + G + T + Wg + I(Year - 1931) + K.lag + P.lag + X.lag
>
>Residuals:
>   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
> -1.890  -0.616  -0.246   0.000   0.885   2.000 
>
>            Estimate Std. Error t value  Pr(>|t|)
>(Intercept)  16.5548    1.46798 11.2772 2.587e-09
>P             0.0173    0.13120  0.1319 8.966e-01
>P.lag         0.2162    0.11922  1.8137 8.741e-02
>I(Wp + Wg)    0.8102    0.04474 18.1107 1.505e-12
>
>Residual standard error: 1.1357 on 17 degrees of freedom

tsls関数の「instruments=」の部分で操作変数のモデルを指定します。


この関数は結果変数が連続変数になっています。

2値結果変数の場合は直接tsls関数を直接利用できませんが、論文ベースでこのような提案がされています。

http://aje.oxfordjournals.org/content/169/3/273.abstract


ちなみに、未調整の線形モデルではこのような結果となります。

lm(C ~ P + P.lag + I(Wp + Wg), data = Klein)

>Coefficients:
>(Intercept)            P        P.lag   I(Wp + Wg)  
>   16.23660      0.19293      0.08988      0.79622  

調整することでPの効果が小さくなり、P.lagの効果が大きくなっているようです。