データ分析の書記録

読んだ本の内容って忘れてしまいますよね。このブログは分析関係の読んだ本、勉強したことの記録です。

調査観察データの統計科学 3章:傾向スコア

はじめに

今までの話
- 調査観察データの統計科学 1~2章

データから因果推論を行う基本戦略3つのうち、 3章は「戦略3:セミパラメトリックな手法を使う」についてです。

3.1 傾向スコア

傾向スコアe

複数の共変量を1次元に集約できれば、マッチングや層別の問題が起こらない(にくい)、という考え方です。

傾向スコアeの推定

傾向スコアeは、割り当てz{0,1}を目的変数、共変量xを説明変数としたロジスティック回帰やプロビット回帰の推定値(z=1の割り当て確率)です。 つまり0~1の値を取ります。

3.2 傾向スコアを用いた具体的な解析方法

共変量を用いたとき同様に、傾向スコアを用いた場合の因果効果の推定も同じです。

マッチング(プロペンシティスコアマッチング)
  • 傾向スコアeが近い人ペアのy1-y0を計算し、平均を取ってE(y_1-y_0)とする
層別解析
  1. 傾向スコアの大小で5つ程度の層に分け
  2. k層ごとに処置群と対照群の平均\bar{y}1\bar{y}0から、効果の推定量\bar{y}1-\bar{y}0を計算し、 層のサンプル数の重みN_k/Nをかける
  3. k層分を合計Σする
線形回帰モデル(共分散分析)
  • 傾向スコアe、割り当て変数zを説明変数とした線形回帰
  • 基本的にはこの回帰モデルが推奨

傾向スコアの使用に際して

前提条件
  • 擬似決定係数や c統計量でモデルによる割付の正判別率が高いかどうかを確かめる
  • 近年の応用例では c統計量が0.8以上であるということが医学系の論文誌でのスタンダードになっているとのこと
傾向スコアの利点
  • 最大の利点:結果変数y1, y0と共変量xとの関係(y~x)を仮定する必要がない➡もし確認したい結果変数が複数あっても、割り当てzは一つであるので、そういう意味でも楽
  • 傾向スコアは1次元に縮約しており、マッチング相手が見つからないという問題にも対処できる
  • 共変量を使用した回帰モデルは統計モデルの選択が重要だがそれも解決する(モデルの誤設定にも強い)
  • 傾向スコア推定のモデル(ロジスティック等)が真のモデルと異なっても、割とロバストな結果が得られるから使いやすい
傾向スコアの問題点
  • 3群以上だったら?・・・結局2群ずつ傾向スコアを推定するしかない
  • 1:1マッチングでは対象者が少ない群にマッチング数を合わせるため多数群データの多くが無駄になる
  • 傾向スコアを用いてもマッチングがある程度恣意的なのは変わらず
  • 傾向スコアを用いた線形回帰:y~ z + eでは、y ~ eの線形関係が前提になるが、傾向スコアeはシグモイド等曲線関係で推定されており0~1の値を取る。そのためy ~ eの関係は曲線関係になり直線関係にはならず無理がある

これらの問題点から、より良い方法が提案されています。

IPW定量…傾向スコアの逆数による重み付け平均推定法


\begin{align}
\hat{E}(y_1) = \sum_{i=1}^{N} \frac{z_iy_i}{e_i} / \sum_{i=1}^{N} \frac{z_i}{e_i} , \hat{E}(y_0) = \sum_{i=1}^{N} \frac{(1-z_i)y_i}{1-e_i} / \sum_{i=1}^{N}\frac{1-z_i}{1-e_i}
\end{align}
hat^がつく\hat{E}()E()の推定量



IPW定量は各々のサンプル i に傾向スコアの逆数で重みをかけて補正したものですね(eはz=1の割り当て確率なのでz=0の時は1-eになっている)。
この方法であれば、平均処置効果ATEだけでなく、各群の周辺期待値\hat{E}(y_1), \hat{E}(y_0)も一緒に求められますね。

3.5 二重にロバストな推定(DR推定)

DR推定量IPW定量のいわば改良版です。 IPW定量の周辺期待値\hat{𝐸}(y_1), \hat{𝐸}(𝑦_0) は実際計算してみればわかりますが、z=1,0各々のサンプル i のe(=それぞれのx)しか使えていません。また、そもそも傾向スコアe:z~xの回帰モデルが正しくない場合(c統計量で確認したとしても)は、推定量が過小/過大評価される可能性があります。
そこでさらにロバストな方法である二重にロバストな推定量(DR推定量)を用いることが提案されています。

y1の二重にロバストな推定量

\begin{align} \hat{E}^{DR}(y_1) =\frac{1}{N}\sum_{i=1}^{N}[ \frac{z_{i}y_{i1}}{e_i} + (1 - \frac{z_{i}}{e_i})g(x_i,\hat{β_1})] \end{align}



この工程はz=1のみのy1~x回帰モデルからy1(z=0,1)の再予測値を出し(g(x_i,\hat{β_1})に該当)、”再予測値に傾向スコアで補正したy1|z=1”と”再予測したy0|z=1(欠測値)”の平均値になります。

同様にy0の二重にロバストな推定量

\begin{align} \hat{E}^{DR}(y0) = \frac{1}{N}\sum_{i=1}^{N}[\frac{(1-z_{i})y_{i0}}{(1-e_i)} + (1 - \frac{(1-z_{i})}{(1-{e_i})})g(x_i,\hat{β_0})] \end{align}

共変量XによるYの再当てはめ値g(x_i,β)を各サンプルの共変量xを用いて推定していることがわかりますね。
またそれにより、E(y)推定のための条件が緩和されていることがポイントです。つまりここでのパラメータ「e:z~xが正しい」or「g(x,β):y1~x, y0~xが正しい」というどちらかが満たされていればDR推定が使用できる、ということになります。

3.6 TET (ATT)の推定

(独立変数を条件づけたときの結果変数の分布の母数推定)

さてさて、実務ですと平均処置効果ATEより実際に知りたいことはz=1のときの効果TET (ATT)であることが多いです。

TETのおさらい

TET = E(y1-y0|z=1) =  E(y_1|z=1) - E(y_0|z=1) *
*:ただしz=1群のy0は観測されない

では傾向スコアで補正したTETはどうなるでしょうか。 \begin{align} \hat{E}(y_1|z=1) &= \bar{y}_1 ※z=1群のy_1標本平均\\ \hat{E}(y_0|z=1) &= \frac{\sum_{i=1}^{N}\frac{(1-z_i)e_i}{1-e_i}y_i}{\sum_{i=1}^{N}\frac{(1-z_i)e_i}{1-e_i}} ※z=1群のy_0期待値(推定)\end{align} TETの推定量は引き算して \begin{align} TET(ATT) = \hat{E}(y_1-y_0|z=1) = \bar{y}_1 - \frac{\sum_{i=1}^{N}\frac{(1-z_i)e_i}{1-e_i}y_i}{\sum_{i=1}^{N}\frac{(1-z_i)e_i}{1-e_i}} \end{align} となります。

※ ちなみにz=0の時の効果TEU(ATU)は、 \begin{align} \hat{E}(y_0|z=0) &= \bar{y}_0 ※z=0群のy_0標本平均\\ \hat{E}(y_1|z=0) &= \frac{\sum_{i=1}^{N}\frac{z_i(1-e_i)}{e_i}y_i}{\sum_{i=1}^{N}\frac{z_i(1-e_i)}{e_i}} ※z=0群のy_1期待値(推定)\end{align} \begin{align} TEU(ATU)= \frac{\sum_{i=1}^{N}\frac{z_i(1-e_i)}{e_i}y_i}{\sum_{i=1}^{N}\frac{z_i(1-e_i)}{e_i}} - \bar{y}_0 \end{align}

二重にロバストな推定量で計算する場合

"強く無視できる割り当て条件”="割り当てZが共変量Xのみに依存し、結果Yには依存しない"から
E(y_0|z=0,x) = E(y_0|z=1,x)が成り立つので、z=0の観測値y0を利用して
\begin{align} \hat{E}^{DR}(y_0|z=1) =\frac{1}{N_1}\sum_{i=1}^{N}[ \frac{(1-z_{i})(e_i)y_{i}}{(1-e_i)} + (1 - \frac{(1-z_{i})}{(1-e_i)})g(x_i,\hat{β_0})] \end{align}

z=0のみの人を対象にしたy0回帰モデルから、全y0の推定値を求めるg(x_i,\hat{β}_0)


E(y_1|z=0,x) = E(y_1|z=1,x)が成り立つので、z=1の観測値y1を利用して
\begin{align} \hat{E}^{DR}(y_1|z=0) =\frac{1}{N_0}\sum_{i=1}^{N}[ \frac{z_{i}(1-e_i)y_{i}}{e_i} + (1 - \frac{z_{i}}{e_i})g(x_i,\hat{β_1})] \end{align}

z=1のみの人を対象にしたy1回帰モデルから、全y1の推定値を求めるg(x_i,\hat{β}_1)


再掲:おさらい
xを条件づけて欠測部分の期待値を傾向スコアで補正推定。

  z=1 z=0
介入を受けた場合:y1 y_1|z=1 y_1|z=0
介入を受けない場合:y0 y_0|z=1 y_0|z=0
共変量X x x

灰文字 は観測できない値

appendix:3群以上の時どうする? 一般化傾向スコアe

2群での効果を求めるために、いろいろと勉強してきましたが、3群以上だった場合にはどう対応すれば良いでしょう。

  • ABC群なら、A vs BC、B vs AC、C vs ABと考えてそれぞれでロジやプロビット回帰を行って傾向スコアを各群で算出する
  • 名義ロジスティック回帰を使う