Experiment Design and Execution Questions

import math

def _norm_cdf(x):
    """Standard normal CDF using erf."""
    return 0.5 * (1 + math.erf(x / math.sqrt(2)))

def _norm_ppf(p):
    """Inverse standard normal CDF (probit) - Acklam's approximation."""
    # Coefficients for approximation
    a = [ -3.969683028665376e+01,  2.209460984245205e+02,
          -2.759285104469687e+02,  1.383577518672690e+02,
          -3.066479806614716e+01,  2.506628277459239e+00 ]
    b = [ -5.447609879822406e+01,  1.615858368580409e+02,
          -1.556989798598866e+02,  6.680131188771972e+01,
          -1.328068155288572e+01 ]
    c = [ -7.784894002430293e-03, -3.223964580411365e-01,
          -2.400758277161838e+00, -2.549732539343734e+00,
           4.374664141464968e+00,  2.938163982698783e+00 ]
    d = [ 7.784695709041462e-03,  3.224671290700398e-01,
          2.445134137142996e+00,  3.754408661907416e+00 ]
    # Define break-points
    plow = 0.02425
    phigh = 1 - plow
    if p <= 0 or p >= 1:
        raise ValueError("p must be in (0,1)")
    if p < plow:
        q = math.sqrt(-2*math.log(p))
        return (((((c[0]*q + c[1])*q + c[2])*q + c[3])*q + c[4])*q + c[5]) / \
               ((((d[0]*q + d[1])*q + d[2])*q + d[3])*q + 1)
    if p > phigh:
        q = math.sqrt(-2*math.log(1-p))
        return -(((((c[0]*q + c[1])*q + c[2])*q + c[3])*q + c[4])*q + c[5]) / \
                ((((d[0]*q + d[1])*q + d[2])*q + d[3])*q + 1)
    q = p - 0.5
    r = q*q
    return (((((a[0]*r + a[1])*r + a[2])*r + a[3])*r + a[4])*r + a[5])*q / \
           (((((b[0]*r + b[1])*r + b[2])*r + b[3])*r + b[4])*r + 1)

def compute_power_two_proportions(p0, p1, n_a, n_b, alpha=0.05):
    """
    Estimate power for a two-sided two-proportion z-test using normal approximation.
    Assumptions:
    - Large-sample normal approximation is valid (np and n(1-p) >= ~5).
    - Two-sided test by default; for one-sided use alpha accordingly and adjust formula.
    - No continuity correction applied.
    Limitations:
    - Inaccurate for small samples or proportions near 0/1.
    - Does not pool variance under H0 (uses alternative variances), which is common for power calc.
    Returns:
    - estimated power (float between 0 and 1)
    """
    # basic validation
    for val, name in [(p0,"p0"), (p1,"p1")]:
        if not (0 <= val <= 1):
            raise ValueError(f"{name} must be in [0,1]")
    if n_a <= 0 or n_b <= 0:
        raise ValueError("n_a and n_b must be positive")

    # standard error under alternative (unpooled)
    se = math.sqrt(p0*(1-p0)/n_a + p1*(1-p1)/n_b)
    if se == 0:
        # no variability: power is 1 if proportions differ, else alpha (type I) irrelevant
        return 1.0 if p0 != p1 else 0.0

    delta = (p1 - p0) / se  # noncentrality (standardized effect)
    z_crit = _norm_ppf(1 - alpha/2)

    # Power = P(Z < -z_crit | alt) + P(Z > z_crit | alt),
    # with Z ~ N(delta, 1)
    power = _norm_cdf(-z_crit - delta) + (1 - _norm_cdf(z_crit - delta))
    return power

WITH first_exposure AS (
  -- pick each user's first exposure and its variant
  SELECT exposure_id, user_id, variant, exposure_time
  FROM (
    SELECT *,
           ROW_NUMBER() OVER (PARTITION BY user_id ORDER BY exposure_time, exposure_id) AS rn
    FROM exposures
  ) t
  WHERE rn = 1
),
purch_in_window AS (
  -- join purchases within 7 days (7*24 hours) after first exposure
  SELECT fe.variant,
         fe.user_id,
         p.amount
  FROM first_exposure fe
  LEFT JOIN purchases p
    ON p.user_id = fe.user_id
   AND p.purchase_time >= fe.exposure_time
   AND p.purchase_time < fe.exposure_time + INTERVAL '7' DAY
)
SELECT
  variant,
  COUNT(DISTINCT user_id) AS exposed_users,
  COALESCE(SUM(amount), 0) AS revenue_7d,
  CASE WHEN COUNT(DISTINCT user_id) = 0 THEN 0
       ELSE COALESCE(SUM(amount), 0) / COUNT(DISTINCT user_id)
  END AS revenue_per_exposed_user_7d
FROM purch_in_window
GROUP BY variant
ORDER BY variant;