class: center, middle, inverse, title-slide .title[ # Module 4: Power and Sample Size ] .subtitle[ ## How Big Does Your Experiment Need to Be? ] --- <style type="text/css"> .remark-code, .remark-inline-code { font-size: 80%; } .remark-slide-content { padding: 1em 2em; } .small { font-size: 80%; } .tiny { font-size: 65%; } .highlight-box { background: #fff3e0; border-left: 4px solid #e65100; padding: 0.5em 1em; margin: 0.5em 0; } .blue-box { background: #e3f2fd; border-left: 4px solid #1565c0; padding: 0.5em 1em; margin: 0.5em 0; } </style> # Course Map <table> <tr><th>#</th><th>Module</th><th>Status</th></tr> <tr><td>1</td><td><a href="../module-01/slides.html">The Experimental Ideal</a></td><td>✓ done</td></tr> <tr><td>2</td><td><a href="../module-02/slides.html">SUTVA and When It Breaks</a></td><td>✓ done</td></tr> <tr><td>3</td><td><a href="../module-03/slides.html">Designing Around Interference</a></td><td>✓ done</td></tr> <tr><td><b>4</b></td><td><b>Power and Sample Size</b> <i>(you are here)</i></td><td>current</td></tr> <tr><td>5</td><td><a href="../module-05/slides.html">Analyzing Experiments</a></td><td>✓ done</td></tr> <tr><td>6</td><td><a href="../module-06/slides.html">Multiple Testing & Subgroups</a></td><td>✓ done</td></tr> <tr><td>7</td><td><a href="../module-07/slides.html">External Validity</a></td><td>✓ done</td></tr> <tr><td>8</td><td><a href="../module-08/slides.html">Beyond the A/B Test</a></td><td>✓ done</td></tr> </table> --- # The Problem: An Underpowered Experiment A ride-sharing platform runs an A/B test: does a new splash screen increase rider sign-ups? They test on 200 users per arm for one week. -- Result: p = 0.38. "No significant effect. Kill the feature." -- But the true effect **is** 2 percentage points (a 20% relative lift!). With only 200 per arm, the experiment had **9% power** -- it almost never would have detected the real effect. -- .highlight-box[ **The most common mistake in experimentation is running underpowered tests.** You don't learn "the feature doesn't work" -- you learn nothing. ] --- name: power-definition # What Is Statistical Power? .blue-box[ **Power** = `\(P(\text{reject } H_0 \mid H_1 \text{ is true})\)` -- the probability of detecting a real effect when it exists. ] -- <img src="slides_files/figure-html/power-diagram-1.png" style="display: block; margin: auto;" /> --- name: four-knobs # The Four Knobs Power depends on exactly **four** things: | Knob | Turn it this way | Power goes... | |------|-----------------|---------------| | Effect size `\((\delta)\)` | Larger | Up | | Sample size `\((n)\)` | Larger | Up | | Variance `\((\sigma^2)\)` | Smaller | Up | | Significance level `\((\alpha)\)` | Larger | Up | -- The formula for a two-sample test with equal arms: `$$n_{\text{per arm}} = \frac{(z_{1-\alpha/2} + z_\beta)^2 \cdot 2\sigma^2}{\delta^2}$$` -- .highlight-box[ **But you should rarely use this formula.** Simulation-based power analysis is more flexible, more intuitive, and handles real-world complications (clustering, non-normal outcomes, covariates). ] --- name: sim-based-power # Simulation-Based Power Analysis The idea is dead simple: 1. **Simulate** data under `\(H_1\)` (with a known effect) 2. **Analyze** it the way you would analyze the real experiment 3. **Check** if you reject `\(H_0\)` 4. **Repeat** 1,000+ times. Power = fraction of rejections. -- .small[ ```r sim_power <- function(n_per_arm, effect, baseline = 0.10, nsim = 2000) { # map_lgl: run the block below once per i in 1:nsim (i is just the run # counter) and collect the nsim TRUE/FALSEs into one logical vector rejections <- map_lgl(1:nsim, function(i) { control <- rbinom(n_per_arm, 1, prob = baseline) treated <- rbinom(n_per_arm, 1, prob = baseline + effect) t.test(treated, control)$p.value < 0.05 # TRUE if this run rejects H0 }) mean(rejections) } ``` ] -- ```r sim_power(n_per_arm = 200, effect = 0.02) # too small ``` ``` ## [1] 0.104 ``` ```r sim_power(n_per_arm = 5000, effect = 0.02) # much better ``` ``` ## [1] 0.8875 ``` --- # Power Curves: Seeing the Trade-offs <img src="slides_files/figure-html/power-curve-1.png" style="display: block; margin: auto;" /> -- You need roughly **5,000-10,000 users per arm** to reliably detect a 2 percentage point lift on a 10% baseline conversion rate. --- name: mde-intro # Minimum Detectable Effect (MDE) .blue-box[ **MDE** = the smallest effect your experiment can reliably detect (at 80% power, alpha = 0.05). ] -- The business doesn't ask "what's our power?" They ask: > "We can run this test for 2 weeks and get 10,000 users per arm. **What's the smallest lift we can detect?**" -- <img src="slides_files/figure-html/mde-curve-1.png" style="display: block; margin: auto;" /> --- # Application: Online Ads An ad platform wants to test a new ad format. Baseline click-through rate (CTR) is 3%. They want to detect a 2 percentage point lift (to 5%). -- .pull-left[ ```r # How many impressions per arm? sizes <- c(100, 200, 500, 1000, 2000) ads_power <- tibble( n = sizes, power = map_dbl(sizes, ~ sim_power( .x, effect = 0.02, baseline = 0.03, nsim = 1000 )) ) ads_power |> knitr::kable(digits = 2) ``` | n| power| |----:|-----:| | 100| 0.10| | 200| 0.17| | 500| 0.37| | 1000| 0.62| | 2000| 0.91| ] -- .pull-right[ At CTR = 3%, the variance is very low ( `\(p(1-p) \approx 0.03\)` ), so detecting a 2pp absolute lift is easier than at a 50% baseline. You need about **500-1,000 impressions per arm** -- quite feasible for an ad platform. But what if you want to detect a **0.1pp** lift? That's ~500,000 per arm. ] --- name: clustering-problem # The Clustering Problem What if randomization is at the **city level**, not the individual level? -- ```r sim_power_cluster <- function(n_clusters, n_per_cluster, effect, icc = 0.05, nsim = 1000) { rejections <- map_lgl(1:nsim, function(i) { cluster_id <- rep(1:n_clusters, each = n_per_cluster) treatment <- rep(sample(rep(c(0, 1), each = n_clusters / 2)), each = n_per_cluster) # Cluster random effect cluster_effect <- rep(rnorm(n_clusters, 0, sd = sqrt(icc / (1 - icc))), each = n_per_cluster) y <- rnorm(n_clusters * n_per_cluster, mean = treatment * effect + cluster_effect, sd = 1) # Cluster-robust inference df <- data.frame(y = y, treatment = treatment, cluster = cluster_id) coef(summary(lm(y ~ treatment, data = df)))["treatment", "Pr(>|t|)"] < 0.05 }) mean(rejections) } ``` --- name: clustering-kills-power # Clustering Kills Power <img src="slides_files/figure-html/cluster-comparison-1.png" style="display: block; margin: auto;" /> --- name: design-effect-icc-m4 # The Design Effect and ICC .blue-box[ **ICC** (Intra-cluster Correlation) = fraction of total variance that is *between* clusters. **Design effect** `\(= 1 + (m - 1) \times \text{ICC}\)`, where `\(m\)` = cluster size. ] -- | Cluster size (m)| ICC = 0.01| ICC = 0.05| ICC = 0.10| |----------------:|----------:|----------:|----------:| | 10| 1.1| 1.4| 1.9| | 50| 1.5| 3.5| 5.9| | 100| 2.0| 6.0| 10.9| | 500| 6.0| 26.0| 50.9| -- .highlight-box[ **Design effect** = how many times more samples you need compared to individual-level randomization. With 500 users per cluster and ICC = 0.05, you need **25x** the sample size. Adding more users per cluster barely helps -- you need **more clusters**. ] --- # Application: Ride-Sharing City-Level Experiment A ride-sharing company tests the **zone-notification** feature (Modules 1–2) in 20 cities (10 treatment, 10 control), ~500 drivers per city. The ICC for driver accept rate within a city is 0.08. -- .pull-left[ **Design effect:** `\(1 + (500 - 1) \times 0.08 = 40.9\)` **Effective sample size:** `\(\frac{10{,}000}{40.9} \approx 245\)` per arm Even with 10,000 drivers total, you effectively have 245. ] -- .pull-right[ | Cities per arm| Power| |--------------:|-----:| | 10| 0.80| | 20| 0.86| | 40| 0.92| | 60| 0.96| | 100| 0.98| With only 10 cities per arm, power is low even for a large effect (0.15 SD). You need **40+ cities per arm** to get adequate power. ] --- # The Same Problem at an Online Retailer **The membership program.** - Member spend is heavy-tailed: a few high-frequency "whale" customers dominate the variance, so `\(\sigma^2\)` is large and the required sample for a given MDE is steep. - Winsorizing the spend outcome or including last year's spend as a covariate sharply reduces the effective variance and cuts the required `\(n\)`. **The delivery rollout.** - Randomization is at the metro level; effective `\(n\)` = number of metros, not customers. - Intra-metro ICC on weekly orders is high: common demand shocks within a market produce a large design effect and a tiny effective sample. - A simulation-based power calculation on the metro-week panel is the honest approach; the standard formula with customer-level `\(\sigma\)` drastically understates the required number of metros. --- name: more-vs-bigger-clusters # More Clusters vs. Bigger Clusters <img src="slides_files/figure-html/more-vs-bigger-1.png" style="display: block; margin: auto;" /> --- name: covariates-for-power # Reducing Variance with Pre-Treatment Covariates You can boost power without more users by **reducing residual variance**: `$$n \propto \frac{\sigma^2}{\delta^2} \quad \Rightarrow \quad \text{cut } \sigma^2 \text{ in half} \Leftrightarrow \text{double } n$$` -- .tiny[ ```r sim_power_covariate <- function(n_per_arm, effect, r_squared, nsim = 1000) { rejections <- map_lgl(1:nsim, function(i) { x_pre <- rnorm(2 * n_per_arm) treat <- rep(c(0, 1), each = n_per_arm) y <- sqrt(r_squared) * x_pre + treat * effect + rnorm(2 * n_per_arm, sd = sqrt(1 - r_squared)) coef(summary(lm(y ~ treat + x_pre)))[2, 4] < 0.05 }) mean(rejections) } cat("Without covariate:", sim_power(2000, effect = 0.1, baseline = 0.5, nsim = 1000), "\n") ``` ``` ## Without covariate: 1 ``` ```r cat("With covariate (R^2=0.3):", sim_power_covariate(2000, effect = 0.1, r_squared = 0.3), "\n") ``` ``` ## With covariate (R^2=0.3): 0.964 ``` ] -- .highlight-box[ **CUPED preview:** Using pre-experiment data (last week's clicks, last month's rides) as a covariate can reduce variance by 30-50%, equivalent to doubling your sample size for free. Module 5 covers the mechanics. ] --- name: baseline-panel-power # Using a Baseline Measurement: DiD vs. ANCOVA A ride-sharing platform is about to test a driver-app change. For every driver you also have **last week's earnings** as a baseline, correlated with this week's earnings at `\(\rho\)` (driver skill and location persist week to week). -- Three ways to test for a lift, from the same data: | Estimator | Uses the baseline how | |---|---| | **Post-only** | Ignores it. Compare this week's earnings, treated vs. control. | | **DiD / change score** | Subtracts it. Compare `\(\Delta = Y - X\)`, treated vs. control. | | **ANCOVA** | Regresses on it. `\(Y \sim \text{treat} + X\)`. | -- Relative to post-only variance `\(2\sigma^2/n\)`, the two baseline-adjusted estimators' variance scales by: `$$\text{DiD: } 2(1-\rho) \qquad \text{ANCOVA: } 1-\rho^2$$` -- .highlight-box[ **DiD is not always the free lunch it looks like.** It only beats post-only once `\(\rho > 0.5\)` -- below that, subtracting a weakly correlated baseline just adds independent noise on top. ANCOVA never has this problem: `\(1-\rho^2 \le 2(1-\rho)\)` for every `\(\rho\)`, so it dominates both, always. ] --- name: baseline-panel-sim-setup # Simulating the Three Estimators Same dataset, three tests. Vary `\(\rho\)` (how strongly baseline predicts follow-up) and see who wins. .tiny[ ```r sim_power_panel <- function(n_per_arm, effect, rho, mu = 800, sigma = 150, nsim = 800) { results <- map(1:nsim, function(i) { treat <- rep(c(0, 1), each = n_per_arm) z1 <- rnorm(2 * n_per_arm); z2 <- rnorm(2 * n_per_arm) baseline <- mu + sigma * z1 # last week followup0 <- mu + sigma * (rho * z1 + sqrt(1 - rho^2) * z2) # this week, untreated followup <- followup0 + treat * effect df <- tibble(baseline, followup, treat) c(post = t.test(followup ~ treat, data = df)$p.value < 0.05, did = t.test(I(followup - baseline) ~ treat, data = df)$p.value < 0.05, ancova = coef(summary(lm(followup ~ treat + baseline, data = df)))["treat", 4] < 0.05) }) colMeans(do.call(rbind, results)) } ``` ] -- .tiny[ ```r sim_power_panel(n_per_arm = 200, effect = 30, rho = 0.2) # weak baseline: DiD hurts ``` ``` ## post did ancova ## 0.53875 0.34875 0.55000 ``` ```r sim_power_panel(n_per_arm = 200, effect = 30, rho = 0.8) # strong baseline: DiD and ANCOVA win big ``` ``` ## post did ancova ## 0.51125 0.88875 0.91500 ``` ] --- name: baseline-panel-results # When Does the Baseline Pay Off? <img src="slides_files/figure-html/panel-power-curve-1.png" style="display: block; margin: auto;" /> -- At `\(\rho = 0.7\)` (earnings persist fairly strongly week to week): post-only gets **~50% power**, DiD gets **~73%**, ANCOVA gets **~80%** -- from the exact same 200 drivers per arm. --- name: pcpanel-concept-main # Where Does Rho Actually Come From? The DiD/ANCOVA formulas need a value for `\(\rho\)`. Last section we just picked one. In practice, when you have real historical data, you don't have to. -- A widely-used approach for panel RCTs (the `pcpanel` calculator behind Burlig, Preonas and Woerman 2020): instead of assuming a distribution and a `\(\rho\)`, **resample your actual historical panel**. -- 1. Take real baseline/follow-up pairs from before any experiment ran. 2. For each simulation: draw units **with replacement**, keeping each unit's own pair intact, so whatever correlation, skew, or clustering really exists comes along for free. 3. Randomly label half "placebo-treated," add a candidate effect, run your planned estimator, record reject/no-reject. 4. Power(effect) = fraction of rejections. Sweep effect sizes for the full curve. -- .highlight-box[ The only assumption left is that history repeats: the correlation in your historical data resembles what the experiment will see. No parametric form, no hand-picked `\(\rho\)`. ] --- name: pcpanel-resample-demo # Live: Power From Resampling, Not From Rho A fixed batch of 30 "already observed" cities (what you'd actually have in hand before running the experiment). We never tell the power function their correlation. .tiny[ ```r set.seed(7) n_hist <- 30 quality <- rnorm(n_hist) # unobserved city factor historical <- tibble( baseline = 800 + 150 * quality + rnorm(n_hist, 0, 90), followup0 = 800 + 150 * quality + rnorm(n_hist, 0, 90) # same factor: real, un-assumed rho ) cat("actual historical correlation:", round(cor(historical$baseline, historical$followup0), 2), "\n") ``` ``` ## actual historical correlation: 0.81 ``` ```r power_resample <- function(effect, panel, n_per_arm = 15, nsim = 1000) { rejections <- map_lgl(1:nsim, function(i) { draw <- panel[sample(nrow(panel), 2 * n_per_arm, replace = TRUE), ] treat <- rep(c(0, 1), each = n_per_arm) followup <- draw$followup0 + treat * effect t.test(I(followup - draw$baseline) ~ treat)$p.value < 0.05 }) mean(rejections) } power_resample(effect = 30, panel = historical) ``` ``` ## [1] 0.106 ``` ] -- The correlation never appears in `power_resample()`; it rides along inside `historical` every time a unit is redrawn. Same logic as the DiD formula, no formula required. --- # Summary: The Power Checklist Before running any experiment, answer these: -- 1. **What is the MDE?** What's the smallest effect that matters for the business? -- 2. **What is the baseline rate and variance?** Look at historical data. -- 3. **Is randomization at the individual level or clustered?** If clustered, how many clusters do you have? What is the ICC? -- 4. **Can you reduce variance?** Pre-treatment covariates, stratification, CUPED? -- 5. **How long do you need to run?** More time = more users, but watch out for novelty effects and contamination. -- .highlight-box[ **Rules of thumb:** (1) 80% power is the standard minimum. (2) You almost always need more samples than you think. (3) For clustered designs, the number of *clusters* matters far more than the number of units per cluster. ] --- # Key Takeaways 1. **Power** = `\(P(\text{reject } H_0 \mid H_1 \text{ true})\)`. Low power means you learn nothing from a null result. 2. **Four knobs:** effect size, sample size, variance, significance level. You usually control `\(n\)` and can reduce `\(\sigma^2\)`. 3. **Simulation-based power** is more flexible than formulas. Simulate your analysis pipeline directly. 4. **MDE** is the question the business actually asks: "Given our traffic, what can we detect?" 5. **Clustering** destroys power via the design effect `\(1 + (m-1) \times \text{ICC}\)`. Add more clusters, not more units per cluster. 6. **Pre-treatment covariates** can reduce variance dramatically -- a free lunch for power. (Details in Module 5.) --- # Exercise Preview In the exercise you will: 1. Build a simulation-based power calculator for proportions 2. Generate power curves for different sample sizes and effect sizes 3. Demonstrate how clustering kills power (compare individual vs. cluster randomization) 4. Compute the design effect and effective sample size 5. Show how adding a covariate recovers power lost to noise See `exercise.R` for the starter code. --- name: course-map-end # Course Map <table> <tr><th>#</th><th>Module</th><th>Status</th></tr> <tr><td>1</td><td><a href="../module-01/slides.html">The Experimental Ideal</a></td><td>✓ done</td></tr> <tr><td>2</td><td><a href="../module-02/slides.html">SUTVA and When It Breaks</a></td><td>✓ done</td></tr> <tr><td>3</td><td><a href="../module-03/slides.html">Designing Around Interference</a></td><td>✓ done</td></tr> <tr><td><b>4</b></td><td><b>Power and Sample Size</b> <i>(just finished)</i></td><td>✓ done</td></tr> <tr><td>5</td><td><a href="../module-05/slides.html">Analyzing Experiments</a></td><td>✓ done</td></tr> <tr><td>6</td><td><a href="../module-06/slides.html">Multiple Testing & Subgroups</a></td><td>✓ done</td></tr> <tr><td>7</td><td><a href="../module-07/slides.html">External Validity</a></td><td>✓ done</td></tr> <tr><td>8</td><td><a href="../module-08/slides.html">Beyond the A/B Test</a></td><td>✓ done</td></tr> </table> --- class: center, middle, inverse # Backup Slides --- name: power-formula-derivation # Backup: Power Formula Derivation .tiny[(for two-sample z-test)] Under `\(H_0\)`: `\(\bar{Y}_T - \bar{Y}_C \sim N\left(0, \frac{2\sigma^2}{n}\right)\)` Under `\(H_1\)`: `\(\bar{Y}_T - \bar{Y}_C \sim N\left(\delta, \frac{2\sigma^2}{n}\right)\)` -- Reject `\(H_0\)` when `\(\frac{\bar{Y}_T - \bar{Y}_C}{\sqrt{2\sigma^2/n}} > z_{1-\alpha/2}\)` -- Power = `\(P\left(\frac{\bar{Y}_T - \bar{Y}_C}{\sqrt{2\sigma^2/n}} > z_{1-\alpha/2} \;\middle|\; H_1\right)\)` -- `\(= P\left(Z + \frac{\delta}{\sqrt{2\sigma^2/n}} > z_{1-\alpha/2}\right)\)` where `\(Z \sim N(0,1)\)` `\(= P\left(Z > z_{1-\alpha/2} - \frac{\delta\sqrt{n}}{\sqrt{2}\sigma}\right)\)` -- Setting power `\(= 1 - \beta\)` and solving for `\(n\)`: `$$n = \frac{2\sigma^2 (z_{1-\alpha/2} + z_{1-\beta})^2}{\delta^2}$$` --- name: icc-derivation # Backup: ICC and Design Effect **ICC** = `\(\rho = \frac{\sigma_b^2}{\sigma_b^2 + \sigma_w^2}\)`, where `\(\sigma_b^2\)` = between-cluster variance and `\(\sigma_w^2\)` = within-cluster variance. -- The variance of the cluster mean is: `$$\text{Var}(\bar{Y}_{c}) = \frac{\sigma_w^2}{m} + \sigma_b^2$$` -- The variance of the treatment effect estimator with `\(J\)` clusters per arm of size `\(m\)`: `$$\text{Var}(\hat{\tau}) = \frac{2(\sigma_w^2 + m\sigma_b^2)}{Jm} = \frac{2\sigma^2}{Jm} \cdot [1 + (m-1)\rho]$$` -- The term `\(\text{DEFF} = 1 + (m-1)\rho\)` is the **design effect**. It inflates the variance relative to individual-level randomization. Effective sample size: `\(n_{\text{eff}} = \frac{Jm}{\text{DEFF}}\)` --- name: mde-formula # Backup: MDE Formula From the power formula, solving for `\(\delta\)` instead of `\(n\)`: `$$\text{MDE} = (z_{1-\alpha/2} + z_{1-\beta}) \cdot \sqrt{\frac{2\sigma^2}{n}}$$` -- For a binary outcome with baseline `\(p\)`: `$$\text{MDE} = (z_{0.975} + z_{0.80}) \cdot \sqrt{\frac{2p(1-p)}{n}} \approx 2.8 \cdot \sqrt{\frac{2p(1-p)}{n}}$$` -- For clustered designs, multiply the standard error by `\(\sqrt{\text{DEFF}}\)`: `$$\text{MDE}_{\text{clustered}} = (z_{1-\alpha/2} + z_{1-\beta}) \cdot \sqrt{\frac{2\sigma^2 \cdot [1 + (m-1)\rho]}{Jm}}$$`