Topic models

What this chapter covers

Here is a problem with no labels in sight. You have a pile of documents — emails, abstracts, product reviews, whatever — and you suspect they cluster into a handful of themes. Some are about animals, some about computers, some about food. Nobody tagged them. Nobody told you how many themes there are or which words belong to which. You just have the words. Can an algorithm read the pile and hand back the themes, along with which document belongs to which?

That is what a topic model does, and this chapter builds the simplest one that works. We model each document as generated by a single hidden topic, and each topic as a bias over the vocabulary — a probability that leans hard on some words and away from others. Then we fit it with EM, the same E-step/M-step loop from the EM chapter, only now the hidden thing isn't which Gaussian drew a number but which topic wrote a document. We build the whole loop by hand in NumPy on a synthetic corpus with themes we planted ourselves, so we can check at the end whether EM found them. Then we hold it up against scikit-learn's LDA, the famous Bayesian descendant of the model we build.

The payoff to watch for is the animation. It starts every topic as flat noise, each word equally likely, no theme at all, and you get to watch the topics sharpen, iteration by iteration, until the animals topic has pulled all its mass onto dog and cat and tiger and the food topic onto pizza and bread. Themes emerging from noise, with nobody supervising. That is the whole idea of the page.

A bit of history

The lineage starts with a 1999 paper by Thomas Hofmann, "Probabilistic Latent Semantic Indexing." Text retrieval at the time leaned on latent semantic analysis, a purely linear-algebraic trick — take the document-term matrix and factor it with an SVD. It worked, but it had no probabilities in it and nothing you could call a generative story. Hofmann gave it one. His probabilistic LSA, pLSA, said a document is a mixture over latent topics and each topic is a distribution over words, and he fit it with EM. For the first time the latent "concepts" were honest probability distributions you could reason about. That model — a mixture of multinomials over words, fit by EM — is exactly the one we build here.

pLSA had a hole, and it was a deep one. It learned a topic mixture for each document in the training set but had no principled way to assign one to a document it had never seen, because the mixtures were parameters, not draws from anything. The fix came in 2003 from David Blei, Andrew Ng, and Michael Jordan in "Latent Dirichlet Allocation." They put a Dirichlet prior over the topic mixtures, which turned pLSA into a fully generative Bayesian model: now a document's mixture is itself sampled, so a fresh document is just another draw from the same story. LDA became one of the most cited papers in machine learning and the workhorse of topic modeling for a decade — every "what are the themes in this corpus" pipeline ran on some flavor of it. We build Hofmann's simpler model because it lays the mechanism bare; LDA is what you reach for in practice, and we compare against it at the end.

The intuition

Think about how you would write a document if you were a very simple author. First you pick a topic — say you roll a weighted die and it comes up "animals." Then you reach into the animals bag, which is stuffed mostly with the words dog, cat, tiger, wolf, rabbit and only a few of everything else, and you pull out words one at a time until the document is long enough. The document that comes out is a bag of words heavy on animal words. Someone reading it can't see the die roll or the bag — they see only the words — but the theme leaks through plainly in which words showed up and how often.

Now run that backwards. You are handed a stack of these bags and told nothing. You don't know the die weights, you don't know what's in each bag, and you don't know which bag each document was drawn from. But notice the same chicken-and-egg structure as every EM problem. If someone told you which topic wrote each document, filling the bags would be trivial — just pool all the animals documents and count their words. And if someone told you what's in each bag, guessing which bag wrote a given document would be trivial — score the document under each bag and pick the best fit. You have neither, so you bootstrap: start with random bags, use them to softly guess which topic wrote each document, refill the bags from those soft guesses, and loop. That loop is EM, and here is the structure it's trying to find, the three planted topics as a heatmap over the 15-word vocabulary:

Three bright blocks down the diagonal — each topic concentrates on its own five words and stays dim on the other ten. That block structure is the signal. EM's job is to recover it from the counts alone, without ever being shown this picture.

The math

Write each document as a row of word counts $x_i = (x_{i1}, \dots, x_{iV})$, where $x_{iw}$ is how many times word $w$ appeared and $V$ is the vocabulary size. We model the corpus with $K$ topics. Topic $k$ has a prior weight $\pi_k$ (the die, with $\sum_k \pi_k = 1$) and a word distribution $\phi_k = (\phi_{k1}, \dots, \phi_{kV})$ over the vocabulary (the bag, with $\sum_w \phi_{kw} = 1$). A single hidden label $z_i \in \{1, \dots, K\}$ says which topic wrote document $i$. The generative story — pick a topic, then draw its words — gives the probability of a document as a mixture over that hidden choice:

The inner product is the multinomial likelihood of the bag under topic $k$: each word contributes its topic probability raised to the number of times it occurred. The E-step asks the reverse question — given the current topics, what's the posterior probability that topic $k$ wrote document $i$? It's Bayes' rule, the topic's weighted likelihood over the total, and we call it the responsibility $W_{ik}$:

Each row $W_{i\cdot}$ sums to 1 — a soft assignment of document $i$ across the topics. The M-step then treats those responsibilities as fractional memberships and refits each topic by pooling responsibility-weighted word counts. The updated probability of word $w$ in topic $k$ is the weighted count of that word across all documents, normalized so the topic's distribution sums to one:

The prior update is just as plain — the new weight of topic $k$ is its average responsibility across the corpus, $\pi_k = \frac{1}{N}\sum_i W_{ik}$. And the quantity both steps quietly improve is the data log-likelihood, the log-probability of all the documents with the hidden topics summed out:

This is the number EM climbs, and by the same argument as the last chapter, where the E-step tightens a lower bound and the M-step raises it, it never goes down from one iteration to the next. One practical note that shapes the code: those products over the vocabulary are products of many small probabilities, which underflow to zero in floating point almost immediately. So we do everything in log-space, where the product becomes a sum $\sum_w x_{iw} \log \phi_{kw}$ and the normalization becomes a log-sum-exp. Every function in the build carries logs, not raw probabilities.

What it's good at, what it isn't

The appeal of a topic model is that it finds structure in text with no labels at all and hands it back in a form you can read. The output isn't an opaque embedding or a cluster id. It's a list of words per topic, ranked by probability, that a person can look at and name. "This topic is dog, cat, tiger, wolf: it's about animals." That interpretability is rare and valuable, and it comes for free from the multinomial-over-words structure. On top of that you inherit everything EM gives you: soft assignments, which are the honest answer when a document straddles two themes, and a monotone climb with no learning rate to babysit.

The costs are the ones every latent-variable model carries, sharpened by text. EM finds a local optimum, and with topics that's a real hazard. A bad start can merge two genuine themes into one blurry topic or split one theme across two, and the log-likelihood will happily plateau there looking converged. You have to tell it how many topics to look for, and picking $K$ is more art than science. And the whole model rests on bag-of-words: it throws away word order entirely, so "the dog bit the man" and "the man bit the dog" are identical to it, and it has no notion that "cat" and "kitten" are related unless both happen to co-occur with the same company. That last limitation is exactly the gap word embeddings and, later, large language models were built to close, which is where we land at the end.

The data

I built a synthetic corpus so we can grade the recovery against a truth we control. The vocabulary is 15 words in three blocks of five: an animals block (dog, cat, tiger, wolf, rabbit), a tech block (cpu, laptop, server, code, network), and a food block (pizza, bread, cheese, coffee, sugar). Each block is one planted topic. Every document is generated by the story from the intuition section: draw a topic from the prior $[0.40, 0.33, 0.27]$ — deliberately unequal, so the topics aren't the same size — then draw a bag of 24 to 60 words from that topic's multinomial, which puts heavy varied weight on its own five words and only a small leak on the other ten. The result is 150 documents, 6,038 word tokens in all, split 61 / 44 / 45 across the three topics. The true_topic of each document is recorded and then sealed: EM never sees it, and we open the envelope only at the end to score how well the recovered topics match the planted ones.

Here is the raw document-term count matrix — one row per document, one column per word, sorted so the documents cluster by their true topic. The color is the word count, and the block-diagonal glow is the same planted structure from before, except now it's the actual data EM will read, noise and all:

Read it top to bottom and you can see the three bands: the top group of documents lights up the leftmost five columns, the middle group the next five, the bottom group the last five. That's the ground truth. But remember EM sees these rows in no particular order and with no band labels — to the algorithm this is 150 rows of 15 numbers, and it has to discover that they fall into three groups and what each group is about.

Build it, one function at a time

Six functions, all in log-space. Where we start, the two EM steps, the score they climb, and the loop. Everything begins with the initialization, and for a topic model the init is the whole reason the animation is worth watching:

def init_params(X, k, rng):
    """The starting point: nearly flat topics and a uniform prior.

    Each topic-word distribution is uniform over the vocabulary jittered by a
    little noise, so no two topics are identical (EM can't split them if they
    are) but every one starts almost featureless. This is what makes the
    animation worth watching — the topics emerge from that flat noise.
    """
    d = X.shape[1]
    P = rng.uniform(0.95, 1.05, size=(k, d))
    P /= P.sum(axis=1, keepdims=True)
    log_pi = np.log(np.full(k, 1.0 / k))
    return np.log(P), log_pi

Each topic starts as very nearly the uniform distribution — every word about $1/15$ — jittered by a hair of noise. The jitter matters: if two topics started byte-identical, EM could never pull them apart, because the responsibilities would be identical too and the symmetry would never break. So we start them almost flat but not quite, and let EM find the asymmetry. The prior starts uniform. From this featureless start the topics have to grow their themes.

The E-step computes the responsibilities. For each document and each topic, the prior plus the log-likelihood of the bag under that topic, the log-joint, and normalizing each row over topics with a log-sum-exp turns it into the log posterior:

def e_step(X, log_P, log_pi):
    """E-step: the responsibility of each topic for each document.

    X is (N, d) integer counts, log_P is (k, d) log topic-word probabilities,
    log_pi is (k,) log topic priors. For document i and topic j the unnormalized
    log-joint is log(pi_j) + sum_w x_iw * log(P_jw) — the prior plus the
    log-likelihood of the bag under that topic's multinomial. Normalizing each
    row over topics (subtract the log-sum-exp) turns it into log W_ij, the log
    posterior probability that topic j generated document i. Returns (N, k).
    """
    log_joint = X @ log_P.T + log_pi[None, :]              # (N, k)
    log_evidence = logsumexp(log_joint, axis=1, keepdims=True)  # (N, 1)
    return log_joint - log_evidence

The whole thing is one matrix multiply: X @ log_P.T computes, for every document-topic pair at once, the sum $\sum_w x_{iw} \log \phi_{kw}$ — the log-multinomial-likelihood of every bag under every topic. Add the log prior, subtract the row-wise log-sum-exp, and you have log $W_{ik}$. The M-step then splits into two pieces. First the topics: weight every document's counts by the responsibility each topic takes for it, pool them, and normalize each topic to sum to one:

def m_step_topics(X, log_W, eps=1e-12):
    """M-step, part one: re-estimate each topic's word distribution.

    Weight every document's word counts by the responsibility the topic takes for
    it and pool them: E_jw = sum_i W_ij * x_iw (a tiny eps keeps the log finite
    for words no document credited to the topic). Normalizing each topic's row to
    sum to one gives the updated multinomial. Returns log_P, shape (k, d).
    """
    W = np.exp(log_W)                                     # (N, k) responsibilities
    weighted = W.T @ X + eps                              # (k, d) pooled counts
    P = weighted / weighted.sum(axis=1, keepdims=True)
    return np.log(P)

W.T @ X is the pooled, responsibility-weighted count of every word in every topic , exactly the numerator of the $\phi_{kw}$ update, and dividing each row by its sum does the normalization. The tiny eps keeps the log finite for a word no document credited to a topic. Second the priors, the average responsibility each topic carried, done in log-space:

def m_step_priors(log_W):
    """M-step, part two: re-estimate the topic priors.

    The updated prior for topic j is the average responsibility it carried across
    all documents: pi_j = (1/N) * sum_i W_ij. Done in log-space with logsumexp
    down the document axis. Returns log_pi, shape (k,).
    """
    N = log_W.shape[0]
    return logsumexp(log_W, axis=0) - np.log(N)           # (k,)

Now the score both steps are climbing, the data log-likelihood — marginalize the topic out of each document with a log-sum-exp, then sum over documents:

def log_likelihood(X, log_P, log_pi):
    """Data log-likelihood under the current parameters.

    Marginalize the latent topic out of each document — sum_j pi_j * P(doc|topic
    j) — in log-space with logsumexp, then sum over documents. This is the number
    EM is climbing; it must never go down from one full iteration to the next.
    """
    log_joint = X @ log_P.T + log_pi[None, :]              # (N, k)
    return float(logsumexp(log_joint, axis=1).sum())

This is the number that must never fall; if it ever does, there's a bug in one of the steps. Finally the loop that ties them together:

def em(X, k, rng, max_iter=100, tol=1e-4):
    """The EM loop: E-step, M-step, repeat until the log-likelihood stops rising.

    Records a snapshot every iteration — the full topic-word table, the priors,
    and the data log-likelihood — which the chapter replays frame by frame as the
    topics sharpen out of the flat start. Returns the hard document assignments
    (argmax responsibility), the recovered topics P (k, d), the priors pi (k,),
    the final log-likelihood, and the snapshot history.
    """
    log_P, log_pi = init_params(X, k, rng)
    history = [_snapshot(X, log_P, log_pi)]
    prev_ll = history[0]["ll"]
    for _ in range(max_iter):
        log_W = e_step(X, log_P, log_pi)          # responsibilities
        log_P = m_step_topics(X, log_W)           # re-estimate topics
        log_pi = m_step_priors(log_W)             # re-estimate priors
        snap = _snapshot(X, log_P, log_pi)
        history.append(snap)
        if snap["ll"] - prev_ll < tol:            # converged — no more to gain
            break
        prev_ll = snap["ll"]
    log_W = e_step(X, log_P, log_pi)
    labels = np.argmax(log_W, axis=1)
    return labels, np.exp(log_P), np.exp(log_pi), history[-1]["ll"], history

E-step, M-step topics, M-step priors, score, repeat until the log-likelihood stops rising by more than a whisker. It snapshots the full topic table, the priors, and the log-likelihood every iteration — that history is exactly what the animation replays. At the end it returns the hard document assignments (each document's most-responsible topic), the recovered topics, the priors, the final log-likelihood, and the history.

Watch it work

This is the centerpiece. Below is a real run of the em function above, one frame per iteration, starting from the nearly-flat random init. Three panels, one per recovered topic, arranged so that after convergence panel one reads as the animals topic, panel two as tech, panel three as food. Within each panel every bar is a word, and the bar is colored by the topic that truly planted that word — cyan for animals words, purple for tech words, orange for food words. So a panel that has locked onto the animals theme should end up all tall cyan bars and short everything else.

Watch the first frame: every panel is a flat picket fence of short bars in mixed colors, because at the random start no topic favors anything, each word sitting near the chance rate of $1/15$. Then press play and watch the bars sharpen. The colors sort themselves: one panel's cyan bars grow while its purple and orange bars shrink away, and the same for the other two panels with their colors. By the last frame each panel is a clean spike of one color — the topic has concentrated its whole probability mass onto its own five words. The caption tracks the iteration and the data log-likelihood climbing underneath it all.

The emergence is fast because the planted topics are well separated, and that's worth seeing plainly. At the random start the log-likelihood is −16,365 and every topic is noise. One E/M pass takes it to −15,015 and you can already see the bars starting to tilt toward their colors. The second pass jumps to −12,688, the frame where the themes become unmistakable, each panel now visibly dominated by one color. The third reaches −11,648, and the fourth confirms it: no more change, the topics have locked. Four iterations from noise to the answer.

Now find a frame where the log-likelihood goes down. There isn't one. The same guarantee as every EM, made visible again. Here is that climb on its own:

The familiar EM shape: a steep early climb, then a knee, then flat. It flattens after four iterations here rather than the twenty of the Gaussian chapter, because the topics are cleanly separated and the responsibilities snap to near-certain almost at once. On a messier corpus, with themes that share vocabulary, you'd see the long slow tail instead.

The full implementation

The whole file, no library, top to bottom. This is exactly what the animation ran:

"""Mixture-of-multinomials topic model, fit by EM — built from scratch.

Each document is a bag of words generated by ONE latent topic. A topic is a
multinomial distribution over the vocabulary. Given the count matrix X we don't
know which topic wrote which document, nor what the topics are — so we alternate:

  E-step: given the current topics and priors, compute the posterior probability
          that each document came from each topic (its responsibilities).
  M-step: given those responsibilities, re-estimate each topic's word
          distribution and the topic priors from responsibility-weighted counts.

Everything runs in log-space with scipy's logsumexp for numerical stability —
the per-word products become sums of log-probabilities and never underflow.
Pure NumPy for the model; the `# region:` markers are what the chapter includes.
"""

import numpy as np
from scipy.special import logsumexp


# region: e_step
def e_step(X, log_P, log_pi):
    """E-step: the responsibility of each topic for each document.

    X is (N, d) integer counts, log_P is (k, d) log topic-word probabilities,
    log_pi is (k,) log topic priors. For document i and topic j the unnormalized
    log-joint is log(pi_j) + sum_w x_iw * log(P_jw) — the prior plus the
    log-likelihood of the bag under that topic's multinomial. Normalizing each
    row over topics (subtract the log-sum-exp) turns it into log W_ij, the log
    posterior probability that topic j generated document i. Returns (N, k).
    """
    log_joint = X @ log_P.T + log_pi[None, :]              # (N, k)
    log_evidence = logsumexp(log_joint, axis=1, keepdims=True)  # (N, 1)
    return log_joint - log_evidence
# endregion


# region: log_likelihood
def log_likelihood(X, log_P, log_pi):
    """Data log-likelihood under the current parameters.

    Marginalize the latent topic out of each document — sum_j pi_j * P(doc|topic
    j) — in log-space with logsumexp, then sum over documents. This is the number
    EM is climbing; it must never go down from one full iteration to the next.
    """
    log_joint = X @ log_P.T + log_pi[None, :]              # (N, k)
    return float(logsumexp(log_joint, axis=1).sum())
# endregion


# region: m_step_topics
def m_step_topics(X, log_W, eps=1e-12):
    """M-step, part one: re-estimate each topic's word distribution.

    Weight every document's word counts by the responsibility the topic takes for
    it and pool them: E_jw = sum_i W_ij * x_iw (a tiny eps keeps the log finite
    for words no document credited to the topic). Normalizing each topic's row to
    sum to one gives the updated multinomial. Returns log_P, shape (k, d).
    """
    W = np.exp(log_W)                                     # (N, k) responsibilities
    weighted = W.T @ X + eps                              # (k, d) pooled counts
    P = weighted / weighted.sum(axis=1, keepdims=True)
    return np.log(P)
# endregion


# region: m_step_priors
def m_step_priors(log_W):
    """M-step, part two: re-estimate the topic priors.

    The updated prior for topic j is the average responsibility it carried across
    all documents: pi_j = (1/N) * sum_i W_ij. Done in log-space with logsumexp
    down the document axis. Returns log_pi, shape (k,).
    """
    N = log_W.shape[0]
    return logsumexp(log_W, axis=0) - np.log(N)           # (k,)
# endregion


# region: init_params
def init_params(X, k, rng):
    """The starting point: nearly flat topics and a uniform prior.

    Each topic-word distribution is uniform over the vocabulary jittered by a
    little noise, so no two topics are identical (EM can't split them if they
    are) but every one starts almost featureless. This is what makes the
    animation worth watching — the topics emerge from that flat noise.
    """
    d = X.shape[1]
    P = rng.uniform(0.95, 1.05, size=(k, d))
    P /= P.sum(axis=1, keepdims=True)
    log_pi = np.log(np.full(k, 1.0 / k))
    return np.log(P), log_pi
# endregion


# region: em
def em(X, k, rng, max_iter=100, tol=1e-4):
    """The EM loop: E-step, M-step, repeat until the log-likelihood stops rising.

    Records a snapshot every iteration — the full topic-word table, the priors,
    and the data log-likelihood — which the chapter replays frame by frame as the
    topics sharpen out of the flat start. Returns the hard document assignments
    (argmax responsibility), the recovered topics P (k, d), the priors pi (k,),
    the final log-likelihood, and the snapshot history.
    """
    log_P, log_pi = init_params(X, k, rng)
    history = [_snapshot(X, log_P, log_pi)]
    prev_ll = history[0]["ll"]
    for _ in range(max_iter):
        log_W = e_step(X, log_P, log_pi)          # responsibilities
        log_P = m_step_topics(X, log_W)           # re-estimate topics
        log_pi = m_step_priors(log_W)             # re-estimate priors
        snap = _snapshot(X, log_P, log_pi)
        history.append(snap)
        if snap["ll"] - prev_ll < tol:            # converged — no more to gain
            break
        prev_ll = snap["ll"]
    log_W = e_step(X, log_P, log_pi)
    labels = np.argmax(log_W, axis=1)
    return labels, np.exp(log_P), np.exp(log_pi), history[-1]["ll"], history
# endregion


def _snapshot(X, log_P, log_pi):
    """One frame of the run: the topics, the priors, and the log-likelihood."""
    return {
        "P": np.exp(log_P),
        "pi": np.exp(log_pi),
        "ll": log_likelihood(X, log_P, log_pi),
    }


def load_counts(path="../data/counts.csv"):
    """The committed document-term matrix plus the true topic per document.

    Returns (X, true_topic, vocab). The true_topic column is ground truth we
    NEVER fit on — EM sees only the counts. It exists to grade recovery.
    """
    import pandas as pd
    df = pd.read_csv(path)
    vocab = [c for c in df.columns if c != "true_topic"]
    X = df[vocab].to_numpy(float)
    return X, df["true_topic"].to_numpy(int), vocab

The library version

Nobody hand-rolls a topic model in production, and the standard tool isn't quite the model we built. scikit-learn's LatentDirichletAllocation is LDA, the Bayesian descendant from the history section. The difference is real and worth stating: our model assigns each document to a single latent topic, while LDA lets a document be a mixture of topics, each topic drawn from a Dirichlet prior, and fits it by variational inference rather than plain EM. It's the richer model. On a corpus like ours, where each document really was written by one topic, both should recover the same planted structure — but don't expect identical numbers, because they're optimizing different objectives under different assumptions.

def sklearn_lda(counts, k, seed=0):
    """Fit LDA with k topics on the integer count matrix.

    Returns the hard document assignments (the topic each document loads on most)
    and the normalized topic-word distributions (k, d) — the same two things we
    grade our own topics on.
    """
    lda = LatentDirichletAllocation(
        n_components=k, learning_method="batch", max_iter=50, random_state=seed
    )
    doc_topic = lda.fit_transform(counts)                 # (N, k) topic mixtures
    topic_word = lda.components_ / lda.components_.sum(axis=1, keepdims=True)
    labels = doc_topic.argmax(axis=1)
    return labels, topic_word

We ask it for the same two things we grade our own model on: a hard topic per document (whichever topic the document loads on most) and the normalized topic-word distributions. Everything else, the Dirichlet priors and the variational bound, is LDA's business, and we treat it as a black-box reference.

Scratch versus library

Both models see the same 150 documents and are asked for three topics. To score them we open the sealed envelope of true labels and use two measures. The first is the adjusted Rand index between each model's hard document assignment and the true topics — agreement between two labelings, corrected for chance, where 1.0 is perfect. The second is how close each recovered topic-word distribution is to the one that planted it, measured by cosine similarity after matching recovered topics to planted ones (both models return topics in arbitrary order, so we pair them up by maximum similarity first).

Both models score a perfect adjusted Rand index of 1.00: every document sorted into the right topic, all 150 of them. And both recover the word distributions almost exactly: mean cosine 0.998 for our from-scratch EM, 0.998 for LDA. Our model even recovers the unequal prior it was never told about, landing on topic weights of 0.41, 0.29, 0.30 against the planted 0.40, 0.33, 0.27. On this clean corpus the two models are indistinguishable at the level of who's right, which is the point of a clean corpus. It confirms both implementations are correct before you trust either on data where the answer is unknown. To see the recovery word by word, here is our model's topics laid over the planted truth, one panel per topic:

The recovered bars sit right on top of the planted ones in every panel. EM, given nothing but 150 bags of words, rebuilt the three topics we hid inside them: the right words, the right weights, the right documents in each, with no supervision at any step. That's the whole trick working.

Takeaways

A topic model is EM wearing a text costume. Swap the Gaussians of the last chapter for multinomials over words and the exact same E-step/M-step loop discovers themes in an unlabeled corpus: the E-step softly guesses which topic wrote each document, the M-step refills the topics from those guesses, and the data log-likelihood climbs monotonically to a stationary point. Build it once and you have the mechanism behind pLSA and, with a Bayesian prior bolted on top, LDA, the model that ran corpus exploration for a decade. The recovery on our planted corpus was perfect because we made it clean; on real text you fight the usual EM demons, local optima that merge or split genuine themes, and the eternal question of how many topics to ask for.

The honest place to end is on where this sits now. The two families didn't stay even. pLSA is the maximum-likelihood point estimate, fast and simple, exactly what we built, but it's really only a training-set model and it overfits without a prior. LDA answered that with full Bayesian treatment, and for years it was the default. Then the ground shifted. Both models share the one crippling assumption we flagged: bag-of-words, which knows nothing about word order or that two words might mean the same thing. Word embeddings started closing that gap, and large language models closed it completely. Ask a modern model to cluster or summarize a corpus by theme and it reads the documents, meaning and syntax and all, in a way no multinomial ever could. Classical topic models haven't vanished; they still earn their keep where you need speed, interpretable word lists, or a model small enough to reason about, and short-text clustering pipelines still run them. But the frontier of "what is this pile of text about" moved to representations that understand words rather than just count them. Knowing the EM version is how you understand what those representations replaced, and why the replacement was such a leap.