Zenn
🎓

A/Bテストにおける効果の確からしさを導出する

2020/11/23に公開

qiitaにも同じ記事があるのですが、試しにzennにも移植してみました。

はじめに

A/Bテストを行う際、施策を打ったことによる効果を推定するだけでなく、推定値がどれだけ信頼出来るものなのかを計算することはとても重要です。
そこで、効果の推定誤差と信頼区間を導出します。

例) ATEの推定値と上下95%信頼区間

問題設定

あるキャンペーンを打つことで、ECサイトの商品購入率がどれだけ変化するかを推定します。
キャンペーンの対象となる介入群をTT, キャンペーン対象外となる統制群をCC, サンプルサイズをそれぞれnn, mm、母平均をpTp_{T}, pCp_{C}とすると、介入群と統制群におけるユーザーの行動XTX_{T}, XCX_{C}はベルヌーイ分布に従います。
ユーザーが商品を購入したなら1, 購入しなかったら0となります。
なお、ユーザーはA/Bを独立して割り付けられ、ユーザーの行動は独立とします。

XT1,XT2,XTnBern(pT)XC1,XC2,XCmBern(pC) X_{T1}, X_{T2}, \ldots X_{Tn} \sim Bern(p_{T}) \\ X_{C1}, X_{C2}, \ldots X_{Cm} \sim Bern(p_{C})

効果の分布の推定

キャンペーンを打つことが購入率に与える施策の効果をATEATEと表記します。
ATEATEの分布、平均、分散を導出します。

中心極限定理

平均、分散パラメータにμ,σ2\mu, \sigma^2を持つ確率分布FFから独立にサンプリングされた標本があるとすると、

X1,X2,,iidF(μ,σ2) X_1, X_2, \ldots, iid \sim F(\mu, \sigma^2)

中心極限定理より、標本平均Xˉ\bar{X}について次の分布収束が成り立ちます。

limnP(n(Xˉμ)σx)=Φ(x) \lim_{n\rightarrow\infty}P \left( \frac{\sqrt n (\bar{X}-\mu)}{\sigma} \leq x \right) = \Phi(x)

つまり、平均と分散パラメータを持つ任意の確率分布から独立に発生した標本の標本平均が従う分布は、漸近的に正規分布に近づいていきます。
ベルヌーイ分布の母平均をppとすれば母分散はp(1p)p(1-p)で表せるので、サンプルサイズが大きい場合、介入群と統制群の標本平均XTˉ,XCˉ\bar{X_T}, \bar{X_C}は以下の正規分布に従います。

XTˉ=1nnXTiXCˉ=1mmXCiXTˉN(pT,1npT(1pT))XCˉN(pC,1mpC(1pC)) \bar{X_T} = \frac{1}{n}\sum^n X_{Ti} \\ \bar{X_C} = \frac{1}{m}\sum^m X_{Ci} \\ \bar{X_T} \sim N(p_T, \frac{1}{n}p_{T}(1-p_{T})) \\ \bar{X_C} \sim N(p_C, \frac{1}{m}p_{C}(1-p_{C})) \\

効果の定義

効果ATEATEは、2つの潜在的な結果変数の期待値の差と定義します。

ATE=E[XTXC]=E[XT]E[XC]=pTpC ATE = E[X_T - X_C] = E[X_T] - E[X_C] = p_T - p_C

「ユーザーの母集団全てが介入群に割り当てられた際の結果」と「ユーザーの母集団全てが統制群に割り当てられた際の結果」の差の平均です。

実際は、ユーザーは介入群か統制群どちらかにしか割り当てられません。
例えば、あるユーザーがキャンペーン対象者になったとき、その後商品を購入したかは観測できます。
しかし、そのユーザーは既にキャンペーン対象に割り当ててしまっているので、そのユーザーがキャンペーン対象にならなかったとき商品を購入したかどうかはわかりません。
なので、ATEATEは直接データから値を計算することはできません。

そこで、通常A/Bテストは割り当てを完全にランダムに行います。
割り当てがランダムに行われていれば、観測された各群の平均値の差を用いることで効果をバイアスなく(不偏)推定できます。
ATEATEの推定量をATE^\hat{ATE}とおきます。

ATE^=1nnXT,i+1mmXC,i=XTˉXCˉE[ATE^]=E[XTˉ]E[XCˉ]=pTpC=ATE \hat{ATE} = \frac{1}{n}\sum^nX_{T,i} + \frac{1}{m}\sum^mX_{C,i} = \bar{X_T} - \bar{X_C} \\ E[\hat{ATE}] = E[\bar{X_T}] - E[\bar{X_C}] = p_T - p_C = ATE \\

効果の推定量が従う確率分布の導出

効果の推定量は介入群の購入率 - 統制群の購入率で計算します。
介入群の購入率統制群の購入率はそれぞれ正規分布に従うので、それらを元に計算した統計量も正規分布に従います。

ATE^=XTˉXCˉN(pTpC,1npT(1pT)+1mpC(1pC)) \hat{ATE} = \bar{X_T} - \bar{X_C} \sim N(p_T - p_C, \frac{1}{n}p_{T}(1-p_{T})+\frac{1}{m}p_{C}(1-p_{C}))

施策を打ったことによる効果の統計量ATE^\hat{ATE}が従う確率分布の母数は

  • 平均pTpCp_T - p_C
  • 分散1npT(1pT)+1mpC(1pC)\frac{1}{n}p_{T}(1-p_{T})+\frac{1}{m}p_{C}(1-p_{C})
  • 標準偏差1npT(1pT)+1mpC(1pC)\sqrt{\frac{1}{n}p_{T}(1-p_{T})+\frac{1}{m}p_{C}(1-p_{C})}

となります。

推測統計における単語の整理

ここまで様々な記号を前振りなく使ってきました。
一度様々な用語を整理します。

統計調査と推測統計

統計調査は、母集団の特性を調べるためにデータを抽出し、そのデータを加工して母集団に関する有益な情報を導いて判断や決定の材料にするために行います。
母集団の全員に調査するのが難しい場合、母集団から標本を抽出して母集団に対する推測を行います。

推測統計では、母集団に確率モデルを想定し、その確率分布に従う確率変数の実現値としてデータを捉えます。
ある確率分布に従う標本X1,X2,,XnX_1, X_2, \ldots, X_nは確率変数で、その実現値x1,x2,,xnx_1, x_2, \ldots, x_nがデータです。nnを標本の大きさ(sample size)と言います。

母数と統計量・推定量

母集団の平均μ\muや分散σ2\sigma^2など、母集団の特性を表すものを母数と言います。
標本X1,X2,,XnX_1, X_2, \ldots, X_nに基づいた関数で母数を含んでいないものを統計量といい、その確率分布を標本分布と言います。
例えば、標本X1X_1や標本平均XTˉ\bar{X_T}などが統計量です。また、統計量を元に計算するATE^\hat{ATE}も統計量となります。
推測統計においては、統計量を計算することで母数を推測します。
標本から得た統計量をもとに、母集団のパラメータを推定するとき、標本の統計量のことを、推定量といいます。

推定量と推定値について

推定量は確率変数です。確定的に何かしらの値をとるわけではありません。
実際にデータを観測して得られた数字をつかって計算したものが、推定値です。

推定量と推定値の違い

効果の確率分布の母数の推定

実際に観測されたデータを用いて母数を推定します。

母平均の推定

母集団の平均μ\muは母数ですが、母集団の全てのデータを元に直接計算するのは容易ではありません。
そこで、n個のデータを観測したのち、標本平均Xˉ=1nnXi\bar{X}=\frac{1}{n}\sum^n X_iを計算します。
標本平均はその平均が母平均と一致します。E[Xˉ]=μE[\bar{X}]=\muです。
なので、標本平均を通じて母平均を推測することが可能になります。

データを(x1,x2,,xn)(x_1, x_2, \ldots, x_n)とすると、推定値はxˉ=nxi\bar{x}=\sum^n x_iとなります。

母平均の不確実性

ただし、統計量(=推定量)は確率変数なので、不確実性がつきまといます。
10個のデータを元に計算した標本平均と10000個のデータを元にした標本平均では、その値に対する不確実性が大きく変わってきます。
そこで、標準誤差を計算して不確実性を明らかにします。

標準偏差と標準誤差の違い

確率変数XXがあるとき、平均μ=E[X]\mu = E[X]、分散σ2=Var(X)=E[(Xμ)2]\sigma^2 = Var(X) = E[(X-\mu)^2]が定義出来ます。
分散のルートを取ったものを標準偏差σ=Var(X)\sigma = \sqrt{Var(X)}といいます。

統計量も確率変数なので、同様に標準偏差を定義出来ます。
統計量はその平均が母平均になることを通じて、母数の値を推定することによく使われます。
なので、統計量の標準偏差、すなわち真の値に対する推定値の不確実性は「標準誤差」と呼ばれます。

例えば、標本平均Xˉ=1ni=1nXi\bar{X}=\frac{1}{n}\sum_{i=1}^n X_iの標準偏差はσn\frac{\sigma}{\sqrt n}となりますが、Xˉ\bar{X}は母平均μ\muを推定するための統計量であるという観点で、その不確実性を表すσn\frac{\sigma}{\sqrt n}は標準誤差です。

標準偏差の推定

推定量の期待値が母数と等しい場合、不偏推定量と言います。
例えば、標本平均は母平均の不偏推定量です。

母分散σ2\sigma^2の不偏推定量V2V^2は以下のようになります。

V2=1n1n(XiXˉ)2E[V2]=σ2 V^2 = \frac{1}{n-1} \sum^n (X_i - \bar{X})^2 \\ E[V^2] = \sigma^2 \\

V2V^2は不偏分散と呼ばれます。
標本平均は個々のデータの合計をnnで割りましたが、不偏分散はn1n-1で割るのがポイントです。
Pandasの関数df.var()df.std()ではデフォルトで不偏分散が計算されるようになっています。

上の定義に従って不偏分散を計算することも出来ますが、ベルヌーイ分布の分散がp(1p)p(1-p)で計算できることを利用して、標本平均Xˉ\bar{X}や標本分散Xˉ(1Xˉ)\bar{X}(1-\bar{X})を使って不偏分散V2V^2を計算することも出来ます。

E[Xˉ(1Xˉ)]=E[XˉXˉ2]=E[Xˉ]E[Xˉ2]=E[Xˉ](E[Xˉ]2+Var(Xˉ))=pp21n2nVar(Xi)=p(1p)1np(1p)=n1np(1p)E[nn1Xˉ(1Xˉ)]=σ2 \begin{align} E[\bar{X}(1-\bar{X})] &= E[\bar{X} - \bar{X}^2] \\ &= E[\bar{X}] - E[\bar{X}^2] \\ &= E[\bar{X}] - (E[\bar{X}]^2 + Var(\bar{X})) \\ &= p - p^2 - \frac{1}{n^2} \sum^n Var(X_i) \\ &= p(1 - p) - \frac{1}{n} p(1 - p) \\ &= \frac{n-1}{n} p (1 - p) \\ E \left[ \frac{n}{n-1}\bar{X}(1-\bar{X}) \right] &= \sigma^2 \\ \end{align} \\

従って、XXの不偏分散はV2=nn1Xˉ(1Xˉ)V^2 = \frac{n}{n-1}\bar{X}(1-\bar{X})となります。

介入群と統制群の不偏分散は以下のようになります。

VT2=nn1XTˉ(1XTˉ)VC2=mm1XCˉ(1XCˉ) V_T^2 = \frac{n}{n-1} \bar{X_T} (1 - \bar{X_T}) \\ V_C^2 = \frac{m}{m-1} \bar{X_C} (1 - \bar{X_C}) \\

XXの母平均μ\muの推定値として標本平均xˉ\bar{x}を計算します。
XXの母分散σ2=p(1p)\sigma^2=p(1-p)の推定値v2v^2は以下のようになります。

xˉ=1nnxiv2=nn1xˉ(1xˉ) \bar{x} = \frac{1}{n} \sum^n x_i \\ v^2 = \frac{n}{n-1} \bar{x} (1 - \bar{x}) \\

標準偏差の推定値はvvです。

標準誤差の推定

標本平均Xˉ\bar{X}の分散(標準誤差の2乗)は以下のようになります。

SE2=V2n SE^2 = \frac{V^2}{n}

介入群と統制群の標本平均の分散は以下のようになります。

SET2=VT2nSEC2=VC2m SE_T^2 = \frac{V_T^2}{n} \\ SE_C^2 = \frac{V_C^2}{m} \\

従って、効果の分散は以下のようになります。

VATE2=SET2+SEC2=VT2n+VC2m V_{ATE}^2 = SE_T^2 + SE_C^2 = \frac{V_T^2}{n} + \frac{V_C^2}{m}

プールした分散

介入群と統制群が従う確率分布の分散がVTC2V_{TC}^2で等しいと仮定すると、ATE^\hat{ATE}の分散を小さく(=推定精度を高める)ことができます。
等分散性が成り立つか予めチェックが必要です。

VT2=nn1XTˉ(1XTˉ)VC2=mm1XCˉ(1XCˉ)VTC2=(n1)VT2+(m1)VC2(n1)+(m1)=(n1)VT2+(m1)VC2n+m2VATE,pool2=VTC2n+VTC2m V_T^2 = \frac{n}{n-1}\bar{X_T}(1-\bar{X_T}) \\ V_C^2 = \frac{m}{m-1}\bar{X_C}(1-\bar{X_C}) \\ V_{TC}^2 = \frac{(n-1)V_T^2 + (m-1)V_C^2}{(n-1)+(m-1)} \\ = \frac{(n-1)V_T^2 + (m-1)V_C^2}{n+m-2} \\ V_{ATE,pool}^2 = \frac{V_{TC}^2}{n} + \frac{V_{TC}^2}{m} \\

効果の確率分布の母数の推定まとめ

今までに出てきた式をまとめます。

介入群と処置群の確率分布について

XT1,XT2,XTnBern(pT)XC1,XC2,XCmBern(pC)XTˉN(pT,1npT(1pT))XCˉN(pC,1mpC(1pC))ATE^N(pTpC,1npT(1pT)+1mpC(1pC)) X_{T1}, X_{T2}, \ldots X_{Tn} \sim Bern(p_{T}) \\ X_{C1}, X_{C2}, \ldots X_{Cm} \sim Bern(p_{C}) \\ \bar{X_T} \sim N(p_T, \frac{1}{n}p_{T}(1-p_{T})) \\ \bar{X_C} \sim N(p_C, \frac{1}{m}p_{C}(1-p_{C})) \\ \hat{ATE} \sim N(p_T - p_C, \frac{1}{n}p_{T}(1-p_{T})+\frac{1}{m}p_{C}(1-p_{C}))

効果の母数の推定量

平均パラメータの推定量

XTˉ=1nnXTiXCˉ=1mnXCiATE^=XTˉXCˉ \bar{X_T} = \frac{1}{n}\sum^n X_{Ti} \\ \bar{X_C} = \frac{1}{m}\sum^n X_{Ci} \\ \hat{ATE} = \bar{X_T} - \bar{X_C}

分散パラメータの推定量

VT2=nn1XTˉ(1XTˉ)VC2=mm1XCˉ(1XCˉ)VATE2=VT2n+VC2mVTC2=(n1)VT2+(m1)VC2n+m2VATE,pool2=VTC2n+VTC2m V_T^2 = \frac{n}{n-1}\bar{X_T}(1-\bar{X_T}) \\ V_C^2 = \frac{m}{m-1}\bar{X_C}(1-\bar{X_C}) \\ V_{ATE}^2 = \frac{V_T^2}{n} + \frac{V_C^2}{m} \\ V_{TC}^2 = \frac{(n-1)V_T^2 + (m-1)V_C^2}{n+m-2} \\ V_{ATE,pool}^2 = \frac{V_{TC}^2}{n} + \frac{V_{TC}^2}{m} \\

効果の母数の推定値について

観測されたデータを用いて、統計量の実現値である推定値を計算します。
ほぼ推定量と一緒なので一部省略します。

ate^=xTˉxCˉvATE2=vT2n+vC2mvATE,pool2=vTC2n+vTC2m \hat{ate} = \bar{x_T} - \bar{x_C} \\ v_{ATE}^2 = \frac{v_T^2}{n} + \frac{v_C^2}{m} \\ v_{ATE,pool}^2 = \frac{v_{TC}^2}{n} + \frac{v_{TC}^2}{m} \\

効果の推定値の信頼区間の計算

最後に、効果の推定値の不確実性を測るため、ATE^\hat{ATE}の信頼区間を計算します。

仮説検定と信頼区間

仮説検定では、帰無仮説と対立仮説という相対する2つの仮説を立てます。
通常なら母平均μ\muはある値μ0\mu_0だ、というのが帰無仮説です。
それに対し、母平均の値はμ0\mu_0ではないのではないか、という仮説を対立仮説と言います。

確率変数XXの実現値xxを観測することで、帰無仮説が正しい場合にその実現値が観測される確率がどれだけ珍しいのかを計算し、計算結果に基づいて仮説を受容したり棄却します。

ある確率変数が正規分布に従うとします。
分散σ02\sigma_0^2は既知とします。

XN(μ,σ02) X \sim N(\mu, \sigma_0^2)

帰無仮説と対立仮説を以下のように設定します。

H0:μ=μ0H1:μμ0 H_0: \mu = \mu_0 \\ H_1: \mu \neq \mu_0 \\

帰無仮説の棄却域、受容域は以下のようになります。

R={xR|xμ0σ0>zα/2}A={xR|xμ0σ0zα/2} R = \left\{ x \in \mathbb{R} \middle| \frac{|x - \mu_0|}{\sigma_0} > z_{\alpha/2}\right\} \\ A = \left\{ x \in \mathbb{R} \middle| \frac{|x - \mu_0|}{\sigma_0} \leq z_{\alpha/2}\right\}

zα/2z_{\alpha/2}とは標準正規分布のz値です。
標準正規分布Φ=N(0,1)\Phi = N(0,1)からzαz_{\alpha}という値以上が取り出される確率がα\alphaとなります。

標準正規分布表の使い方

ここで、検定方式をXXからμ0\mu_0に反転させると

A={μ0R|xμ0σ0zα/2}={μ0R|xμ0σ0zα/2}={μ0R|xσ0zα/2μ0x+σ0zα/2} A = \left\{ \mu_0 \in \mathbb{R} \middle| \frac{|x - \mu_0|}{\sigma_0} \leq z_{\alpha/2}\right\} \\ = \left\{ \mu_0 \in \mathbb{R} \middle| |x - \mu_0| \leq \sigma_0 z_{\alpha/2}\right\} \\ = \left\{ \mu_0 \in \mathbb{R} \middle| x - \sigma_0 z_{\alpha/2} \leq \mu_0 \leq x + \sigma_0 z_{\alpha/2} \right\}

となります。
帰無仮説が正しい場合、標本XXの実現値がxxとして観測される確率が1α1 - \alpha%以内に収まる区間は[xσ0zα/2,x+σ0zα/2][x - \sigma_0 z_{\alpha/2}, x + \sigma_0 z_{\alpha/2}]となります。
これを、信頼係数1α1 - \alphaの信頼区間と言います。

効果の推定量の信頼区間の計算

信頼区間の計算方法がわかったので、効果の信頼区間を実際に計算します。

効果が従う確率分布の母数

施策を打ったことによる効果の推定量ATE^\hat{ATE}が従う確率分布の母数は

μ=ATE=pTpCσ=1npT(1pT)+1mpC(1pC) \mu = ATE = p_T - p_C \\ \sigma = \sqrt{\frac{1}{n}p_{T}(1-p_{T})+\frac{1}{m}p_{C}(1-p_{C})} \\

でした。

データから母数の推定値を計算

平均の推定値を計算します。
標準偏差の真の値はわからないので、代わりに不偏分散のルートで代用します。

μ^=ate^=xTˉxCˉσ0vATE=1n1xTˉ(1xTˉ)+1m1xCˉ(1xCˉ) \hat{\mu} = \hat{ate} = \bar{x_T} - \bar{x_C} \\ \sigma_0 \approx v_{ATE} = \sqrt{\frac{1}{n-1}\bar{x_T}(1-\bar{x_T})+\frac{1}{m-1}\bar{x_C}(1-\bar{x_C})} \\

ATEATEが従う分布が正規分布の場合(※後述)、施策を打つことによる効果ATEATE1α1 - \alpha信頼区間は

[μ^vATEzα/2,μ^+vATEzα/2] [\hat{\mu} - v_{ATE} \, z_{\alpha/2}, \hat{\mu} + v_{ATE} \, z_{\alpha/2}]

となります。

α\alphazα/2z_{\alpha/2}の関係は

α=0.05zα/2=1.96α=0.1zα/2=1.65 \alpha = 0.05 \leftrightarrow z_{\alpha/2} = 1.96 \\ \alpha = 0.1 \leftrightarrow z_{\alpha/2} = 1.65 \\

なので、例えばα=0.05\alpha=0.05の場合の両側95%信頼区間は

[μ^1.96vATE,μ^+1.96vATE] [\hat{\mu} - 1.96v_{ATE}, \hat{\mu} + 1.96v_{ATE}]

です。

「標本平均 ± 2 × 標準誤差」が両側95%信頼区間の近似になります。
実用上はまずこれを使うところから始めるといいと思います。

コード

import numpy as np
import pandas as pd
import plotly.express as px

## 問題設定
t_imp = 10000
c_imp = 50000
t_click = 100
c_click = 350

t_ctr = t_click / t_imp
c_ctr = c_click / c_imp


## ATEの母平均と母分散の推定値を計算
# 平均
ate_mean = t_ctr - c_ctr

# 分散(等分散性を仮定せず)
ate_var = (t_ctr * (1 - t_ctr) / (t_imp - 1)) + (c_ctr * (1 - c_ctr) / (c_imp - 1)) # 等分散性の仮定なし
ate_se = np.sqrt(ate_var)
# 分散(等分散を仮定)
t_var = (t_ctr * (1 - t_ctr)) * (t_imp / (t_imp - 1)) # 介入群の不偏分散
c_var = (c_ctr * (1 - c_ctr)) * (c_imp / (c_imp - 1)) # 統制群の不偏分散
tc_var_pool = ((t_imp - 1) * t_var + (c_imp - 1) * c_var) / (t_imp + c_imp - 2) # プールした分散
ate_var_pool = tc_var_pool * (1 / t_imp + 1 / c_imp) # ateの分散
ate_se_pool = np.sqrt(ate_var_pool)


## 可視化
df_plot = pd.DataFrame({
    "mean": [ate_mean, ate_mean],
    "2se": [2 * ate_se, 2 * ate_se_pool], # 上下95%信頼区間の近似
    "type": ["equal_var=False", "equal_var=True"]
})
px.bar(df_plot, x="type", y="mean", error_y="2se", height=400, width=400)

等分散性を仮定したほうが、推定誤差がわずかに小さいことがわかります。

感想など

色々調べたお陰で推論統計の概念を整理できました。

  1. ある事象の発生過程に確率分布を仮定し、その形状を決める母数を知りたい。
  2. 母数は観測不可なので、標本をもとに(不偏)推定量を導出する。
  3. 推定量は確率変数なので、確からしさを標準偏差や標準誤差として導出する。
  4. 実際に観測されたデータをもとに平均や標準誤差などの推定値を計算する。

用語がこんがらがってややこしい。。

参考文献

(久保川, 2017) 「現代数理統計学の基礎」
(星野, 2009) 「調査観察データの統計科学」
(安井, 2020) 「効果検証入門〜正しい比較のための因果推論/計量経済学の基礎」

(余談) 正規分布とt分布

信頼区間に正規分布を仮定しましたが、分散が未知の場合ATE^\hat{ATE}は正規分布ではなくt分布に従います。
今回の問題設定ではECサイトの商品購入率なので、十分なデータ数が集まり、t分布と正規分布の形状はほぼ一致します。そのため説明を省きました。
t分布を用いてより厳密な検定や信頼区間を推定したい場合は、zα/2z_{\alpha/2}tα/2t_{\alpha/2}に置き換えてください。
実際に計算する際は、t値を手計算するのではなく、データから直接stats.ttest_indなどのパッケージを用いた方が楽です。

Discussion

ログインするとコメントできます