効果検証の２章では回帰モデルを使い、介入変数の説明力から因果関係を解釈します。
概要
回帰モデルをつくり、介入変数の説明力から因果関係の有無を確かめます。
実践
回帰モデルをつくり、介入変数の説明力を確認してみます。
データ項目
RCTを試してみます。男性向けに配信されたメールの効果を検証してみます。
# install package library("tidyverse") # read data email_data <- read_csv("http://www.minethatdata.com/Kevin_Hillstrom_MineThatData_E-MailAnalytics_DataMiningChallenge_2008.03.20.csv")
データ項目を確認してみます。
- recency：最後の購入からの経過月数
- history_segment：昨年の購入額の階層
- history：昨年の購入額
- mens：昨年に男物の商品を購入しているか
- womens：昨年に女物の商品を購入しているか
- zip_code ：地区分類コード
- newbie：過去12ヶ月に新しくユーザになったか
- channel：昨年に置いてどのチャンネルから購入したか
- segment：どのメールが配信されたか
- visit：メールが配信されてから２週間いないにサイトに訪問したか
- conversion：メールが配信されてから２週間いないに購入したか
- spend：購入した際の購入額
head(email_data,5) # A tibble: 5 x 12 recency history_segment history mens womens zip_code newbie channel segment visit conversion spend <dbl> <chr> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <chr> <dbl> <dbl> <dbl> 1 10 2) $100 - $200 142. 1 0 Surburban 0 Phone Womens E-Mail 0 0 0 2 6 3) $200 - $350 329. 1 1 Rural 1 Web No E-Mail 0 0 0 3 7 2) $100 - $200 181. 0 1 Surburban 1 Web Womens E-Mail 0 0 0 4 9 5) $500 - $750 676. 1 0 Rural 1 Web Mens E-Mail 0 0 0 5 2 1) $0 - $100 45.3 1 0 Urban 0 Web Womens E-Mail 0 0 0
RCTでデータを絞り込む。
介入は男性向けメールと女性向けメールがあるので、まず女性向けメールのレコードを取り除きます。
男性向けのメールが配信されたサンプルと配信されなかったサンプルに絞り込みます。(後から色々と検証できるのでログは残しておいた方が良さそうです)
# remove delivery message for women and add columns of delivery message or not male_df <- email_data %>% filter(segment != "Womens E-Mail") %>% mutate(treatment= if_else(segment == "Mens E-Mail", 1, 0))
女性に対して配信されたメールが削除されて、施策の介入有無を意味するカラムが追加されていることを確認できます。
head(male_df, n=5) # A tibble: 5 x 13 recency history_segment history mens womens zip_code newbie channel segment visit conversion spend treatment <dbl> <chr> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <chr> <dbl> <dbl> <dbl> <dbl> 1 6 3) $200 - $350 329. 1 1 Rural 1 Web No E-Mail 0 0 0 0 2 9 5) $500 - $750 676. 1 0 Rural 1 Web Mens E-Mail 0 0 0 1 3 9 5) $500 - $750 675. 1 1 Rural 1 Phone Mens E-Mail 0 0 0 1 4 2 2) $100 - $200 102. 0 1 Urban 0 Web Mens E-Mail 1 0 0 1 5 4 3) $200 - $350 241. 0 1 Rural 1 Multichannel No E-Mail 0 0 0 0
バイアスのあるデータを作る
意図的にバイアスのあるデータを用意します。
set.seed(1) obs_rate_c <- 0.5 obs_rate_t <- 0.5 bias_data <- male_df %>% mutate(obs_rate_c=if_else( (history > 300) | (recency < 6) | (channel == "Multichannel"), obs_rate_c,1), obs_rate_t=if_else( (history > 300) | (recency < 6) | (channel == "Multichannel"), 1, obs_rate_t), random_number=runif(n=NROW(male_df))) %>% filter((treatment==0 & random_number < obs_rate_c)| (treatment==1 & random_number < obs_rate_t))
作られたデータは次のようになっています。
bias_data # A tibble: 31,863 x 16 recency history_segment history mens womens zip_code newbie channel segment visit conversion spend treatment obs_rate_c obs_rate_t <dbl> <chr> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> 1 6 3) $200 - $350 329. 1 1 Rural 1 Web No E-M… 0 0 0 0 0.5 1 2 9 5) $500 - $750 676. 1 0 Rural 1 Web Mens E… 0 0 0 1 0.5 1 3 9 5) $500 - $750 675. 1 1 Rural 1 Phone Mens E… 0 0 0 1 0.5 1 4 2 2) $100 - $200 102. 0 1 Urban 0 Web Mens E… 1 0 0 1 0.5 1 5 4 3) $200 - $350 241. 0 1 Rural 1 Multic… No E-M… 0 0 0 0 0.5 1 6 5 1) $0 - $100 30.0 1 0 Surburb… 0 Phone Mens E… 0 0 0 1 0.5 1 7 5 6) $750 - $1,0… 828. 1 0 Surburb… 1 Multic… Mens E… 0 0 0 1 0.5 1 8 9 1) $0 - $100 30.0 0 1 Surburb… 1 Phone No E-M… 0 0 0 0 1 0.5 9 11 2) $100 - $200 182. 1 0 Surburb… 0 Phone Mens E… 0 0 0 1 1 0.5 10 2 2) $100 - $200 118. 1 0 Surburb… 0 Web Mens E… 1 0 0 1 0.5 1
回帰分析
バイアスのあるデータで回帰分析を行ます。
biased_reg <- lm(data=bias_data, formula=spend ~treatment+history)
分析結果のレポート
説明変数の説明力を確認したいので、「Coefficients」に着目し分析結果の解釈をします。今回知りたいのは介入効果なので「treatment」に着目します。treatmentの推定結果は 0.9026109であり、その検定におけるP値の値も2.25e-07と非常に小さい値です。よってこの介入によって売上平均が0.9ほど増加すると解釈することができます。
summary(biased_reg) Call: lm(formula = spend ~ treatment + history, data = bias_data) Residuals: Min 1Q Median 3Q Max -4.74 -1.46 -1.26 -0.48 497.74 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 0.3241996 0.1444390 2.245 0.02480 * treatment 0.9026109 0.1743057 5.178 2.25e-07 *** history 0.0010927 0.0003366 3.246 0.00117 ** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 15.36 on 31860 degrees of freedom Multiple R-squared: 0.001339, Adjusted R-squared: 0.001276 F-statistic: 21.35 on 2 and 31860 DF, p-value: 5.406e-10
今回のように「Coefficients」にのみ着目したい場合は以下のようにするとわかりやすくなります。
library("broom") biased_reg_coef <- tidy(biased_reg) biased_reg_coef # A tibble: 3 x 5 term estimate std.error statistic p.value <chr> <dbl> <dbl> <dbl> <dbl> 1 (Intercept) 0.324 0.144 2.24 0.0248 2 treatment 0.903 0.174 5.18 0.000000225 3 history 0.00109 0.000337 3.25 0.00117
回帰分析におけるセレクションバイアスと解決策
セレクションバイアスの問題は解決されてるのか？
回帰分析において、セレクションバイアスが起きているか確認してみます。RCTとバイアスのあるデータで介入変数の説明力を確認します。
rct_reg <- lm(data=male_df, formula=spend~treatment) rct_reg_coef <- summary(rct_reg) %>% tidy() nonrct_reg <- lm(data=bias_data, formula=spend~treatment) nonrct_reg_coef <- summary(nonrct_reg) %>% tidy()
RCTの説明力を確認。
rct_reg_coef # A tibble: 2 x 5 term estimate std.error statistic p.value <chr> <dbl> <dbl> <dbl> <dbl> 1 (Intercept) 0.653 0.103 6.36 2.09e-10 2 treatment 0.770 0.145 5.30 1.16e- 7
バイアスのあるデータの説明力を確認。RCTと比較すると介入効果を過大評価しており、セレクションバイアスの問題は解決されていないことが確認できます。
nonrct_reg_coef # A tibble: 2 x 5 term estimate std.error statistic p.value <chr> <dbl> <dbl> <dbl> <dbl> 1 (Intercept) 0.548 0.127 4.32 0.0000156 2 treatment 0.979 0.173 5.67 0.0000000143
回帰分析におけるバイアスを取り除く方法
セレクションバイアスは人為的に引き起こされます。今回のセレクションバイアスは次の３つによって引き起こされているため、これらを共変量として回帰モデルに追加します。
- recency
- channel
- history
バイアスのあるデータで重回帰分析を行います。
nonrct_mreg <- lm(data=bias_data, formula=spend ~ treatment+recency+channel+history) nonrct_mreg_coef <- tidy(nonrct_mreg)
重回帰分析の結果を確認すると介入(treatment )の説明力は0.847とRCTデータの結果に近づいたことが確認できます。
nonrct_mreg_coef # A tibble: 6 x 5 term estimate std.error statistic p.value <chr> <dbl> <dbl> <dbl> <dbl> 1 (Intercept) 0.502 0.379 1.32 0.185 2 treatment 0.847 0.178 4.74 0.00000211 3 recency -0.0403 0.0259 -1.55 0.121 4 channelPhone -0.00178 0.304 -0.00585 0.995 5 channelWeb 0.226 0.303 0.745 0.456 6 history 0.00103 0.000375 2.74 0.00608
追加した共変量のように本来必要だがモデルから抜け落ちている変数を脱落変数と呼びます。このように脱落変数が抜け落ちたバイアスを「脱落変数バイアス(OVB)」と呼びます。
これは直感的には脱落変数が目的変数Yに対して与える影響が、脱落変数と介入変数との相関を通してZの効果として現れています。実際に回帰モデルで脱落変数バイアスのあるモデルとないモデルで推定し、OVBの値を求めると、OVBの差がtreatment の説明力の差になっていることを確認できます。この脱落変数のように目的変数と介入変数両方に相関のある変数のことを交路因子と呼びます。
CIA (Conditional Independence Assumption)
効果検証のための回帰分析における共変量の選択は理想的にはモデルに含まれていない変数によるOVBが全て0になるような状態を目指す。これはモデルに含めた共変量で条件付けしたときに、介入変数が目的変数と独立している状態になり、CIAと呼ばれ(conditional independence assumption) RCTのように介入Zがランダムに割り振られた状態を作り出します。
CIAを満たす上での問題点
- バイアスの評価ができない問題
OVBはモデル間でのバイアスの変化を示す指標。残りのバイアスを評価できるものではない。なので分析者はセレクションバイアスがどのような理由で発生しどのような変数を加えることでコントロール可能であり、どのような変数を組み合わせることでCIAが満たされるのか仮定する必要があります。
- 必要な共変量がデータにはないと言う問題
つまり完全に影響を取り除くことは難しい場合がある。
詳細は記載されてませんが、対処法として次のような手法が考えられます。
- 操作変数法
- 固定効果モデル
- 差分の差分法(４章で説明。必要な共変量がデータとして入手できない場合の対応)
Sensitivity Analysis
手持ちのデータに含まれない変数がセレクションバイアスを起こしているか評価するための手法。これは重要だと分析者が認識している共変量以外の共変量をモデルから抜くことで、効果の推定値が大きく変動しないか確認する分析手法です。
Post treatment bias
介入によって影響を受けた変数によって起きるバイアスをPost treatment biasといいます。このように介入より後のタイミングで値が決まるような変数はモデルから除外する必要があり、セレクションバイアスが減る可能性があるからといって、OVBの値が0でない値を全てモデルに入れて良いわけではありません。
予測精度と推定効果の関係
本書籍ではモデルの予測精度は考慮していません。それは「モデルのデータに対する説明能力や、道のサンプルに対する予測能力を高めることが効果検証において有用である保証にはならないため」です。例えばRCTを行ったデータセットに置いて介入変数のみを入れて回帰した場合、ほとんどの場合、バイアスが非常に小さい値を得ることができます。しかし、モデルの目的変数の予測能力は非常に低くなります。脱落変数に着目した場合、介入変数と目的変数に相関する変数を加えれば、より正確な介入効果を得ることができますが、単にYの予測を改善するような変数を追加しても、その変数が介入効果と相関がない場合には効果の推定は改善しません。このような理由からYに対する予測能力が一定以上なければ効果検証としてのモデルに価値がないとは言えません。
多重共線性
回帰モデルに含まれている変数のうち2つ以上が強い相関をもつ状況のことです。予測問題に置いて問題視される多重共線性ですが、効果検証に置いても標準誤差が変化することで検定の結果が歪んだものになってしまいます。介入変数における多重共線性が確認できる場合対応が必要ですが、それ以外の変数に置いてはもし起きていても介入効果を検証する上では考慮する必要はありません。
所感
モデルの説明変数を活用しすれば全量を変数に入れOVBを可能なかぎり０に近づけることで、施策のユーザ属性を考慮せずにABテストの評価自動化システムをつくり出すことができそうに思う。評価の敷居値はOVBでやれば良さそう。ただし、リークのように入れてはいけない項目があることを認識して設計する必要がありそうです。
