Web における A/B テストでは、ランダム化単位と分析単位が異なるということがよくあります。 例えば、A か B かのランダム割付がユーザーごとに行われるのに対して、評価指標として分析したいのはページビューごとのクリック率だったりします。 この場合、ランダム化単位はユーザー、分析単位はページビューということになります。
一般に、Web の A/B テストはサンプルサイズが非常に大きいため、中心極限定理により平均値は正規分布に従うと仮定できるので、Z検定がよく使われます。 しかし、ランダム化単位と分析単位が異なると、この Z検定にまずいことが起きます。 具体的には、一人のユーザーが複数のページビューを発生させることができ、それぞれのユーザーは異なるクリック率を持つため、Z検定が仮定する独立同分布 (i.i.d.) に違反してしまいます。
データがこのような性質を持つと、A と B で評価指標に差が無い場合でも、Z検定の p値は小さめに出やすくなります。 （シミュレーションのコードはこの記事の最後にあります）
したがって、A と B に差が無い場合に、誤って有意差があると判定するリスクが高くなってしまいます。
デルタ法を使った統計的仮説検定
この問題に対して、Deng et al. (2018) では、デルタ法を使った解決策を提示しています。 相関のあるデータにおいて、Z検定の p値が小さくなりやすいのは、分散（ひいては標準誤差）が過小推定されることが原因です。 そこで、デルタ法を使って、データ内相関を考慮した形で分散推定式を導出します。 Z検定の分散推定式をこれに置き換えることで、分散の過小推定を防ぐことができます。
今回、この検定手法を簡単に実行できるパッケージ deltatest を CRAN に登録しました。 インストールするには次を実行します。
install.packages("deltatest")
利用するには、ユーザー（ランダム化単位）ごとに評価指標の分子と分母を集計したデータを準備します。
# A tibble: 200 × 4
user group click pageview
<int> <chr> <int> <int>
1 1 B 7 8
2 2 A 4 8
3 3 B 5 11
4 4 A 3 9
5 5 A 3 12
6 6 B 3 9
7 7 A 6 14
8 8 A 3 8
9 9 A 3 9
10 10 A 5 11
# ℹ 190 more rows
このデータに対して、デルタ法を使った Z検定を実行するには、次のようにします。
library(deltatest) deltatest(data, click / pageview, by = group)
Two Sample Z-test Using the Delta Method
data: click/pageview by group
Z = 1.0601, p-value = 0.2891
alternative hypothesis: true difference in means between control and treatment is not equal to 0
95 percent confidence interval:
-0.02525010 0.08474725
sample estimates:
mean in control mean in treatment difference
0.47876448 0.50851305 0.02974857
この手法を使うと、A と B に差が無い場合、検定の p値は一様分布します。
したがって、A と B に差が無い場合に、有意差を誤って検出するリスク（タイプ I エラー率）を正しい値（有意水準）に抑えることができます。
詳しい使い方については次をご参照ください。
参考
- deltatest: Statistical Hypothesis Testing Using the Delta Method
- Deng, A., Knoblich, U., & Lu, J. (2018). Applying the Delta Method in Metric Analytics: A Practical Guide with Novel Ideas. Proceedings of the 24th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining. doi:10.1145/3219819.3219919
- 確率変数の比の分布における平均と分散をデルタ法で求める - 人間だったら考えて
コード
generate_data <- function(n_user) { data <- data.frame() for (i in 1:n_user) { group <- sample(c("A", "B"), size = 1) N <- rpois(1, 9) + 1 prob <- rnorm(1, mean = 0.5, sd = 0.1) if (prob < 0 || prob > 1) next click <- rbinom(N, size = 1, prob) d <- data.frame(user = i, group, click) data <- rbind(data, d) } data } ztest <- function(data) { data_A <- subset(data, group == "A") data_B <- subset(data, group == "B") mean_A <- mean(data_A$click) mean_B <- mean(data_B$click) diff <- mean_B - mean_A se2_A <- var(data_A$click) se2_B <- var(data_B$click) se <- sqrt(se2_A / nrow(data_A) + se2_B / nrow(data_B)) z <- diff / se p <- 2 * pnorm(-abs(z)) list(p.value = p) } library(mirai) daemons(0) daemons(10, seed = 314) m <- mirai_map(1:10, function(x) { p_values <- double(500) for (i in 1:500) { data <- generate_data(n_user = 200) p_value <- ztest(data)$p.value p_values[i] <- p_value } p_values }, generate_data = generate_data, ztest = ztest) p_value <- m[.flat] library(dplyr) df <- data.frame(p_value) |> mutate(range = cut(p_value, breaks = seq(0, 1, by = 0.05))) |> group_by(range) |> summarise(p = factor(ceiling(max(p_value) * 20) / 20), n = n()) |> mutate(prop = n / sum(n)) library(ggplot2) ggplot(df, aes(p, prop)) + geom_col() + geom_hline(yintercept = 0.05, color = "red") + scale_y_continuous(breaks = seq(0, 1, by = 0.01), minor_breaks = NULL) + xlab("p-value") + ylab("proportion") library(mirai) daemons(0) daemons(10, seed = 314) everywhere(library(dplyr)) everywhere(library(deltatest)) m <- mirai_map(1:10, function(x) { p_values <- double(500) for (i in 1:500) { data <- generate_data(n_user = 200) |> group_by(user, group) |> summarise(click = sum(click), pageview = n(), .groups = "drop") p_value <- deltatest(data, click / pageview, by = group, quiet = TRUE)$p.value p_values[i] <- p_value } p_values }, generate_data = generate_data, ztest = ztest) p_value <- m[.flat] library(dplyr) df <- data.frame(p_value) |> mutate(range = cut(p_value, breaks = seq(0, 1, by = 0.05))) |> group_by(range) |> summarise(p = factor(ceiling(max(p_value) * 20) / 20), n = n()) |> mutate(prop = n / sum(n)) library(ggplot2) ggplot(df, aes(p, prop)) + geom_col() + geom_hline(yintercept = 0.05, color = "red") + scale_y_continuous(breaks = seq(0, 1, by = 0.01), minor_breaks = NULL, limits = c(0, 0.1)) + xlab("p-value") + ylab("proportion")
