Bayesian Thinking for Models

Two engineers estimate a new checkout button's conversion rate. Engineer A, who has shipped dozen of buttons, starts convinced the rate is near 5% and keeps a tight belief there. Engineer B has no experience and starts by saying "somewhere between 0% and 20%, honestly." The first ten users arrive and two of them convert (2/10 = 20%). Engineer A concludes the rate is still near 5%; Engineer B, updating a flat prior, lands near 21%. Both did the same arithmetic correctly. This is Bayesian updating, and the disagreement between them — which shrinks as data arrives — is the feature that frequentist estimation does not have.

Bayes' theorem again, now as a procedure rather than an identity. You start with a prior P(H), observe data D, and compute:

posterior P(H | D) = likelihood P(D | H) × prior P(H) / evidence P(D)

The computational heart is the numerator: the posterior is proportional to likelihood times prior. The evidence term P(D) is the same for every hypothesis you compare, so for finding the argmax hypothesis you can skip it — a fact Naive Bayes exploits later in this topic.

A fully worked micro-example, numbers included. Suppose the only two hypotheses are "conversion rate = 0.05" and "conversion rate = 0.20", with prior probabilities 0.5 each. You observe 2 conversions in 10 independent visitors. The binomial likelihood from Key Distributions for ML gives:

The evidence normalizes: 0.0373 + 0.1510 = 0.1883, and each row divides by it. Observing 20% conversion in ten users made the p = 0.20 world roughly 4× more plausible than the p = 0.05 world — but not certain, because ten users is thin evidence and the prior still holds 20% of the mass.

from math import comb

hypotheses = [0.05, 0.20]
prior = {0.05: 0.5, 0.20: 0.5}
k, n = 2, 10   # 2 conversions in 10 visitors

like = {p: comb(n, k) * p**k * (1 - p)**(n - k) for p in hypotheses}
unnorm = {p: prior[p] * like[p] for p in hypotheses}
z = sum(unnorm.values())
posterior = {p: unnorm[p] / z for p in hypotheses}
print(posterior)   # {0.05: 0.198, 0.20: 0.802} — wait: recompute below

Running it: 0.5×0.0746 = 0.0373 and 0.5×0.3020 = 0.1510, total 0.1883, so the posteriors are 0.0373/0.1883 ≈ 0.20 for p = 0.05 and 0.1510/0.1883 ≈ 0.80 for p = 0.20. (The table above rounds the same computation.) The hypothesis "p = 0.20" is now four times more plausible than its competitor, starting from even odds.

The interview question this topic really tests: "where can the prior and the data disagree?" Answer with the shrinkage view — the posterior compromise is weighted by how much information each side carries. A tight prior from 10,000 historical visitors outweighs 10 new ones; a flat prior lets 10 new ones speak almost unopposed. When data is scarce the prior matters most, and when data is overwhelming the prior washes out no matter where it started. That is the whole Bayesian contract, and the next chapter makes it computable in closed form.

The previous chapter's two-hypothesis table was honest but impractical: real conversion rates live on [0, 1], not two points. Conjugacy solves this. A prior is conjugate to a likelihood when updating produces a posterior in the same family — and for Bernoulli/binomial data the conjugate prior is the Beta(α, β) distribution over p ∈ [0, 1].

The Beta distribution is parameterized by two pseudo-counts: α can be read as "imaginary successes" and β as "imaginary failures." Its mean is α/(α+β), and it is flexible enough to be nearly flat or sharply peaked:

The update rule is the reason to care — add, never retrain:

Beta(α, β) + k successes in n trials → Beta(α + k, β + n − k)

Each observation is one more pseudo-count. With prior Beta(2, 8) — mean 0.20, representing two imaginary conversions in ten — and fresh data of 3 successes in 20 trials, the posterior is Beta(5, 25), whose mean is 5/30 ≈ 0.167. Note how the posterior mean is exactly a weighted average: (2 + 3) / (10 + 20) = 0.167, blending prior and data by their pseudo-counts. Increasing the prior's strength (Beta(20, 80) + same data → Beta(23, 105), mean ≈ 0.180) pulls the answer back toward the prior — prior strength is the regularization knob on a probability estimate. The analogy to ridge regression's L2 penalty is deliberate: both shrink estimates toward a prior when data is thin.

import numpy as np
from scipy.stats import beta

def update(alpha, beta_, k, n):
    """Beta-Binomial conjugate update; returns posterior parameters."""
    return alpha + k, beta_ + n - k

a, b = update(2, 8, 3, 20)                 # (5, 25)
print(a / (a + b))                          # 0.1667 — posterior mean
lo, hi = beta.interval(0.95, a, b)          # 95% credible interval
print(np.round([lo, hi], 3))                # [0.058 0.317]

# Day-by-day updating: the posterior is reusable as tomorrow's prior
alpha, beta_ = 2, 8
for day_k, day_n in [(1, 10), (2, 15), (0, 12), (4, 20)]:
    alpha, beta_ = update(alpha, beta_, day_k, day_n)
    print(f"posterior mean {alpha / (alpha + beta_):.3f}")

Three consequences worth saying in an interview. First, Laplace's rule of succession — adding one pseudo-success and one pseudo-failure (the Beta(1,1) prior) — is exactly the add-one smoothing you use in Naive Bayes and language models; smoothing is Bayesian updating wearing a product costume. Second, the interval from beta.interval is a credible interval: "95% of the posterior mass lies here," a direct statement about the parameter. Compare that to the confidence intervals you will build in Hypothesis Tests and Confidence Intervals, which mean something subtly different and weaker. Third, sequential updating is free — a production conversion tracker maintains two counters instead of recomputing over all history.

A candidate in an ML interview is asked to classify an email from its words. The principled approach is Bayes' rule over the whole word sequence:

P(spam | w₁…wₙ) ∝ P(spam) · P(w₁…wₙ | spam)

The prior P(spam) is the base rate (spam ≈ 40% of a typical inbox). The likelihood requires the joint distribution of every word combination given the class — a table with |Vocabulary|ⁿ rows. Estimating it from data is impossible; that combinatorial explosion is the problem Naive Bayes dissolves with one deliberately false assumption: features are conditionally independent given the class label.

Under that assumption the joint factorizes into per-feature terms:

P(w₁…wₙ | spam) = P(w₁ | spam) · P(w₂ | spam) · … · P(wₙ | spam)

Now every feature's likelihood is a count you can estimate from training data, and the classifier is three lines of arithmetic. For text, the standard variant is Multinomial Naive Bayes with word counts; Gaussian Naive Bayes instead assumes each numeric feature is normally distributed within a class, with per-class mean and variance estimated from the data — the direct application of the distributions from the previous topics.

Two implementation details are not optional. First, work in log space: sum log-likelihoods instead of multiplying probabilities, because the product of 200 small probabilities underflows float64 and a single zero-count word (a word unseen in spam training data) would zero the entire product. Second, apply Laplace smoothing: add one pseudo-count to every word's count, which is exactly the Beta(1,1)/add-one updating from the conjugate chapter and eliminates both the zero-probability and underflow failure modes at once.

from sklearn.naive_bayes import MultinomialNB, GaussianNB
from sklearn.feature_extraction.text import CountVectorizer

texts = ["free money now", "claim your prize", "meeting at nine",
         "quarterly report attached", "project deadline friday"]
labels = [1, 1, 0, 0, 0]   # 1 = spam

X = CountVectorizer().fit_transform(texts)
clf = MultinomialNB(alpha=1.0).fit(X, labels)   # alpha=1.0 → Laplace smoothing
print(clf.predict(CountVectorizer().fit(texts).transform(["free prize now"])))

# Gaussian variant for numeric features: per-class N(mu_c, sigma_c^2)
GaussianNB().fit(X_numeric, labels)

Now the question every interviewer follows with: "the independence assumption is obviously false — words like 'free' and 'money' are correlated. Why does it work?" The answer has two parts, and you want both. First, Naive Bayes is estimating a ranking of class posteriors for an argmax decision, not calibrated probabilities; correlated features double-count evidence and distort the values, but the ordering of classes usually survives the distortion. Second, the alternative — estimating dense joint distributions over hundreds of features — needs vastly more data than a real training set has; the naive factorization is a bias-variance trade that wins in the small-to-medium-data regime. The failure mode to name: when features are strongly correlated and their combined evidence pushes the overconfident ranking past a rival class, logistic regression (which learns feature weights jointly and does not assume independence) pulls ahead. Prior strength matters here too: with 12 labeled training emails, the class prior and smoothed likelihoods carry most of the decision; with 12 million, the likelihood term dominates and the prior is nearly irrelevant — the same scarcity principle from the first two chapters, now inside a production classifier.