岩波書店の確率と情報の科学シリーズの星野崇宏 著「調査観測データの統計科学」を読んでいき、まとめていってます。
時折、サンプルデータを利用して実装検証していきます。
今回はやっと第２章で欠測と因果効果の推定についての基礎です。
第２章 欠測データと因果推論
欠測の分類
データが欠損する状況は、主に４種に大別できます。
1.記入漏れや未回答
2.閾値を境とする打ち切りや切断
3.経時データやパネルデータの脱落や磨耗
4.調査や観測への不参加、測定不能
これらを不完全データと呼び、不完全データ出ないデータを完全データと呼びます。
欠測のメカニズム
3種類で考えるのが一般的
- 完全にランダムな欠測(missing completely at random)
欠損はモデリングに用いる変数に依存しない - ランダムな欠測(missing at random)
欠測は観測値に依存する - ランダムでない欠測（not missing at random）
欠測は欠測値や観測していない他の変数にも依存する
関心のある観測ベクトルを$y=(y^t_{obs},y^t_{mis})$とし、$y_{obs}$と$y_{mis}$をそれぞれ観測されているデータのベクトル、観測されていないデータのベクトルとする。
また、欠測するかどうかをインディケータ変数$m$とする。この時の、$y$と$m$の同時分布は選択モデルと呼ばれ次のように表します。
P(y, m | \theta, \phi )=p(y|\theta)p(m|y,\phi)\\
欠測しているデータを積分消去すると$\theta$と$\phi$の完全尤度を求めることができる。
P(y_{obs}, m | \theta, \phi )= \int p(y|\theta)p(m|y,\phi) dy_{mis}\\
この完全尤度を最大化する推定量が、最尤推定量です。
完全にランダムな欠測では、
P( m | y, \theta)= p(m|\theta)\\
であり、
P( y_{obs}, m | \theta, \phi) = p(y_{obs}|\theta)p(m|\phi)\\
と言えます。
P( y_{obs}, m | \theta, \phi) = p(y_{obs}|\theta)p(m|\phi)\\
ランダムな欠測では、
P( m | y, \theta)= p(m|y_{obs},\theta)\\
であり、
P( y_{obs}, m | \theta, \phi) = p(y_{obs}|\theta)p(m|y_{obs}, \phi)\\
と言えます。
欠測モデルからみた調査観測データと因果効果の定義
英語の早期学習の効果を考える場合、早期教育の処理群(z=1)と対照群(z=0)はそれぞれ従属変数$y_1$と$y_0$が得られます。
ここで、観測することが不可能なそれぞれの従属変数$y_0$と$y_1$も存在すると考えます。
ルービンの因果モデルでは、この観測できなかった従属変数を仮想的な従属変数と呼びます。
また、独立変数が取り得る値の数と同数だけ存在する「仮想的な従属変数」を潜在的な効果変量と呼びます。
実際は、観測することができなかった結果も考えることがルービンの因果モデルの特徴で有り、この考えを反実仮想モデル/アプローチと呼ばれます。
ルービンの因果効果とは、次のように表せれる２つの潜在的な結果変量の期待値の差です。
E(y_1 - y_0) = E(y_1)-E(y_0)
また、対象者iにおける２つの潜在的結果変数の差$y_{i1}-y_{i0}$を対象者iに対する因果効果と定義されます。
「もし英語教育を受けた場合の中学校の成績$y_1$」と「もし英語教育を受けていなかった場合の中学校での成績$y_0$」についての個人内の差$y_1-y_0$の期待値が因果効果です。
\hat E(y_1 - y_0) = \frac{1}{N} \sum_{i=1}^N (y_{i1}-y_{i0})
しかし、潜在的な結果のうち必ず一方は欠測するため、推定値を算出することになります。
割り当て$z$がランダム出会った場合の因果効果は、次のように表します。
E(y_1)-E(y_0) = E(y_1|z=1)-E(y_0|z=0)
E(y_1|z=1)=E(y|z=1),E(y_0|z=0)=E(y|z=1)
この時、欠測値は無視すれば良いので、観測された各軍の平均値の差を用いることで因果効果をバイアスなく（不偏）推定することができます。
また、処置群での平均介入効果または因果効果(average trearment dffect on the treared:TET)と、対照群での因果効果(average trearment effect on the untreared:TEU)を考えることができる。
TET=E(y_1 - y_0|z=1)\\
TEU=E(y_1 - y_0|z=0)
一般的には一方を正しく推定することでできる可能性が高いです。
また、TETとTEUを用いると因果効果は,処置群と対照群の構成比率$p(z=1),p(z=0)$を用いると次のように表すことができます。
E(y_1-y_0)=TET \times p(z-1)+TEU \times p(z=0)
Rでルービンの因果効果を推定してみる
文章ばっかりなのは嫌なので、気休めに...
６変数の回帰モデルのある変数dの効果をルービンの因果効果を簡単に推定します。
変数dの割り当ては完全にランダムな状態です。
# モデルの設定
regression<-function(x1,x2,x,x4,x5,d){
e <- rnorm(length(x1),0,0.3)
res <- 0.5*x1+0.3*x2+0.6*x3-0.8*x4-0.3*x5-0.5*d
return(res)
}
set.seed(1234)
N=500
Effect<-0
for (i in 1:100){
x1 <- rnorm(N)
x2 <- rnorm(N)
x3 <- rnorm(N)
x4 <- rnorm(N)
x5 <- rnorm(N)
# 0.5の確率で1と0を生成する
d <- sample(x=c(1,0), size=N, replace = T, prob=c(0.5,0.5))
y <- regression(x1,x2,x3,x4,x5,d)
# 効果量の推定
effect <- mean(y[d==1]) - mean(y[d==0])
Effect[i] <- effect
}
print(Effect)
[1] -0.5677094 -0.4360929 -0.3850071 -0.5420914 -0.5684106 -0.4920748 -0.5211063 -0.2749452 -0.4011044 -0.5830069
mean(Effect)
[1] -0.4771549
sd(Effect)
[1] 0.1005282
10回、同様のサンプリングと推定を行いましたが、分散が大きそうです。
サンプル数が小さい影響ですね。
変数dの効果の平均は-0.5と設定していたので、だいたい合っています。
まぁ、当たり前ですね。
共変量調整による因果効果推定のための条件
因果推定を行う際に、潜在的な結果変数と割り当ての両者に影響を与える様々な今日変量の影響を考慮しなければならないです。
社会科学や医学等では無作為割り当てを強制的に行うことはできない、そのため、共変量の影響を除くことが必要になります。
英語教育の場合であれば、「親の教育意欲が高いほど、子供には早期教育を行う」「収入が高い親ほど、学校外での英語教育を行いやすい」といった間接的な関係が生じる可能性があります。
このように、本来独立変数による結果変量への因果効果がないにも関わらず、見せかけ上の相関が生じる疑似相関が現れる可能性があります。
共変量を除くことで、本来の早期の英語教育を行なった場合の中学校での英語テストの成績への本来の効果を求めることができるようになります。
ここで、潜在的は結果変数$(y_1,y_0)$、割り当て$z$、共変量$x$の同時分布を考えます。
p(y_1,y_0,x,z)=p(z|y_1,y_0,x)p(y_1,y_0|x)p(x)
ここから本来求めたい$p(y_1)$を求めるには次のような操作を行います。
p(y_1)= \int p(z|y_1,y_0,x)p(y_1,y_0|x)p(x)dy_0dzdx
この操作を共変量についての周辺化や共変量調整と呼びます。
割り当てはあくまでも共変量にのみ依存し、結果変数には依存しないと仮定します。この仮定は強く無視できる割り当て条件と呼び、次のようにあわらすことができます。
p(z|y_1,y_0,x)=p(z|x)
一方、ランダムな欠測であるとは、次にように表すことができる。
処置群に割り当てられる確率は処置群で観測される$y_1$に、対象群に割り当てられる確率は対象群で観測される$y_0$に、依存することを意味する。
p(z=j|y_1,y_0,x)=p(z=j|y_j,x) (j=1,0)
これらことから、"強く無視できる割り当て"条件が成立すれば、"ランダムな欠測"条件が成立すると言える。
共変量調整による因果効果の推定法
「平均での特立性」をうまく利用する因果効果を推定する方法は主に４種類存在してます。
-
マッチング
処置群と対照群の値が同じになる対象者のペアを作り、ペアごとの差のん期待値から効果を推定する。
欠点は、一般に完璧なペアを作ることができないこと。特に、医療研究のように健常群と疾患群の比較のように一方が少ない場合はペアを作成することができない。
また、共変量が多くなればなるほど、必要なサンプル数が多く必要になる次元問題、、共変量の分布の重なりがない部分についてはマッチングできないサポート問題が存在する。 -
恒常化・限定
共変量の特定の値のみの対象に限定して、条件付きの効果を推定する。
欠点は、条件付き効果のみしか推定できないこと。 -
層別解析
今日変量量の値をいくつかの層に分けて、各層ごとに共変量が等しくなるように比較する。
欠点は、どのように層別を行うかの恣意性が残ること。 -
回帰モデルを用いる方法
各郡ごとに回帰関数$E(y_1|z=1,x)$,$E(y_1|z=0,x)$をデータから推定、その差の標本平均から効果を推定する。
つまり、潜在的な結果変量の共変量への回帰関数を推定し、その回帰関数を共変量ごろに期待値を求めることで、潜在的な結果変数の周辺期待値を算出する方法である。
これを式で表すと、
\begin{align}
E(y_1-y_0) &= E_x(E(y_1-y_0|x))\\
&= \int (E(y_1|x)-E(y_0|x))p(x)dx \\
&= \int (E(y_1|z=1,x)-E(y_0|z=0,x))p(x)dx
\end{align}
となる。「$z＝1$での$y_1$の$x$への回帰関数」$E(y_1|x)$を$g(x|\beta_1)$、「$z＝0$での$y_0$の$x$への回帰関数」$E(y_0|x)$を$g(x|\beta_0)$とし、$\beta_1$と$\beta_0$のぞれぞれを回帰関数のパラメータとすると、因果効果の推定値は次のように表すことができます。
\frac{1}{N} \sum_{i=1}^N (g(x|\beta_1)-g(x|\beta_0))
回帰モデルを用いた因果効果の推定の問題点
(1) 結果変数と共変量のモデリングが必要
結果変量と共変量の関係を正しくモデル化する必要です。もし、正しいモデリングできていない場合、推定値に大きなバイアスが存在する事になります。
(2) 直接因果効果の推定値は得られない
回帰モデルで推定するものは切片や回帰係数であり、因果効果自体は求めることができません。求められたとしても条件付き効果であり、調整が必要です。
これらを回避するために、ノンパラメトリックな手法、セミパラメトリックな手法が存在します。
特に回帰モデルの誤設定を避けるためにカーネル回帰モデルを利用するのが一般的です。
今回はここまでです。