Naive Bayes
Bayes' theorem, conditional independence, MLE estimates, and Laplace smoothing.
On this page
Beginner: Naive Bayes classifies by playing out a little story for each possible label and seeing which story fits the evidence best. For an email, the story goes: "if this were spam, how likely would I be to see exactly these words? And if it were not spam, how likely would I be to see exactly these words?" Whichever story makes the observed words more likely — weighted by how common that label is overall — wins. That comparison, done properly, is nothing more than Bayes' theorem. The "naive" part is a simplifying assumption bolted on to make the comparison tractable: pretend each word's presence is generated independently of every other word, once you already know the label. That's obviously not true of real language — but as the derivations below show, it turns out not to matter nearly as much as it should.
Intermediate: here is the comparison worked through with real numbers. Suppose an email contains the words "free", "win", and "meeting", and training data says: P(free|spam)=0.34, P(free|ham)=0.02; P(win|spam)=0.32, P(win|ham)=0.01; P(meeting|spam)=0.03, P(meeting|ham)=0.30. With priors P(spam)=0.4, P(ham)=0.6, the naive factorization gives: P(spam)·0.34·0.32·0.03 ≈ 0.00131 versus P(ham)·0.02·0.01·0.30 ≈ 0.000036. Spam wins by more than 30-to-1, even though "meeting" alone is ten times more typical of ham than of spam — because "free" and "win" together supply far more evidence than "meeting" can outweigh. This is exactly what the naive factorization does: it treats each word as an independent vote, and lets the votes accumulate multiplicatively.
Advanced: "features are conditionally independent given the class" says nothing about what distribution each individual factor P(xⱼ|y=k) should follow — that is a second, separate modeling choice, and it is exactly what distinguishes the three variants this lesson covers. Continuous, real-valued features (pixel intensities, sensor readings, a patient's blood pressure) get a per-class, per-feature Gaussian — Gaussian Naive Bayes. Nonnegative integer counts (how many times each vocabulary word appears in a document) get a per-class multinomial — Multinomial Naive Bayes, the classic bag-of-words model. Pure binary presence/absence features (did this word appear at all, ignoring how many times) get a per-class Bernoulli — Bernoulli Naive Bayes. All three share the identical Bayes'-rule combination step and the identical independence factorization; they differ only in which one-dimensional family is plugged in for each P(xⱼ|y=k) factor, each fit by its own maximum-likelihood estimate (section 2.1.6), each vulnerable to the same zero-probability failure that Laplace smoothing fixes, and each implemented, in practice, entirely in log-space.
Bayes' theorem, applied to classification:
the naive conditional-independence factorization of the class-conditional density, giving the decision rule actually used at prediction time:
and the log-space form every real implementation actually computes:
Every symbol here is derived from scratch below: where Bayes' theorem itself comes from, why the product factorizes at all, how each variant's P(xⱼ|y=k) is actually fit, why it needs smoothing, and why the sum replaces the product in every real implementation.
Start from the definition of conditional probability, applied twice to the same joint event A \cap B:
The second equation rearranges to P(A \cap B) = P(B\mid A)\,P(A). Substitute that directly into the numerator of the first:
That is Bayes' theorem — a purely algebraic consequence of conditional probability, no extra assumption required. Now set A to "the label is class k" and B to "the observed features are x":
Three names are worth fixing: P(y=k) is the prior — how common class k is before looking at any features at all; P(x|y=k) is the class-conditional likelihood — how plausible this exact feature vector is under class k's story; P(y=k|x) is the posterior — what we actually want, updated belief in class k after seeing x. The denominator P(x), the evidence, expands by the law of total probability to \sum_j P(x\mid y=j)P(y=j) — but note it does not depend on k at all. Since classification only needs \arg\max_k P(y=k\mid x), and dividing every candidate by the same positive constant never changes which one is largest, the denominator can simply be dropped for the purpose of choosing a label:
Where this is used: this proportionality, not the full normalized posterior, is what every generative classifier in this module (section 2.5's overview) actually computes at prediction time — compute an unnormalized score per class, and take the argmax. Recovering the actual normalized probability, when it's needed, is what the log-sum-exp derivation at the end of this lesson is for.
P(x|y=k) is a joint density over all d features at once. In general that joint requires estimating an object with exponentially many parameters — for example, d binary features have 2^d − 1 free joint probabilities to estimate per class, which is hopeless with any realistic amount of training data once d reaches even a few dozen. The naive assumption sidesteps this entirely by declaring the features conditionally independent given the class:
Read this assumption precisely — it is easy to misstate. It says xᵢ and xⱼ are independent within a fixed class k; it says nothing about whether they are independent overall, unconditionally, across the whole population. Those are entirely different claims. The words "prescription" and "pharmacy" are strongly correlated across all emails taken together — emails that mention one tend to mention the other — but the model only needs the weaker, conditional claim: once you already know an email is spam, does seeing "prescription" change how likely "pharmacy" is within that class? The naive model says no. Marginally the two words are still correlated in the data as a whole (both are more common in spam than in ham, so knowing one is present raises your belief the email is spam, which in turn raises your belief in the other) — the independence claim is scoped strictly to "given the class label is already fixed."
Geometrically, this is exactly why the diagram below shows an axis-aligned classification region rather than a tilted one: each class's story is a product of d independent one-dimensional distributions, one per feature, with no cross term between any pair of features ever estimated. A class's data cloud is modeled as though its contours ran parallel to the coordinate axes, even when the real data is visibly correlated and elliptical within that class. Section 2.5.2's Gaussian Discriminant Analysis is precisely what you get by dropping this restriction and estimating a full covariance matrix instead.
The assumption is almost always literally false, and yet Naive Bayes is a genuinely strong classifier in practice. The resolution is that classification only ever needs \arg\max_k P(y=k\mid x) to be right — it never needs the absolute value P(x|y=k) to be an accurate density estimate. When features are correlated within a class, the independence assumption multiplies in the same correlated evidence once per correlated feature instead of once total — every class's score gets distorted, and typically distorted in a similar direction (favoring whichever class the real evidence favors, just more emphatically than it should). As long as that distortion doesn't flip the *ordering* between classes more often than it preserves it, the decision boundary the naive model draws can still track the true optimal boundary closely, even though the numeric value of P(x|y=k) it reports is a poor estimate of the true class-conditional density. The Expert note below returns to exactly this point from the calibration side.
For a continuous feature xⱼ, Gaussian NB assumes x_j \mid y=k \sim \mathcal{N}(\mu_{jk}, \sigma_{jk}^2) — its own mean and variance for every (feature, class) pair. Collect the nₖ training examples of class k; their feature-j values x_{1j},\ldots,x_{n_k j} are treated as i.i.d. draws from that Gaussian, so the log-likelihood is:
Differentiate with respect to μ and set to zero:
— the ordinary sample mean of feature j, restricted to class k's examples. Now substitute μ̂ⱼₖ back in and differentiate with respect to u = σ² (treating it as a single variable makes the algebra cleaner):
Multiplying through by 2u²/nₖ and solving:
— the ordinary sample variance, again restricted to class k. Both estimates are exactly what you'd compute by hand: split the training data by class, then take the per-feature mean and variance within each split, one pass through the data, no iteration required. The class priors get the same treatment as any categorical MLE: \hat\pi_k = n_k / n, the empirical class frequency.
Where this is used: exactly these closed-form estimates are what fit the diagram below — no gradient descent, no iteration, just per-class sample statistics — and exactly this factorized Gaussian is why the diagram's classification region has no tilt: each feature's contribution to the log-posterior is a separate parabola in that one feature alone.
Multinomial NB is built for count data — the canonical example is bag-of-words, where xⱼ is "how many times word j appears in this document." Within class k, each vocabulary word j has a probability θⱼₖ of being the next word drawn, with \sum_j \theta_{jk} = 1. Let Nⱼₖ be the total number of times word j occurs across all class-k training documents; the multinomial log-likelihood (ignoring the multinomial coefficient, which doesn't involve θ) is \sum_j N_{jk}\log\theta_{jk}. Maximize this subject to the simplex constraint with a Lagrange multiplier λ:
Setting \partial\mathcal L/\partial\theta_{jk} = N_{jk}/\theta_{jk} - \lambda = 0 gives \theta_{jk} = N_{jk}/\lambda; summing both sides over j and using \sum_j\theta_{jk}=1 pins down \lambda = \sum_j N_{jk} = N_k, the total word count in class k. So:
— the raw relative frequency of word j among all word occurrences in class k's documents. This is the textbook MLE for a multinomial parameter, and it's exactly what the Laplace-smoothing derivation below starts from.
Bernoulli NB is built for a genuinely different kind of feature: pure binary presence/absence — "did word j appear in this document at all", regardless of how many times. Within class k, feature j has its own parameter pⱼₖ = P(xⱼ=1|y=k). With mₖ documents of class k and Dⱼₖ of them containing word j at least once, the likelihood over these mₖ independent Bernoulli trials is:
Log, differentiate, set to zero exactly as for any Bernoulli MLE:
— the fraction of class-k documents in which word j appears at least once. So far this looks like Multinomial NB's estimate with a relabeling, but the likelihood a document actually gets scored with is structurally different, and this is the point the content scope for this lesson insists on getting right. Bernoulli NB's per-document likelihood explicitly multiplies in a term for every feature in the vocabulary, present or not:
A word the class strongly expects (pⱼₖ close to 1) that is absent from this particular document contributes a small factor (1-pⱼₖ) — actively penalizing that class, not just failing to reward it. Multinomial NB has no such term at all: its "event" is a sequence of drawn word-tokens, and words that never get drawn simply never enter the product — there is no notion of "expected but missing" to penalize. Practically, this makes Bernoulli NB noticeably more sensitive to short documents with a lot of absent-but-expected vocabulary (and is why it's the traditional choice for short texts), while Multinomial NB, which does care how many times each word appears, is the traditional choice for longer, count-sensitive documents.
Start from the raw multinomial MLE just derived, \hat\theta_{jk} = N_{jk}/N_k. If word j never once occurred in class k's training documents, Nⱼₖ = 0, so θ̂ⱼₖ = 0 — not approximately small, exactly zero. Now look at what that zero does to a document's score. Say a document's naive factorization is a product of five per-word estimates, four of them strong evidence for class k:
Multiplication by zero is absorbing — the whole product is exactly zero regardless of how overwhelmingly the other four factors favor class k. One missing training example for one word is enough to make class k mathematically impossible for every future document containing that word, no matter what else is in it. This is the precise failure the diagram below demonstrates live.
The fix is a MAP estimate (section 2.1.6's MLE-vs-MAP framing) rather than a raw MLE. Place a symmetric Dirichlet prior over the vector (\theta_{1k},\ldots,\theta_{Vk}) — the conjugate prior for a multinomial — with density p(\theta) \propto \prod_j \theta_{jk}^{\alpha} for some α ≥ 0. The posterior is proportional to likelihood times prior:
which has exactly the same algebraic shape as the raw likelihood, with Nⱼₖ replaced by Nⱼₖ + α. Maximizing it — the same Lagrange-multiplier argument as Derivation 4, term for term — gives:
where V is the vocabulary size (so that the estimates still sum to 1 across j). α = 1 is the classic Laplace (add-one) smoothing; general α > 0 is sometimes called Lidstone smoothing, but this lesson (and the diagram below) just calls it add-α smoothing throughout. The intuition is exact, not approximate: it's as if α extra "pseudo-occurrences" of every word had been observed in every class before the real data was ever looked at — a mild prior belief that no word is truly impossible. As α → 0 this recovers the raw MLE (and its zero-probability failure mode); larger α pulls every estimate toward the uniform 1/V, trading a little fit to the observed counts for the guarantee that no estimate is ever exactly zero. The same Beta-prior argument, conjugate to a Bernoulli/Binomial likelihood instead of a multinomial one, gives Bernoulli NB's smoothed estimate \hat p_{jk} = (D_{jk}+\alpha)/(m_k + 2\alpha).
Return to the concrete failure above with α = 1: the zero factor becomes (0+1)/(N_k+V) — small, but strictly positive — so the full product is now a tiny but nonzero number instead of an absolute veto. The other four words' evidence is free to be weighed again, and can still carry class k to the winning posterior if it's strong enough. Smoothing doesn't manufacture evidence for the unseen word; it just stops that one missing training example from silently overriding every other feature in the document.
Even with smoothing guaranteeing every factor is strictly positive, multiplying hundreds of numbers each well below 1 is numerically dangerous. IEEE 754 double-precision floats can represent positive numbers no smaller than about 5\times10^{-324} before rounding to exactly 0.0. A product of just 220 factors averaging 0.03 is already around 10^{-334} — past that floor. This isn't a "very small number" any more; it's identically zero in the machine's arithmetic, and every bit of information about how much more (or less) likely one class was than another is gone. The diagram below shows exactly this happening, term by term.
The fix is to never form the raw product at all. Take the log of the score being maximized before multiplying anything:
Because log turns products into sums, and because log is strictly increasing (so it never changes which k attains the maximum), this sum is exactly as good for classification as the raw product — the argmax is identical either way. But sums behave completely differently in floating point than products of fractions do: adding numbers like -3.1, -2.4, -3.6, \ldots just keeps the running total more negative, with no absorbing floor at zero and no loss of precision as more terms are added, no matter how many hundreds of features there are. This is the entire reason a real implementation computes \log P(y=k) + \sum_j \log P(x_j\mid y=k) per class and takes the argmax over those sums, never forming P(x|y=k) itself.
Log-space is enough to pick the winning class, but sometimes the actual calibrated probabilities P(y=k|x) are wanted too (for a risk score, a confidence threshold, and so on) — and that requires dividing by the evidence P(x) = \sum_j P(x\mid y=j)P(y=j), which in log-space means computing \log\sum_k \exp(\text{score}_k). Naively exponentiating each log-score can itself overflow (the opposite numerical failure — huge numbers instead of vanishing ones), which is precisely what section 1.21's log-sum-exp trick — subtract off the largest score before exponentiating, add it back after taking the log — was built to fix. So the full pipeline in practice is: compute each class's log-score as a sum (this derivation), then, only if normalized probabilities are actually needed, run log-sum-exp (section 1.21) on those log-scores to safely turn them back into a proper distribution over classes.
Where this is used: this is why every from-scratch implementation in the code tabs below accumulates Math.log / std::log terms in a running sum rather than ever calling a running product — it is not a style preference, it is the difference between a classifier that works on realistic feature counts and one that silently returns "0 vs. 0" past a few hundred features.
Each class's Gaussian is fit per feature, independently -- no cross term between x1 and x2 is ever estimated, exactly as Derivation 2 requires. Drag either slider to move class B's feature-1 mean or variance away from its fitted value and watch the classification region and the bottom-panel density curve reshape live.
A tiny bag-of-words example where 'cryptocurrency' never occurred in ham training documents. Toggle smoothing off (alpha = 0) to watch P(doc | ham) collapse to exactly 0 despite three strongly ham-typical words in the same document, then raise alpha and watch the correct prediction come back.
260 realistic per-feature likelihoods multiplied in one at a time for two competing hypotheses. The raw running products underflow to exactly 0.000000 partway through and become an uninformative tie; the running log-sums (plotted above) stay comparable for every single term. Scrub the slider or hit replay.
All three from-scratch implementations follow the same shape: fit per-class parameters with the closed-form MLE (or MAP, once smoothing is added) derived above, then predict by summing log-scores and taking the argmax — never by forming a raw product. The library tab fits the same three models via scikit-learn's GaussianNB, MultinomialNB, and BernoulliNB, which store fitted parameters as theta_/var_ and feature_log_prob_ — already in log-space, for exactly the reason Derivation 6 gives.
#include <cmath>
#include <iostream>
#include <vector>
// Gaussian Naive Bayes -- the closed-form MLE from Derivation 3, and the log-space
// scoring rule from Derivation 6.
struct GaussianClassStats {
std::vector<double> mean, variance;
double prior;
};
std::vector<GaussianClassStats> fitGaussianNB(const std::vector<std::vector<double>>& X,
const std::vector<int>& y, int nClasses) {
int d = X[0].size();
std::vector<GaussianClassStats> stats(nClasses, {std::vector<double>(d, 0.0),
std::vector<double>(d, 0.0), 0.0});
std::vector<int> counts(nClasses, 0);
for (size_t i = 0; i < X.size(); ++i) {
int k = y[i];
counts[k]++;
for (int j = 0; j < d; ++j) stats[k].mean[j] += X[i][j];
}
for (int k = 0; k < nClasses; ++k) {
for (int j = 0; j < d; ++j) stats[k].mean[j] /= counts[k];
stats[k].prior = static_cast<double>(counts[k]) / X.size();
}
for (size_t i = 0; i < X.size(); ++i) {
int k = y[i];
for (int j = 0; j < d; ++j) {
double diff = X[i][j] - stats[k].mean[j];
stats[k].variance[j] += diff * diff;
}
}
for (int k = 0; k < nClasses; ++k)
for (int j = 0; j < d; ++j) stats[k].variance[j] = stats[k].variance[j] / counts[k] + 1e-9;
return stats;
}
double gaussianLogPdf(double x, double mean, double var) {
return -0.5 * std::log(2 * M_PI * var) - (x - mean) * (x - mean) / (2 * var);
}
int predictGaussianNB(const std::vector<double>& x, const std::vector<GaussianClassStats>& stats) {
int best = 0;
double bestScore = -1e300;
for (size_t k = 0; k < stats.size(); ++k) {
// Log-space sum (Derivation 6), never a raw product.
double score = std::log(stats[k].prior);
for (size_t j = 0; j < x.size(); ++j)
score += gaussianLogPdf(x[j], stats[k].mean[j], stats[k].variance[j]);
if (score > bestScore) { bestScore = score; best = static_cast<int>(k); }
}
return best;
}
int main() {
std::vector<std::vector<double>> X = {{-2.1, -1.0}, {-1.4, -0.5}, {-1.8, -1.3}, {-0.9, -0.4},
{1.3, 1.1}, {1.9, 1.4}, {0.8, 1.7}, {2.0, 0.9}};
std::vector<int> y = {0, 0, 0, 0, 1, 1, 1, 1};
auto stats = fitGaussianNB(X, y, 2);
std::vector<double> query = {1.5, 1.0};
std::cout << "predicted class: " << predictGaussianNB(query, stats) << "\n";
return 0;
}- Spam filtering. Naive Bayes' single most famous application, and still a common baseline layer — Multinomial or Bernoulli NB over word features, exactly the setup in the smoothing diagram above.
- Sentiment analysis. Classifying reviews or social posts as positive/negative/neutral from bag-of-words or presence/absence word features, valued for how cheaply it retrains as vocabulary shifts.
- Document and topic classification. Routing news articles, support tickets, or emails into categories from their word-frequency profile — Multinomial NB's core use case at scale.
- Medical diagnosis screening. Estimating disease probability from a set of symptoms or test results, explicitly leaning on the (approximate) assumption that symptoms are independent given the diagnosis — usually false in detail, but often close enough for a fast first-pass screen, echoing the "wrong but useful" theme running through this whole module.
- Real-time and resource-constrained classification. Training is one pass over the data to compute counts, means, and variances (no iterative optimization at all), and the fitted model is just a small table of per-class, per-feature parameters — making Naive Bayes a natural fit for on-device or streaming settings where both training time and memory footprint matter.
- Skipping Laplace (or add-
α) smoothing and being surprised that a well-trained model occasionally makes wildly wrong, overconfident predictions — as Derivation 5 shows, a single unseen feature value for one class doesn't just weaken that class's score, it zeroes it out completely, regardless of every other feature. - Feeding continuous, unnormalized, or negative-valued features to Multinomial NB. Its likelihood is only defined for nonnegative counts; standardized features (which can be negative) or raw real-valued measurements need Gaussian NB instead, or an explicit transformation into counts first.
- Treating the conditional-independence assumption as something that must hold for the model to be trustworthy. As Derivation 2 argues, it is almost always false and the model is frequently excellent anyway — but it's worth checking when it might actually be costly: features that are strongly correlated within a class (not just overall) are the case where the naive factorization's distortion is largest.
- Implementing the posterior as a literal running product instead of a running sum of logs. It works perfectly on toy examples with a handful of features and silently starts returning ties or wrong answers once the feature count grows into the hundreds — exactly the underflow Derivation 6 and the third diagram walk through.
Going deeper
Good decisions, bad probabilities. Naive Bayes is notorious for producing classification decisions that are quite good while its predicted probabilities P(y=k|x) are badly miscalibrated — frequently pushed to extremes near 0 or 1 even when the true confidence should be much more moderate. Derivation 2 explains why: when features are positively correlated within a class but modeled as independent, each correlated feature contributes its own full log-likelihood term as though it were independent new evidence, effectively double- (or triple-, or more-) counting the same underlying signal. That inflation happens on whichever side the real evidence already favors, so the argmax — the classification decision — is often preserved even as the reported probability gets pushed further toward an extreme than the data actually justifies. If a downstream system only needs the label, this rarely matters. If it needs the probability itself — to rank candidates by risk, or to set a decision threshold other than the default 0.5 — Naive Bayes' raw output is a poor input, and it typically needs post-hoc recalibration (Platt scaling or isotonic regression) rather than being trusted directly.
A second thing worth knowing: run the naive-Bayes generative recipe with Bernoulli or multinomial class-conditionals through Bayes' rule and simplify, and the resulting posterior P(y|x) comes out in exactly the same sigmoid/softmax functional form that logistic and softmax regression (sections 2.4.1, 2.4.2) assume directly — the same phenomenon module-5's overview flags for Gaussian Discriminant Analysis and logistic regression, and the same exponential-family machinery underlying GLMs (section 2.3.9). Naive Bayes and logistic regression are, in a real sense, the generative/discriminative pair for discrete features, exactly as GDA and logistic regression are that pair for continuous Gaussian ones — section 2.5.4 makes this precise for the whole module.
Naive Bayes assumes features are conditionally independent given the class label, an assumption that's almost always false for real data. Why can the resulting classifier still perform well despite this?
Classification only ever needs argmax_k P(y=k|x) to pick the right class -- it never needs the absolute value P(x|y=k) to be an accurate density estimate. When features are correlated within a class but modeled as independent, the naive factorization multiplies in the same correlated evidence once per correlated feature instead of once total, distorting every class's score -- but that distortion tends to push in the same direction the real evidence already favors, inflating whichever class's score was already larger rather than flipping the ranking. As long as the distortion doesn't flip which class has the larger score more often than it preserves it, the decision boundary can stay close to optimal even though the reported P(x|y=k) values and the resulting posterior probabilities are a poor, overconfident estimate of the truth -- exactly the calibration-vs-decision distinction the Expert note draws out.
Every generative classifier in this module runs the same recipe: pick a class-conditional density for P(x|y=k), fit it by maximum likelihood (or, with a prior, by MAP), and combine it with the prior through Bayes' theorem. Naive Bayes made the simplest possible choice for that density — features factorize completely, one independent 1-D distribution per feature — which is exactly what bought its computational simplicity and its small-data robustness, and exactly what produced its axis-aligned decision regions. Section 2.5.2, Gaussian Discriminant Analysis, keeps this identical Bayes'-rule recipe and relaxes only that one assumption: instead of independent per-feature Gaussians, each class gets a full multivariate Gaussian with a genuine covariance structure, trading Naive Bayes' simplicity and data-efficiency for boundaries that can finally tilt to match correlated features directly.