『効果検証入門』第3章_傾向スコア算出に入れるべき変数／入れるべきでない変数

# 1.傾向スコア算出のパターンを変える
## ①RCTの結果
```
# A tibble: 2 × 4
  treatment conversion_rate spend_mean count
      <dbl>           <dbl>      <dbl> <int>
1         0         0.00573      0.653 21306
2         1         0.0125       1.42  21307
```
ATE=ATT=0.747
## ②交絡変数だけを加える
```
 ## 傾向スコアを利用したマッチング
 m_near <- matchit(formula = treatment ~ recency + history + channel,
                   data = biased_data,
                   method = "nearest",
                   replace = TRUE)
 ## マッチング後のデータを作成
 matched_data <- match.data(m_near)
 ## マッチング後のデータで効果の推定
 PSM_result <- lm(data = matched_data,
                  formula = spend ~ treatment) %>%
   tidy()
 PSM_result
```

```
# A tibble: 2 × 5
  term        estimate std.error statistic  p.value
  <chr>          <dbl>     <dbl>     <dbl>    <dbl>
1 (Intercept)    0.665     0.187      3.55 0.000387
2 treatment      0.863     0.226      3.81 0.000137
```

ATT=0.863

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1736314/dde782cf-f4bb-6d72-877b-eb62eed6aa86.png)


## ③媒介変数を加える
```
#実験:媒介変数を加える
m_near_withM <- matchit(formula = treatment ~ recency + history + channel+visit+conversion,
                  data = biased_data,
                  method = "nearest",
                  replace = TRUE)
## マッチング後のデータを作成
matched_data_withM <- match.data(m_near_withM)
## マッチング後のデータで効果の推定
PSM_result_withM <- lm(data = matched_data_withM,
                 formula = spend ~ treatment) %>%
  tidy()
PSM_result_withM
```
```
# A tibble: 2 × 5
  term        estimate std.error statistic    p.value
  <chr>          <dbl>     <dbl>     <dbl>      <dbl>
1 (Intercept)    0.946     0.194      4.88 0.00000105
2 treatment      0.581     0.234      2.48 0.0131    
```

ATT=0.581
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1736314/9b6209f5-d2e4-84da-10cf-fcf2e01e56a5.png)


## ④交絡変数の1つを落とす
```
#実験:交絡変数historryを落とす
m_near_ommit <- matchit(formula = treatment ~ recency + channel ,
                        data = biased_data,
                        method = "nearest",
                        replace = TRUE)
## マッチング後のデータを作成
matched_data_ommit <- match.data(m_near_ommit)
## マッチング後のデータで効果の推定
PSM_result_ommit <- lm(data = matched_data_ommit,
                       formula = spend ~ treatment) %>%
  tidy()
PSM_result_ommit
```
```
# A tibble: 2 × 5
  term        estimate std.error statistic   p.value
  <chr>          <dbl>     <dbl>     <dbl>     <dbl>
1 (Intercept)    0.599     0.196      3.06 0.00224  
2 treatment      0.929     0.233      3.98 0.0000700
```
ATT=0.929
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1736314/b9c741cd-18c7-b5cf-a568-478b8bf4a40f.png)



## 結果まとめ
| RCT| 全交絡変数 |媒介変数入り|交絡変数欠損|
| ---- | ---- | ---- | ---- |
| 0.747 | 0.863 | 0.581| 0.929|

## 考察
・一番RCTに推定が近いのは全交絡変数を入れたとき（②）
・バイアスのかかり方が重回帰分析で言われるバイアスのかかり方と同じ
・傾向スコアのときも因果ダイアグラム、コントロール変数の選択法を意識した方がよさそう

# 2.重みづけのパターンを変える
## ①RCT
```
# A tibble: 2 × 4
  treatment conversion_rate spend_mean count
      <dbl>           <dbl>      <dbl> <int>
1         0         0.00573      0.653 21306
2         1         0.0125       1.42  21307
```
ATE=0.747

## ②全ての交絡変数
```
## 重みの推定
weighting <- weightit(treatment ~ recency + history + channel,
              data = biased_data,
              method = "ps",
              estimand = "ATE")

## 重み付きデータでの効果の推定
IPW_result <- lm(data = biased_data,
                 formula = spend ~ treatment,
                 weights = weighting$weights) %>%
  tidy()
IPW_result
```
```
# A tibble: 2 × 5
  term        estimate std.error statistic     p.value
  <chr>          <dbl>     <dbl>     <dbl>       <dbl>
1 (Intercept)    0.580     0.116      4.99 0.000000601
2 treatment      0.870     0.165      5.27 0.000000136
```
ATE=0.870

## ③媒介変数を入れる
```]
#実験：媒介変数を入れる
weighting_withM <- weightit(treatment ~ recency + history + channel+visit+conversion,
                      data = biased_data,
                      method = "ps",
                      estimand = "ATE")
IPW_result_withM <- lm(data = biased_data,
                 formula = spend ~ treatment,
                 weights = weighting_withM$weights) %>%
  tidy()
IPW_result_withM
```
```
# A tibble: 2 × 5
  term        estimate std.error statistic  p.value
  <chr>          <dbl>     <dbl>     <dbl>    <dbl>
1 (Intercept)   1.01       0.115     8.79  1.59e-18
2 treatment     0.0846     0.163     0.519 6.04e- 1
```
ATE=0.0846

## ④交絡変数を1つ落とす
```
#実験：交絡変数historyを落とす
weighting_ommit <- weightit(treatment ~ recency +  channel,
                            data = biased_data,
                            method = "ps",
                            estimand = "ATE")
IPW_result_ommit <- lm(data = biased_data,
                       formula = spend ~ treatment,
                       weights = weighting_ommit$weights) %>%
  tidy()
IPW_result_ommit
```
```
# A tibble: 2 × 5
  term        estimate std.error statistic      p.value
  <chr>          <dbl>     <dbl>     <dbl>        <dbl>
1 (Intercept)    0.575     0.117      4.91 0.000000897 
2 treatment      0.889     0.166      5.37 0.0000000815
```

## 結果まとめ
| RCT| 全交絡変数 |媒介変数入り|交絡変数欠損|
| ---- | ---- | ---- | ---- |
| 0.747 | 0.870 | 0.0846| 0.889|

## 考察
・一番RCTに推定が近いのは全交絡変数を入れたとき（②）
・バイアスのかかり方が重回帰分析で言われるバイアスのかかり方と同じ
・傾向スコアのときも因果ダイアグラム、コントロール変数の選択法を意識した方がよさそう