AI EngineeringZero to ProductionHome·About·Contact
Appendix · Advanced AI Engineering · Part 11

Probability & statistics

AI systems are probabilistic: a classifier outputs P(class | input), an LLM samples the next token from a distribution, and every eval number is an estimate with error bars. This lesson builds the theory — probability axioms, Bayes, random variables and distributions, expectation/variance, the law of large numbers & CLT, maximum likelihood, and hypothesis testing — and shows why temperature sampling works and when an A/B result is real.

⏱️ ~2.5 hours🎓 Intermediate → Expert🎲 axioms → A/B testsrunnable
SetupThe simulations use only Python’s standard library (random, math) so every result is reproducible — each block seeds the RNG, so the printed numbers below are exactly what you will get. The A/B-test section uses numpy for convenience and is labelled numpy-required; the same math is done with the stdlib elsewhere. No fabricated numbers.

Learning objectives

  • State the probability axioms and compute with conditional probability.
  • Apply Bayes’ theorem and connect it to a probabilistic classifier.
  • Work with Bernoulli, Binomial and Normal random variables.
  • Compute expectation and variance and interpret them.
  • Explain why temperature/top-p sampling works in LLM decoding.
  • Demonstrate the law of large numbers and the CLT by simulation.
  • Derive maximum-likelihood estimation and link it to model training.
  • Run a hypothesis test + confidence interval for an A/B experiment and judge significance.

1 · Axioms & conditional probability intermediate

Probability assigns each event a number in [0, 1]. The axioms (Kolmogorov): P(Ω) = 1 for the whole sample space, P(A) ≥ 0, and for disjoint events P(A ∪ B) = P(A) + P(B). Conditional probability P(A | B) = P(A ∩ B) / P(B) is the probability of A once you know B happened — the update that all inference is built on. Events are independent when P(A ∩ B) = P(A)P(B).

Example code for learning — review, test, and adapt it before running against real or production systems. Commands can create, change, or delete resources. See the Terms & Disclaimer.
Try it — conditional probability by counting
python# stdlib only
# A fair pair of dice; sample space = 36 equally likely outcomes.
outcomes = [(i, j) for i in range(1, 7) for j in range(1, 7)]

def prob(event):                       # event: a predicate on (i, j)
    hits = [o for o in outcomes if event(o)]
    return len(hits) / len(outcomes)

pA = prob(lambda o: o[0] + o[1] == 8)          # sum is 8
pB = prob(lambda o: o[0] == 5)                 # first die is 5
pAB = prob(lambda o: o[0] + o[1] == 8 and o[0] == 5)
print("P(sum=8)      ", round(pA, 4))          # 0.1389
print("P(first=5)    ", round(pB, 4))          # 0.1667
print("P(sum=8|first=5)", round(pAB / pB, 4))  # 0.1667
P(sum=8)       0.1389
P(first=5)     0.1667
P(sum=8|first=5) 0.1667

P(sum=8) = 5/36 ≈ 0.1389. Given the first die is 5, the second must be 3 — exactly one of six equally likely values — so P(sum=8 | first=5) = 1/6 ≈ 0.1667. Conditioning shrank the sample space from 36 outcomes to the 6 where the first die is 5. That shrink-the-space move is the engine of the Bayes update next.

2 · Bayes’ theorem — the classifier’s update advanced

Bayes’ theorem inverts a conditional: P(H | E) = P(E | H)·P(H) / P(E), where P(E) = Σh P(E | h)P(h) (the law of total probability). It turns a prior P(H) plus likelihood P(E | H) into a posterior P(H | E). Every probabilistic classifier — naive Bayes, and conceptually any model that outputs P(class | input) — is doing this update. The classic trap: a very accurate test for a rare condition still yields mostly false positives.

Try it — the base-rate / medical-test problem
python# stdlib only
p_D        = 0.01     # prior: 1% of people have the condition
sens       = 0.99     # P(test + | disease)      (true-positive rate)
false_pos  = 0.05     # P(test + | no disease)   (1 - specificity)

# Law of total probability for the denominator P(test +):
p_pos = sens * p_D + false_pos * (1 - p_D)
# Bayes:
p_D_given_pos = sens * p_D / p_pos
print("P(test +)         ", round(p_pos, 5))          # 0.0594
print("P(disease | +)    ", round(p_D_given_pos, 5))  # 0.16667
P(test +)          0.0594
P(disease | +)     0.16667

Despite a 99%-sensitive test, a positive result means only a 16.7% chance of disease — because the condition is rare, the 5% false positives among the 99% healthy people swamp the true positives. This base-rate reasoning is why a classifier tuned only for high recall on a rare class floods you with false alarms, and why you report precision, not just accuracy (see Ch 5).

✓ Knowledge check

A spam filter flags 95% of spam (P(flag | spam)=0.95) and wrongly flags 2% of ham (P(flag | ham)=0.02). If 20% of mail is spam, what fraction of flagged mail is actually spam?

Show answer
P(flag) = 0.95·0.20 + 0.02·0.80 = 0.190 + 0.016 = 0.206. P(spam | flag) = 0.190 / 0.206 ≈ 0.922. Because spam is common here (20%) and the false-positive rate is low, the posterior is high — contrast with the 1%-prevalence medical test, where the same-quality test gave only 16.7%. The prior/base rate dominates.

3 · Random variables & distributions advanced

A random variable maps outcomes to numbers. Three you must know: Bernoulli(p) — a single 0/1 trial (a coin, a click, a token being ‘correct’); Binomial(n, p) — the count of successes in n independent Bernoulli trials, PMF C(n,k)pk(1−p)n−k; and the Normal(μ, σ²) — the bell curve that the CLT (§5) makes universal. Discrete variables have a PMF (probability per value); continuous ones have a PDF (density).

Try it — Binomial PMF and the Normal PDF
python# stdlib only
import math
from math import comb

def binom_pmf(k, n, p):
    return comb(n, k) * p**k * (1 - p)**(n - k)

n, p = 10, 0.3
print("P(3 successes in 10):", round(binom_pmf(3, n, p), 5))   # 0.26683
print("mean = n*p          :", n * p)                          # 3.0
print("var  = n*p*(1-p)    :", round(n * p * (1 - p), 4))      # 2.1

def normal_pdf(x, mu=0.0, sigma=1.0):
    return math.exp(-0.5 * ((x - mu) / sigma)**2) / (sigma * math.sqrt(2*math.pi))

print("N(0,1) density at 0 :", round(normal_pdf(0.0), 5))      # 0.39894
P(3 successes in 10): 0.26683
mean = n*p          : 3.0
var  = n*p*(1-p)    : 2.1
N(0,1) density at 0 : 0.39894

The most likely single count for Binomial(10, 0.3) is around its mean np = 3, with PMF ≈ 0.267. The standard Normal peaks at 1/√(2π) ≈ 0.399. In AI, per-token correctness is Bernoulli, the number of correct answers on an eval set of n items is Binomial, and by the CLT your accuracy estimate is approximately Normal — which is what lets you put error bars on it (§7).

DistributionModelsMeanVarianceAI use
Bernoulli(p)one 0/1 trialpp(1−p)single correct/incorrect, one click
Binomial(n,p)successes in n trialsnpnp(1−p)#correct on an eval set of n
Normal(μ,σ²)sums/means (via CLT)μσ²metric estimates, error bars

4 · Expectation, variance & sampling — why temperature works advanced

Expectation E[X] = Σ x·P(x) is the long-run average; variance Var(X) = E[(X−μ)²] is the spread. LLM decoding is sampling from a distribution: the model produces logits, softmax(logits / T) turns them into a probability vector, and the next token is drawn from it. Temperature T reshapes that distribution — T<1 sharpens it (greedier, less variance), T>1 flattens it (more random). Top-p truncates to the smallest set of tokens whose probability sums to p, then renormalizes. This is applied conditional probability.

Try it — temperature reshapes the sampling distribution
python# stdlib only
import math, random

def softmax(logits, T=1.0):
    z = [x / T for x in logits]
    m = max(z)                              # stability (see A4)
    e = [math.exp(x - m) for x in z]
    s = sum(e)
    return [x / s for x in e]

logits = [2.0, 1.0, 0.1]                    # raw next-token scores
for T in (0.5, 1.0, 2.0):
    probs = softmax(logits, T)
    print(f"T={T}: {[round(x, 3) for x in probs]}")

# Draw one token from a distribution (inverse-CDF sampling):
def sample(probs, r):
    c = 0.0
    for i, pr in enumerate(probs):
        c += pr
        if r <= c:
            return i
    return len(probs) - 1

random.seed(7)
probs = softmax(logits, 1.0)
draws = [sample(probs, random.random()) for _ in range(10000)]
freq0 = draws.count(0) / len(draws)
print("empirical P(token 0):", round(freq0, 4), "vs", round(probs[0], 4))
T=0.5: [0.864, 0.117, 0.019]
T=1.0: [0.659, 0.242, 0.099]
T=2.0: [0.502, 0.304, 0.194]
empirical P(token 0): 0.6584 vs 0.659

At T=0.5 the top token gets 86% of the mass (near-greedy); at T=2.0 the distribution flattens toward uniform (more ‘creative’/random). The sampling check confirms the theory: drawing 10,000 tokens gives token 0 about 66% of the time, matching its 0.659 probability — the law of large numbers (§5) in action. This is the temperature knob the API exposes; see A4 for the decoding loop.

Temperature is variance controlLowering temperature reduces the variance of the token you draw — more deterministic, more repetitive. Raising it increases variance — more diverse, more error-prone. Greedy decoding is the T→0 limit (always the argmax). Choosing T is choosing how much randomness your product tolerates.

5 · Law of large numbers & the CLT expert advanced

The law of large numbers (LLN): the sample mean converges to the true mean as n grows — why more eval examples give a more trustworthy score. The central limit theorem (CLT): the sample mean of any finite-variance distribution is approximately Normal for large n, with variance σ²/n. Together they are why we can estimate metrics by averaging and attach Normal error bars — the foundation of the confidence intervals in §7.

Try it — LLN convergence and CLT emergence by simulation
python# stdlib only
import random

# --- LLN: sample mean of a fair coin -> 0.5 as n grows ---
random.seed(42)
def coin_mean(N):
    return sum(1 for _ in range(N) if random.random() < 0.5) / N
print("LLN n=100    :", round(coin_mean(100), 3))       # 0.5
print("LLN n=10_000 :", round(coin_mean(10_000), 4))    # 0.4986
print("LLN n=1e6    :", round(coin_mean(1_000_000), 5)) # 0.49982

# --- CLT: mean of 30 Uniform(0,1) draws is ~Normal ---
random.seed(0)
def sample_mean(n=30):
    return sum(random.random() for _ in range(n)) / n
means = [sample_mean() for _ in range(5000)]
mu  = sum(means) / len(means)
var = sum((m - mu)**2 for m in means) / len(means)
print("CLT mean ~0.5:", round(mu, 4))                   # 0.4998
print("CLT var  ~sigma^2/n:", round(var, 6),
      "vs theory", round((1/12) / 30, 6))               # 0.002702 vs 0.002778
LLN n=100    : 0.5
LLN n=10_000 : 0.4986
LLN n=1e6    : 0.49982
CLT mean ~0.5: 0.4998
CLT var  ~sigma^2/n: 0.002702 vs theory 0.002778

As n grows the coin’s sample mean tightens onto 0.5 (LLN). And the mean of 30 uniforms — a flat distribution — produces a bell-shaped spread centered at 0.5 with variance ≈ σ²/n (Uniform variance 1/12, divided by 30 = 0.002778; simulated 0.002702). That the average of a non-Normal variable becomes Normal is the CLT, and it is precisely why your eval accuracy — an average of 0/1 correctness — has a Normal-shaped uncertainty you can bound.

Bernoulli 0/1 per item sum over n eval set ÷ n = mean accuracy ≈ Normal with error bars

6 · Maximum likelihood estimation expert expert

Maximum likelihood estimation (MLE) picks the parameter that makes the observed data most probable: maximize the likelihood L(θ) = ∏ P(xi | θ), or equivalently the log-likelihood (turns the product into a sum, avoids underflow). For a Bernoulli, the MLE of p is just the sample mean — provable by setting the derivative of the log-likelihood to zero. This is what training is: minimizing negative log-likelihood (cross-entropy, see A12) is maximizing the likelihood of the training data under the model.

Try it — MLE for a Bernoulli, analytic vs numeric
python# stdlib only
import math

data = [1, 1, 0, 1, 0, 0, 1, 1, 1, 0]        # 6 ones out of 10

def log_likelihood(p, data):
    if p <= 0 or p >= 1:                      # avoid log(0)
        return float("-inf")
    return sum(math.log(p) if x else math.log(1 - p) for x in data)

# Analytic MLE for Bernoulli is the sample mean:
p_hat = sum(data) / len(data)
print("analytic MLE p:", p_hat)               # 0.6

# Numeric grid search should agree:
best = max((round(0.01*i, 2) for i in range(1, 100)),
           key=lambda p: log_likelihood(p, data))
print("numeric  MLE p:", best)                # 0.6
print("log-likelihood:", round(log_likelihood(p_hat, data), 5))  # -6.73012
analytic MLE p: 0.6
numeric  MLE p: 0.6
log-likelihood: -6.73012

Six ones out of ten → the MLE is p̂ = 0.6, and a numeric search over the log-likelihood lands on the same value — confirming the analytic result. Training a neural network is the same principle at scale: adjust millions of parameters to maximize the log-likelihood of the training tokens, i.e. minimize cross-entropy loss. MLE is the bridge from probability to optimization (A12).

7 · Hypothesis testing & confidence intervals — is the A/B win real? expert expert

You ship a new prompt and conversion rises from 10% to 13%. Real, or noise? A hypothesis test asks: under the null hypothesis (no difference), how surprising is a gap this big? The p-value is that probability; below a threshold (commonly 0.05) you call it significant. A confidence interval gives the plausible range for the true difference. This is exactly how you decide whether an eval or A/B improvement (Ch 5) is a genuine gain or sampling luck.

Try it — two-proportion z-test + 95% CI for an A/B test [numpy-required]
python# numpy-required
import numpy as np
from math import erf, sqrt

def norm_cdf(x):                       # standard-normal CDF via erf
    return 0.5 * (1 + erf(x / sqrt(2)))

# Control: 100/1000 convert.  Variant: 130/1000 convert.
n_c, x_c = 1000, 100
n_v, x_v = 1000, 130
p_c, p_v = x_c / n_c, x_v / n_v

# Pooled two-proportion z-test (H0: p_c == p_v):
p_pool = (x_c + x_v) / (n_c + n_v)
se = sqrt(p_pool * (1 - p_pool) * (1/n_c + 1/n_v))
z = (p_v - p_c) / se
p_value = 2 * (1 - norm_cdf(abs(z)))            # two-sided
print("rates:", p_c, p_v)
print("z =", round(z, 4), " p =", round(p_value, 5))   # z=2.1027 p=0.03549

# 95% CI for the difference (unpooled SE):
se_d = sqrt(p_c*(1-p_c)/n_c + p_v*(1-p_v)/n_v)
lo, hi = (p_v - p_c) - 1.96*se_d, (p_v - p_c) + 1.96*se_d
print("95% CI on lift:", (round(lo, 5), round(hi, 5)))  # (0.00207, 0.05793)
rates: 0.1 0.13
z = 2.1027  p = 0.03549
95% CI on lift: (0.00207, 0.05793)

The z-statistic is 2.10, giving a two-sided p-value of 0.035 — below 0.05, so the 3-point lift is statistically significant: unlikely to be noise. The 95% CI on the lift is roughly (0.2%, 5.8%) — it excludes 0, the same conclusion, and additionally tells you the true gain could be as small as a fraction of a point. Report the CI, not just ‘significant’: a significant but tiny effect may not be worth shipping.

Significance is not size — and beware peekingA large enough sample makes any tiny difference ‘significant’. Always pair the p-value with the effect size / CI. And do not repeatedly test as data trickles in and stop at the first p<0.05 — that inflates false positives. Fix the sample size (or use a sequential test) up front, exactly as you would when comparing two model/eval variants.

Checkpoint expert

✓ Checkpoint — you can move on when you can…

  • State the probability axioms and compute a conditional probability.
  • Apply Bayes’ theorem and explain the base-rate trap for rare classes.
  • Give the mean and variance of Bernoulli, Binomial and Normal variables.
  • Explain how temperature reshapes the LLM sampling distribution.
  • Demonstrate LLN and CLT and state Var(mean) = σ²/n.
  • Derive the Bernoulli MLE and connect NLL minimization to training.
  • Run a two-proportion test + CI and judge whether an A/B lift is real.
✓ Knowledge check

Your eval set has 200 items and your agent scores 75% accuracy. Roughly how wide is the 95% confidence interval, and what does that imply about comparing it to a 78% baseline?

Show answer
For a proportion, SE = √(p(1−p)/n) = √(0.75·0.25/200) ≈ 0.0306, so the 95% CI is 0.75 ± 1.96·0.0306 ≈ (0.69, 0.81) — about ±6 points. A 78% baseline sits well inside that interval, so with only 200 items you cannot claim the 3-point gap is real; you need a larger eval set (LLN) or a paired test to tighten the estimate.

🪜 Practice ladder beginner → industry

Six graded exercises, easy to real-world. Try each before opening its solution.

Exercise 1 · Conditional probability by simulationBeginner

Context: Conditioning — restricting attention to the cases where B happened — is the core move behind every model that outputs P(y | x). Simulating it makes the definition concrete.

Your task: By simulating two dice, estimate P(sum = 8 | first die = 5) and compare to the exact value 1/6.

Requirements:

  • Use random with a fixed seed
  • Keep only trials where the first die is 5, then measure the fraction with sum 8
  • Print the estimate and the exact 1/6
  • No numpy

💡 Hint: Filter to the conditioning event first, then compute the fraction of those with the target event.

Show solution

Restrict to the conditioning event, then measure the target within it.

import random
random.seed(1)
N = 200_000
d = [(random.randint(1,6), random.randint(1,6)) for _ in range(N)]
given = [o for o in d if o[0] == 5]           # first die = 5
hits  = [o for o in given if sum(o) == 8]
print(round(len(hits)/len(given), 4), round(1/6, 4))
0.1689 0.1667

The estimate lands within ~0.001 of 1/6 — and tightens further as N grows (LLN).

Exercise 2 · Bayes for a spam filterIntermediate

Context: Turning likelihoods and a prior into a posterior is the naive-Bayes classifier in one line; the base rate can make a ‘good’ detector untrustworthy.

Your task: Compute P(spam | flagged) given P(flag|spam)=0.95, P(flag|ham)=0.02, and prior P(spam)=0.20.

Requirements:

  • Compute P(flag) via the law of total probability
  • Apply Bayes to get the posterior
  • Print the posterior rounded
  • No numpy

💡 Hint: P(flag) = P(flag|spam)P(spam) + P(flag|ham)P(ham); then Bayes.

Show solution

Total probability for the denominator, then Bayes.

p_spam = 0.20
p_flag_spam, p_flag_ham = 0.95, 0.02
p_flag = p_flag_spam*p_spam + p_flag_ham*(1-p_spam)
post = p_flag_spam*p_spam / p_flag
print(round(p_flag, 4), round(post, 4))       # 0.206 0.9223
0.206 0.9223

92% — high, because spam is common here; the same test at 1% prevalence would be far weaker.

Exercise 3 · Sample from a temperature-scaled distributionAdvanced

Context: LLM decoding samples tokens from softmax(logits/T); showing the empirical frequencies match the probabilities links sampling theory to real generation.

Your task: Given logits [2.0, 1.0, 0.1], build the T=1.0 distribution, draw 20,000 tokens, and confirm the empirical frequency of each token matches its probability.

Requirements:

  • Implement a numerically stable softmax with temperature
  • Sample via inverse-CDF from a uniform draw
  • Seed the RNG for reproducibility
  • Compare empirical frequencies to the probabilities
  • No numpy

💡 Hint: Accumulate probabilities and pick the first index whose cumulative sum exceeds a uniform draw.

Show solution

Softmax gives the distribution; inverse-CDF sampling draws from it; LLN makes frequencies converge.

import math, random
def softmax(z, T=1.0):
    z=[x/T for x in z]; m=max(z); e=[math.exp(x-m) for x in z]; s=sum(e)
    return [x/s for x in e]
def sample(probs, r):
    c=0.0
    for i,pr in enumerate(probs):
        c+=pr
        if r <= c: return i
    return len(probs)-1

random.seed(3)
probs = softmax([2.0,1.0,0.1], 1.0)
draws = [sample(probs, random.random()) for _ in range(20000)]
freq = [round(draws.count(i)/len(draws),3) for i in range(3)]
print('probs', [round(x,3) for x in probs])   # [0.659, 0.242, 0.099]
print('freq ', freq)
probs [0.659, 0.242, 0.099]
freq  [0.66, 0.243, 0.098]
Exercise 4 · See the CLT emerge from a skewed distributionExpert

Context: The CLT is why averaging a non-Normal quantity (like 0/1 correctness) yields a Normal estimate you can bound. Watching it emerge from a skewed source drives the point home.

Your task: Draw sample means of n=40 exponential-like variables (use max of two uniforms, which is skewed), 5000 times, and check the distribution of means is roughly Normal with variance σ²/n.

Requirements:

  • Use a clearly non-Normal base distribution
  • Average n=40 draws, repeated 5000 times
  • Report the empirical mean and variance of the sample means
  • Compare variance to σ²/n
  • No numpy

💡 Hint: Var(mean) = Var(base)/n; estimate Var(base) from many single draws first.

Show solution

Even a skewed base yields Normal-shaped sample means with variance shrinking as σ²/n.

import random
random.seed(5)
def base():                      # skewed: max of two uniforms
    return max(random.random(), random.random())
singles = [base() for _ in range(200000)]
mu_b = sum(singles)/len(singles)
var_b = sum((x-mu_b)**2 for x in singles)/len(singles)

n = 40
means = [sum(base() for _ in range(n))/n for _ in range(5000)]
mu = sum(means)/len(means)
var = sum((m-mu)**2 for m in means)/len(means)
print(round(mu,4), round(var,6), round(var_b/n,6))
0.6669 0.001393 0.001395

The sample-mean variance (0.001393) matches σ²/n (0.001395), and the means cluster in a bell around 2/3 — the CLT, from a skewed source.

Exercise 5 · MLE by maximizing the log-likelihoodProfessional

Context: Training minimizes negative log-likelihood; implementing MLE for a simple model shows what that objective computes and why the log turns a fragile product into a stable sum.

Your task: Given Bernoulli data, compute the analytic MLE (sample mean) and confirm it maximizes the log-likelihood via a grid search.

Requirements:

  • Implement the log-likelihood as a sum of logs
  • Compute the analytic MLE (sample mean)
  • Grid-search p and confirm the maximizer matches
  • Print the maximized log-likelihood
  • No numpy

💡 Hint: log L(p) = (#ones) log p + (#zeros) log(1-p); its maximizer is #ones / n.

Show solution

The Bernoulli MLE is the sample mean; a grid search over the log-likelihood confirms it.

import math
data = [1,0,1,1,0,1,0,1,1,1,0,1]        # 8 ones / 12
def ll(p):
    if p <= 0 or p >= 1: return float('-inf')
    return sum(math.log(p) if x else math.log(1-p) for x in data)
p_hat = sum(data)/len(data)
grid = [round(0.01*i,2) for i in range(1,100)]
best = max(grid, key=ll)
print(round(p_hat,4), best, round(ll(p_hat),5))
0.6667 0.67 -7.63817

Analytic 0.6667 and the grid’s 0.67 agree (grid step is 0.01) — minimizing NLL is maximizing this.

Exercise 6 · Decide an A/B test with a z-test and confidence intervalIndustry scenario

Context: This is the daily decision behind shipping a prompt/model change: is the observed lift real, and how big could it be? A two-proportion test plus CI is the standard tool.

Your task: Given control 100/1000 and variant 130/1000 conversions, run a two-proportion z-test, report the p-value, and give the 95% CI on the lift — then state the ship/hold decision.

Requirements:

  • Compute the pooled-SE z-statistic and two-sided p-value
  • Compute the 95% CI on the difference (unpooled SE)
  • State whether p < 0.05 and whether the CI excludes 0
  • Give a one-line decision that mentions effect size, not just significance
  • numpy allowed and labelled

💡 Hint: Pooled SE for the test; unpooled SE for the CI; z=1.96 for 95%.

Show solution

Pooled-SE z-test for significance, unpooled-SE interval for the effect size.

# numpy-required (numpy optional here; math suffices)
from math import erf, sqrt
def ncdf(x): return 0.5*(1+erf(x/sqrt(2)))
n_c,x_c,n_v,x_v = 1000,100,1000,130
p_c,p_v = x_c/n_c, x_v/n_v
pp = (x_c+x_v)/(n_c+n_v)
se = sqrt(pp*(1-pp)*(1/n_c+1/n_v))
z = (p_v-p_c)/se; pval = 2*(1-ncdf(abs(z)))
se_d = sqrt(p_c*(1-p_c)/n_c + p_v*(1-p_v)/n_v)
lo,hi = (p_v-p_c)-1.96*se_d, (p_v-p_c)+1.96*se_d
print(round(z,4), round(pval,5))            # 2.1027 0.03549
print(round(lo,5), round(hi,5))             # 0.00207 0.05793
print('ship' if pval < 0.05 and lo > 0 else 'hold')
2.1027 0.03549
0.00207 0.05793
ship

Significant (p=0.035) and the CI excludes 0, so ship — but the lift could be as small as ~0.2 points, so weigh it against rollout cost.

© 2026 studybydoing.in · AI Engineering: Zero to Production · All rights reserved. · About · Privacy Policy · Terms · Contact
Educational content, provided as-is and without warranty. Code samples are examples — review, test, and adapt them before using in production. See the Terms of Use & Disclaimer. Use at your own risk.
© studybydoing.in