Probability & statistics

1 · Axioms & conditional probability intermediate

Every number your AI system produces is secretly a probability, whether you asked for one or not. The classifier that says “this ticket is billing” is really reporting P(billing | text) = 0.87. The LLM that writes the next word rolled a weighted die. The eval dashboard that reads “91% accuracy” is a coin-flip estimate that would land somewhere else if you re-ran it. If you don’t speak probability fluently, all of this feels like magic or luck — and you’ll ship systems you can’t reason about. This section gives you the grammar.

Here’s the whole idea in one picture: probability is bookkeeping for uncertainty. Imagine every possible thing that could happen laid out as a pie chart that must sum to 1 (100% of the time, something happens). An “event” is a slice of that pie, and its probability is just how big the slice is. That’s all the axioms say, formalized: the whole pie is 1, no slice is negative, and non-overlapping slices add up. Nothing mystical — it’s the arithmetic of shares of a whole.

The move that does the real work — and the one worth slowing down for — is conditioning. When you learn that B happened, you don’t recompute the world from scratch; you throw away every slice where B didn’t happen and rescale what’s left so it sums to 1 again. “Given the first die is a 5” means: ignore the other 30 outcomes, keep only the 6 that start with 5, and re-ask your question inside that smaller world. Conditional probability P(A | B) is nothing more than “how much of the surviving B-world is also A.” That single operation — shrink the space, renormalize — is the engine of Bayes, of every classifier, of the whole lesson.

The classic beginner mistake is confusing P(A | B) with P(B | A). These are wildly different and the gap ruins real systems. P(positive test | disease) is not P(disease | positive test) — the first can be 99% while the second is 17% (you’ll prove exactly that in §2). Flip the bar and you get “prosecutor’s fallacy” territory: confidently wrong. A second trap is assuming independence when it isn’t there — treating two correlated signals as if multiplying their probabilities were valid, which quietly double-counts evidence. Independence is a strong claim (P(A ∩ B) = P(A)P(B) exactly); most real-world events are not independent, so earn that assumption before you use it.

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).

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

A doctor tells you a test that’s 99% accurate came back positive for a rare disease. How scared should you be? Your gut screams “99% — I’m almost certainly sick.” The correct answer is often closer to 1 in 6, and the gap between those two numbers has sent people into needless panic and steered doctors toward needless treatment. Bayes’ theorem is the tool that turns the gut’s wrong answer into the right one, and it’s the exact same machinery your classifier runs millions of times a day.

The intuition before the formula: you always start with a prior — what you believed before seeing evidence (“this disease is rare, so probably not”). Evidence arrives (“the test is positive”). Bayes tells you precisely how much to update that belief — not to throw it away, but to nudge it in proportion to how strong the evidence is. The updated belief is the posterior. Think of it as a courtroom: the prior is the presumption you walk in with, each piece of evidence shifts the odds, and the posterior is the verdict after weighing everything. Crucially, strong evidence against a strong prior can still lose — a great test against a genuinely rare disease still leaves you probably healthy.

Why the answer feels so wrong is worth naming, because it’s the single most important idea in this section: the base rate (how common the thing is) dominates. Picture 10,000 people where only 100 have the disease. The test correctly flags 99 of them — but it also wrongly flags 5% of the 9,900 healthy people, which is ~495 false alarms. So among everyone who tested positive (99 + 495 ≈ 594), only 99 are actually sick: ~17%. The false positives, drawn from the huge healthy majority, drown out the true positives. The test didn’t lie about being accurate; the rarity did the damage.

The mistake to watch for — the base-rate fallacy — is everywhere in AI, not just medicine. A fraud detector, a rare-defect classifier, a “is this prompt an attack” guardrail: any time the positive class is rare, a high-recall model buries you in false alarms, and reporting only “95% accurate” hides it (a model that always says “not fraud” is 99.9% accurate and useless). This is exactly why you report precision alongside recall, and why the same test can be trustworthy for a common condition and near-worthless for a rare one — the very contrast you’ll compute in the lab and the spam-filter check below.

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.

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).

3 · Random variables & distributions advanced

You already own three of these distributions — you just haven’t named them. Every time you check whether one answer is correct, you’re looking at a Bernoulli. Every time you count how many of 500 eval items passed, that’s a Binomial. And the reason you can slap error bars on the resulting accuracy is a Normal. Learn these three shapes and most of applied AI statistics stops being intimidating, because nearly every measurement you make collapses onto one of them.

Start with the mental model. A random variable is just a rule that turns a messy real-world outcome into a number you can do math on — “heads” becomes 1, “wrong answer” becomes 0. A distribution is the fingerprint of that variable: it says which numbers are likely and which are rare. Picture three fingerprints. Bernoulli is a single yes/no lever, weighted p toward yes. Binomial is what you get when you pull that lever n times and count the yeses — a lopsided hill centered near np. Normal is the smooth bell curve that shows up, almost eerily, whenever you average lots of little independent things together (the CLT in §5 is why).

One distinction trips people constantly, so pin it down: discrete vs. continuous. Bernoulli and Binomial are discrete — they land on whole, countable values (0 heads, 1 head, 2 heads…), and each value has an honest probability you can read off a PMF (probability mass function). The Normal is continuous — it can be any real number, so the probability of exactly 0.39894 is zero. That’s not a bug; for continuous variables you ask about ranges, and the PDF (probability density function) gives you height, not probability — you integrate an area under it to get a real probability. Reading a PDF value as if it were a probability (and panicking when it exceeds 1, which densities are allowed to do) is the classic beginner error here.

Why memorize the means and variances in the table below? Because they’re your sanity checks. When your eval harness reports 340 correct out of 500 and you expected p≈0.7, the Binomial mean np = 350 tells you at a glance whether that’s in the right ballpark, and the variance np(1−p) tells you how much wobble to expect between runs. In AI terms: per-token correctness is Bernoulli, the count of correct answers on an n-item eval is Binomial, and — the payoff — your accuracy estimate is approximately Normal, which is the whole reason §7 can hand you a confidence interval instead of a naked point that lies about its own certainty.

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)p^k(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).

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).

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

The temperature slider in every LLM playground is a variance dial in disguise — and once you see that, you’ll never guess at it again. Turn it down and the model gets repetitive and safe; turn it up and it gets creative and occasionally unhinged. People treat it as a mysterious “creativity knob,” but it’s pure probability: you’re literally reshaping the distribution the next token is drawn from. This section connects two ideas — expectation and variance — directly to that knob.

Two summaries describe any distribution. Expectation (the mean) is the center of gravity: if you sampled forever and averaged, where would you land? Variance is the spread: how far do individual draws typically stray from that center? A useful picture is a dartboard. Expectation is where the cloud of darts is centered; variance is how wide the cloud is. A low-variance thrower clusters tightly (predictable, boring); a high-variance thrower sprays (surprising, risky). Neither summary alone tells the story — two systems can share a mean and feel completely different because their variance differs.

Now the connection to decoding. An LLM doesn’t “know” the next word; it produces a vector of raw scores (logits), and softmax squashes them into a probability distribution over the vocabulary. The next token is sampled from that distribution — a weighted die roll, exactly the random-variable draw from §3. Temperature enters as a divisor before the softmax: softmax(logits / T). Divide by a small T (<1) and the gaps between logits blow up, so the top token’s share balloons toward 1 — sharp, greedy, low variance. Divide by a large T (>1) and the gaps shrink, so probability leaks toward the long tail — flat, adventurous, high variance. That’s the whole mechanism, and the lab shows the same three logits reshaped at T=0.5, 1.0, and 2.0.

The mistake to avoid is treating temperature as a quality knob rather than a risk knob. Raising T doesn’t make the model “smarter” or more creative in any deep sense — it just widens the sampling variance, which means more diverse outputs and more chances to draw a low-probability token that’s simply wrong. For a code generator or a JSON-emitting tool call you want low variance (small T, or greedy T→0 which always takes the argmax); for brainstorming you tolerate more. Top-p is a related lever that first discards the improbable tail (keep the smallest set of tokens summing to p, then renormalize) so a high temperature can’t reach truly absurd tokens. All of it is applied conditional probability — you’re choosing how much randomness your product can afford.

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.

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.

5 · Law of large numbers & the CLT expert advanced

“We tested the new model on 20 examples and it scored 85% — ship it!” is how teams ship regressions. Twenty examples is basically a coin you flipped a handful of times: the number you got is real, but it’s luck as much as skill, and the next 20 could read 70%. The two theorems in this section — the law of large numbers and the central limit theorem — are what separate a trustworthy eval number from a hopeful one, and they underwrite every confidence interval in §7.

The law of large numbers (LLN) is the reassuring one, and you already believe it in your bones: the more you measure, the closer your average gets to the truth. Flip a fair coin 10 times and you might see 70% heads; flip it a million times and you’ll see 50.0-something percent. The wobble doesn’t vanish because luck cancels — it shrinks because each new sample is a smaller and smaller fraction of a growing pile. For you this means: a bigger eval set gives a score you can lean on. There’s no shortcut; you buy precision with sample size.

The central limit theorem (CLT) is the almost-magical one, and it’s worth pausing on because it sounds too good to be true. Take any distribution with finite variance — skewed, lumpy, ugly, doesn’t matter — draw a batch and average it, do that many times, and the histogram of those averages is a smooth bell curve. The mess of the original shape washes out. The analogy: individual people’s heights, incomes, reaction times all have odd distributions, but average enough of anything and the averages pile up symmetrically around the true mean. And the spread of those averages is predictable: variance σ²/n, so it tightens as you gather more data — the LLN and CLT holding hands.

Why an AI engineer should care, concretely, and where people go wrong. Your accuracy is an average of 0/1 correctness scores — a mean of Bernoulli draws. The CLT says that average is approximately Normal, which is the only reason you’re allowed to write “85% ± 3%” and mean it. The common mistake is forgetting the “±” exists at all: reporting a bare eval number as if it were exact, then agonizing over a 1-point difference that’s pure sampling noise. A second trap is trusting the CLT with tiny n or heavy tails — “large n” is a real requirement, and a distribution with near-infinite variance (rare giant outliers) can refuse to behave. The simulation below shows both theorems emerging from nothing but coin flips and uniform draws.

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.

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.

6 · Maximum likelihood estimation expert expert

When you train a neural network, what is the loss function actually asking for? The honest answer is: “pick the parameters that make our training data look least surprising.” That one sentence is maximum likelihood estimation, and it’s the hidden objective behind cross-entropy loss, logistic regression, and the trillion-parameter models you use every day. Understand MLE and “minimize the loss” stops being a ritual and becomes a statement about probability.

The intuition first. Suppose you flipped a coin 10 times and saw 6 heads. Someone asks: what’s the coin’s bias p? You could try p=0.1 — but that makes 6 heads wildly improbable, so it’s a poor explanation of what you saw. You could try p=0.9 — also a stretch. The value that makes your actual observations most probable is p=0.6, and that’s the maximum-likelihood estimate. MLE is just running that reasoning backwards: instead of “given p, how likely is the data,” you fix the data you observed and ask “which p explains it best.” It’s the parameter that would have been least astonished to produce your dataset.

Two technical moves make it practical, and both are worth knowing because you’ll see them in every training loop. First, the likelihood of independent observations is a product of their probabilities — and multiplying thousands of numbers below 1 underflows to zero on a real computer. So we take the log, which turns the fragile product into a stable sum (logs are monotonic, so whatever maximizes the sum maximizes the product). Second, optimizers minimize by convention, so we flip the sign and minimize the negative log-likelihood — and for a classifier that quantity is cross-entropy. “Minimize cross-entropy” and “maximize the likelihood of the training data” are the same instruction wearing different clothes (see A12).

The subtlety to respect — and the reason MLE isn’t the end of the story. Maximum likelihood trusts your data completely, which is dangerous when data is scarce: flip a coin 3 times, see 3 heads, and MLE confidently declares p=1.0 (the coin never lands tails), which is almost certainly wrong. That over-confidence on small samples is precisely why regularization, priors, and Bayesian methods exist — they pull the estimate back toward something reasonable when the evidence is thin. So MLE is the workhorse bridge from probability to optimization, but keep in the back of your mind that “most probable given this data” and “most probably true” part ways when the data is small. The lab confirms the clean case: an analytic formula and a brute-force grid search agreeing on the same p.

Maximum likelihood estimation (MLE) picks the parameter that makes the observed data most probable: maximize the likelihood L(θ) = ∏ P(x_i | θ), 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.

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

This is the section that decides whether you ship — and the one most likely to get you fooled. A new prompt lifts conversion from 10% to 13% and everyone wants to celebrate. But you flipped a noisy coin; the question that pays your salary is whether that 3-point bump is a real improvement or the kind of wobble you’d see from shuffling the same users. Hypothesis testing is the disciplined way to answer it, and the p-value is the most misunderstood number in all of statistics — so we’re going to define it carefully.

The logic runs by contradiction, like a courtroom that presumes innocence. You start with the null hypothesis: “there is no real difference; the variant is exactly as good as the control.” Then you ask: if that null were true, how often would random chance alone hand me a gap at least this big? That probability is the p-value. A small p-value (say 0.035) means “a gap this large would be rare if nothing were really different,” so you doubt the null and call the result significant. The analogy: the null is the defendant, your data is the evidence, and the p-value is how easily an innocent defendant could have produced evidence this incriminating. Small p → hard to explain away as chance → convict.

Now the part everyone gets wrong, stated precisely, because it matters. The p-value is not the probability that the null is true, and it is not the probability your result was a fluke. p=0.035 does not mean “3.5% chance there’s no real effect.” It means: “assuming no real effect, there’s a 3.5% chance of data this extreme.” That’s a statement about the data given the hypothesis, never the hypothesis given the data — the same P(A|B)-vs-P(B|A) flip from §1 that Bayes warned you about. Getting this backwards is how “significant” gets oversold. A confidence interval is the healthier companion: instead of a yes/no verdict it gives the plausible range for the true lift, so you learn both “is it real” and “how big could it be.”

Two traps end the section, and both cost real money. First, significance is not effect size: with a huge sample, a 0.1% lift can be “highly significant” and still not worth the engineering to ship it — which is why you report the CI, not just the p-value. Second, peeking: if you watch the dashboard and stop the moment p dips below 0.05, you’ve massively inflated your false-positive rate, because you gave chance many independent shots at crossing the line. Fix your sample size up front (or use a proper sequential test) — the same discipline you’d apply when comparing two model or eval variants (Ch 5). The lab runs the standard two-proportion z-test and the matching interval end to end.

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.

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.

Checkpoint expert

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

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

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]

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

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

# 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