やりたいこと
E[Y^a] = \theta_0 + \theta_1 aというMarginal Structural Modelのパラメータ推定を考える。IPWを重みとしたWeighted Least Squares(WLS)でフィッティングすることでこれらの解は得られる。そのうち、関心のあるAverage Treatment EffectであるE[Y^1] - E[Y^0]= \theta_1の推定量を導出する。
IPWを用いたATEの推定量としての説明には大抵の場合Horvitz-Thompson推定量が記載されている。しかしWhat IfのTechnical Point 12.1でmodified Horvitz-Thompson推定量(Horvitz-Thompson推定量を重みで除した形)が登場し、こちらの方が好まれるとの記載があった。別名Hajek推定量といい、重みをnormalizedすることでより安定するということで、biasedだがconsistentな推定量だそう。そして実は、IPWを重みとしたWLSを行う場合、用いている推定量はHajek推定量らしいということで、導出を試みるに至った。
※OLSの正規方程式を解く方法とやっていることは変わらないので、統計学入門に載っているような内容ではある。基本的なWLSからの解の導出だが、「医学のための因果推論II Rubi因果モデル」のp99に記載のある式の形を直接得たいというモチベーションでやった。
導出
WLSは、wを重みとしてQ=\sum_{i=1}^{n} w_i(y - \theta_0 - \theta_1 a_i)^2の最小化を考える。aは0か1を取る二値の曝露変数である。nはサンプルサイズ。yは従属変数。
まず、\theta_0の解を得る。
\begin{align*}
\frac{\partial Q}{\partial \theta_0} & = -2 \sum_{i=1}^{n} w_i(y_i - \theta_0 - \theta_1 a_i) \equiv 0\\
\\
\Rightarrow & \sum_{i=1}^{n} w_i(y_i - \theta_0 - \theta_1 a_i) = 0\\
\\
\Rightarrow & \theta_0 = \frac{\sum_{i=1}^{n} w_iy_i - (\sum_{i=1}^{n} w_i a_i)\theta_1}{\sum_{i=1}^{n} w_i}\\
\end{align*}
次に、\sum_{i=1}^{n} w_ia_i^2 = \sum_{i=1}^{n} w_ia_iを利用して\theta_1の解を得る。
\begin{align*}
\frac{\partial Q}{\partial \theta_1} & = -2 \sum_{i=1}^{n} w_ia_i(y_i - \theta_0 - \theta_1 a_i) \equiv 0\\
\\
\Rightarrow & \theta_1 = \frac{\sum_{i=1}^{n} w_ia_iy_i - (\sum_{i=1}^{n} w_i a_i)\theta_0}{\sum_{i=1}^{n} w_ia_i}
\end{align*}
\theta_1を\theta_0を代入して式変形すると、
\begin{align*}
\theta_1 & = \frac{\sum_{i=1}^{n} w_i}{\sum_{i=1}^{n} w_i(1-a_i)} \cdot \frac{\sum_{i=1}^{n} w_ia_iy_i - (\sum_{i=1}^{n} w_ia_i)(\sum_{i=1}^{n} w_iy_i)}{\sum_{i=1}^{n} w_ia_i}
\end{align*}
を得る。ここで、
\begin{align*}
\sum_{i=1}^{n} w_i & = \sum_{i=1}^{n} w_ia_i + \sum_{i=1}^{n} w_i(1-a)\\
\sum_{i=1}^{n} w_iy_i & = \sum_{i=1}^{n} w_ia_iy_i + \sum_{i=1}^{n} w_i(1-a)y_i
\end{align*}
を\theta_1に代入して整理すると、
\begin{align*}
\theta_1 & = \frac{ (\sum_{i=1}^{n} w_i(1-a_i)) (\sum_{i=1}^{n} w_ia_iy_i) - (\sum_{i=1}^{n} w_ia_i)(\sum_{i=1}^{n} w_i(1-a_i)y_i) }{ (\sum_{i=1}^{n} w_ia_i) (\sum_{i=1}^{n} w_i(1-a_i)) }\\
\\
& = \frac{\sum_{i=1}^{n} w_ia_iy_i}{\sum_{i=1}^{n} w_ia_i} - \frac{\sum_{i=1}^{n} w_i(1-a)y_i}{\sum_{i=1}^{n} w_i(1-a)}
\end{align*}
と、無事目的の式を得ることができた。めでたし、めでたし。
個人的な補足
分母の重みの期待値が1になるということで(厳密にはbiasedだけど)期待値的にHorvitz-Thompsonと同じになるだろうことはわかるものの(比の期待値は期待値の比にあらず)、どこから出てきたんだ？と思っていたところに、今井先生のこちらのレクチャーノートで、"Weighted least squares gives automatic normalization"と書いてあることを発見してこの記事を書くに至った。
個人的な補足2
ここら辺の推定量の話はサーベイサンプリングに関係する文脈で出てきたものっぽい。そちらの資料はさらに魔境だったのでもっと強くなってきたら再訪するかもしれない。時期は未定。。。
Discussion