U5.4 Comparing the Models

The three previous pages fitted LDA, NMF and BERTopic one at a time. This page puts them side by side on the same posts, with the same number of topics and two seeds each, sums the differences up in one table, and ends with the pitfalls that catch most topic-modeling projects.

import numpy as np
import pandas as pd
from sklearn.feature_extraction.text import CountVectorizer

def top_words(components, vocab, n=4):
    return [' '.join(vocab[row.argsort()[::-1][:n]]) for row in components]

from sklearn.datasets import fetch_20newsgroups

forums = ['comp.graphics', 'rec.sport.hockey', 'sci.space', 'soc.religion.christian']
ng = fetch_20newsgroups(subset='train', categories=forums, random_state=0,
                        remove=('headers', 'footers', 'quotes'))   # keep only the body text
long_enough = np.array([len(text.split()) >= 20 for text in ng.data])
docs = [text for text, keep in zip(ng.data, long_enough) if keep]
y = ng.target[long_enough]                                         # the forum: used only to check

from sklearn.metrics import adjusted_rand_score

settings = dict(stop_words='english', min_df=5, max_df=0.5,
                token_pattern=r'(?u)\b[a-zA-Z]{3,}\b')             # words of 3+ letters, no numbers
cv = CountVectorizer(**settings)
X = cv.fit_transform(docs)
vocab = cv.get_feature_names_out()

from sklearn.decomposition import LatentDirichletAllocation

from sklearn.decomposition import NMF
from sklearn.feature_extraction.text import TfidfVectorizer

tv = TfidfVectorizer(**settings)
T = tv.fit_transform(docs)                                  # same vocabulary as X, TF-IDF weights

B = CountVectorizer(stop_words='english', token_pattern=r'(?u)\b[a-zA-Z]{3,}\b',
                    binary=True).fit(docs)                 # reference: which post contains which word
Xb = B.transform(docs).tocsc()
col = B.vocabulary_
D = Xb.shape[0]

def n_docs(*ws):                                           # number of posts containing all the words
    hit = np.ones(D, dtype=bool)
    for w in ws:
        hit &= Xb[:, col[w]].toarray().ravel() > 0
    return int(hit.sum())

from itertools import combinations

def npmi(w1, w2):
    p1, p2, p12 = n_docs(w1) / D, n_docs(w2) / D, n_docs(w1, w2) / D
    return -1.0 if p12 == 0 else np.log(p12 / (p1 * p2)) / -np.log(p12)

def coherence(topics, n=10):           # mean NPMI over all pairs of each topic's top-n words
    return np.mean([np.mean([npmi(a, b) for a, b in combinations(t.split()[:n], 2)])
                    for t in topics])

def diversity(topics, n=10):           # share of distinct words among all topics' top-n words
    ws = [w for t in topics for w in t.split()[:n]]
    return len(set(ws)) / len(ws)

from sklearn.cluster import KMeans
from sklearn.metrics.pairwise import cosine_similarity
from sentence_transformers import SentenceTransformer

encoder = SentenceTransformer('all-MiniLM-L6-v2')

from sklearn.neighbors import NearestNeighbors

E = encoder.encode(docs, show_progress_bar=False)       # stage 1: one 384-d vector per post

from umap import UMAP

umap_model = UMAP(n_neighbors=15, n_components=5, min_dist=0.0, metric='cosine',
                  random_state=42)                      # BERTopic's defaults, plus a seed
Z = umap_model.fit_transform(E)                         # stage 2: 384 -> 5 dimensions

from hdbscan import HDBSCAN

hdbscan_model = HDBSCAN(min_cluster_size=10, metric='euclidean',
                        cluster_selection_method='eom', prediction_data=True)
labels = hdbscan_model.fit_predict(Z)                   # stage 3: dense regions; -1 = outlier

cv_topics = CountVectorizer(stop_words='english', token_pattern=r'(?u)\b[a-zA-Z]{3,}\b')

from bertopic import BERTopic

topic_model = BERTopic(embedding_model=encoder, umap_model=umap_model,
                       hdbscan_model=hdbscan_model, vectorizer_model=cv_topics)
topics = np.array(topic_model.fit_transform(docs, embeddings=E)[0])   # reuse stage-1 vectors

from bertopic.representation import KeyBERTInspired
topic_model.update_topics(docs, vectorizer_model=cv_topics, representation_model=KeyBERTInspired())

LDA vs NMF vs BERTopic, measured

Stage 5. The same 2,155 posts, the same number of topics for all three (24, the number BERTopic found), two seeds each. The scores are those of the NMF page, plus the homogeneity of the topics with respect to the forums (homogeneity_score: 1.0 when every topic contains posts of a single forum; it does not penalize splitting a forum into several topics). For BERTopic, the outliers are left out of the homogeneity:

from sklearn.metrics import homogeneity_score

rows, assign = [], {}
n_topics = int(topics.max()) + 1
for seed in [42, 43]:
    lda_k = LatentDirichletAllocation(n_components=n_topics, random_state=seed).fit(X)
    nmf_k = NMF(n_components=n_topics, init='nndsvda', random_state=seed, max_iter=500).fit(T)
    tm = BERTopic(embedding_model=encoder, vectorizer_model=cv_topics,
                  umap_model=UMAP(n_neighbors=15, n_components=5, min_dist=0.0,
                                  metric='cosine', random_state=seed),
                  hdbscan_model=HDBSCAN(min_cluster_size=10, prediction_data=True))
    bt_label = np.array(tm.fit_transform(docs, embeddings=E)[0])
    for name, lab, tops in [
            ('LDA', lda_k.transform(X).argmax(axis=1), top_words(lda_k.components_, vocab, 10)),
            ('NMF', nmf_k.transform(T).argmax(axis=1), top_words(nmf_k.components_, vocab, 10)),
            ('BERTopic', bt_label, [' '.join(w for w, _ in tm.get_topic(t))
                                    for t in range(bt_label.max() + 1)])]:
        assign[name, seed] = lab
        inside = lab >= 0
        rows.append(dict(model=name, seed=seed, topics=len(tops), outliers=f"{1 - inside.mean():.0%}",
                         coherence=round(coherence(tops), 3), diversity=round(diversity(tops), 2),
                         homogeneity=round(homogeneity_score(y[inside], lab[inside]), 3)))
print(pd.DataFrame(rows).to_string(index=False))
for name in ['LDA', 'NMF', 'BERTopic']:
    print(f"{name:9s} seed 42 vs seed 43, ARI {adjusted_rand_score(assign[name, 42], assign[name, 43]):.3f}")

   model  seed  topics outliers  coherence  diversity  homogeneity
     LDA    42      24       0%      0.265       0.72        0.526
     NMF    42      24       0%      0.391       0.92        0.646
BERTopic    42      24      20%      0.251       0.88        0.927
     LDA    43      24       0%      0.276       0.70        0.574
     NMF    43      24       0%      0.407       0.91        0.644
BERTopic    43      16      10%      0.295       0.82        0.662
LDA       seed 42 vs seed 43, ARI 0.164
NMF       seed 42 vs seed 43, ARI 0.736
BERTopic  seed 42 vs seed 43, ARI 0.324

How to read it:

• Homogeneity. With seed 42, BERTopic's topics are by far the purest (0.927 against 0.646 for NMF and 0.526 for LDA), but that is measured on the 80% of posts it kept, while NMF and LDA have to place every post; the pitfalls section measures what happens when BERTopic's outliers are forced back in. With seed 43, BERTopic finds 16 broader clusters and its homogeneity falls to 0.662, the level of NMF.
• Coherence and diversity favor NMF (0.391 and 0.407; 0.92 and 0.91). NPMI rewards words that co-occur in many posts, and NMF's objective is precisely to reconstruct those co-occurrences; BERTopic's words are chosen to be distinctive of a cluster, which brings in rarer words (names of software, people, places). A coherence score measures one aspect of quality, not the truth.
• Stability. No model gives the same topics with another seed. LDA's assignments change the most (ARI 0.164 between the two runs), NMF's the least (0.736, thanks to its SVD-based start), and BERTopic's change in a more visible way: the seed of UMAP changes even the number of topics, 24 against 16 (ARI 0.324).
• Cost (not in the table, since timings vary by machine): the embedding step dominates BERTopic's cost, a few seconds per thousand short documents on a GPU and tens of seconds on a CPU; after that, UMAP and HDBSCAN take seconds. LDA and NMF take seconds on a CPU and need no model download.

The overall comparison:

Pitfalls

Stop words in the representation. Seen in step 5: without a stop-word list, BERTopic names topics by the, of, to; and even with the English list, LDA's topics on the posts kept everyday words such as don, think, like, just. A list for the wrong language does not help: in this corpus, BERTopic's topic 6 is named by van, het, een, Dutch function words from a handful of posts written in Dutch, which the English list does not know. Other classic junk topics are signatures, e-mail footers and mailing-list templates: text that repeats verbatim and forms very dense clusters.

Too many outliers. 20% of the posts in topic −1 is common. reduce_outliers reassigns them to the nearest topic, using either the c-TF-IDF of the document or its embedding:

for strategy, extra in [('c-tf-idf', {}), ('embeddings', {'embeddings': E})]:
    new = np.array(topic_model.reduce_outliers(docs, topics, strategy=strategy, **extra))
    print(f"{strategy:10s}  outliers {(new == -1).sum():>3}   homogeneity on all posts "
          f"{homogeneity_score(y, new):.3f}")

c-tf-idf    outliers   0   homogeneity on all posts 0.822
embeddings  outliers   0   homogeneity on all posts 0.882

Every post now has a topic, and homogeneity, measured on all 2,155 posts, falls from 0.927 (on the posts BERTopic kept) to 0.822 or 0.882, depending on the strategy. That is still well above LDA and NMF in the table. The outliers were the ambiguous posts, and forcing them into topics makes the topics less pure: a trade-off, not a free fix. After reassigning, call topic_model.update_topics(docs, topics=new) so that the topic words are recomputed with the new members. A threshold argument keeps as outliers the documents that are not similar enough to any topic.

Topic granularity is a cliff, not a dial. min_cluster_size (called min_topic_size in BERTopic) decides the number of topics, and not smoothly:

for mcs in [5, 10, 15, 25]:
    lab = HDBSCAN(min_cluster_size=mcs, metric='euclidean', cluster_selection_method='eom').fit_predict(Z)
    print(f"min_cluster_size={mcs:<3} -> {lab.max() + 1:>3} topics, {(lab == -1).mean():.0%} outliers")

min_cluster_size=5   ->  96 topics, 28% outliers
min_cluster_size=10  ->  24 topics, 20% outliers
min_cluster_size=15  ->   5 topics, 0% outliers
min_cluster_size=25  ->   3 topics, 0% outliers

From 10 to 15 the model goes from 24 topics to 5, because HDBSCAN's "excess of mass" rule then prefers the large parent clusters to their children. Scan a few values and read the topics. To get fewer, broader topics without changing the clustering, topic_model.reduce_topics(docs, nr_topics=...) merges the most similar topics after fitting (the lab does this); BERTopic(nr_topics=...) does it during fitting.

Reproducibility. Three sources of randomness: UMAP (fix random_state, which also makes it single-threaded and slower), the embeddings (identical on the same hardware and library versions, but not bit-for-bit between a CPU and a GPU) and, for LDA and NMF, the starting point (random_state). The embeddings deserve a test. Add to them a noise far too small to see, keep every seed, and run stages 2 and 3 again:

print(f"original       : {labels.max() + 1} topics, {(labels == -1).mean():.0%} outliers")
for noise_seed in [0, 1, 2]:
    rng = np.random.default_rng(noise_seed)
    E_noisy = E + rng.normal(0, 1e-6, size=E.shape).astype(np.float32)   # changes of ~0.000001
    Z_noisy = UMAP(n_neighbors=15, n_components=5, min_dist=0.0, metric='cosine',
                   random_state=42).fit_transform(E_noisy)              # same seed as stage 2
    lab_noisy = HDBSCAN(min_cluster_size=10, metric='euclidean',
                        cluster_selection_method='eom').fit_predict(Z_noisy)
    print(f"noise, seed {noise_seed}: {lab_noisy.max() + 1} topics, {(lab_noisy == -1).mean():.0%} outliers,"
          f" ARI with the original {adjusted_rand_score(labels, lab_noisy):.3f}")

original       : 24 topics, 20% outliers
noise, seed 0: 24 topics, 19% outliers, ARI with the original 0.908
noise, seed 1: 24 topics, 20% outliers, ARI with the original 0.910
noise, seed 2: 3 topics, 0% outliers, ARI with the original 0.514

Changes of about one millionth, the order of the differences between computing the embeddings on a CPU or on a GPU, change the partition a little twice (ARI 0.908 and 0.910) and, the third time, collapse the 24 topics into 3. min_cluster_size=10 sits at the edge of the cliff measured above, and the smallest push can send HDBSCAN over it. Fixing seeds makes a run repeatable; it does not make the topics true. Report topics that survive several seeds and small perturbations (the table above showed two UMAP seeds giving 24 and 16 topics), prefer settings away from a cliff (the lab uses min_cluster_size=20 for that reason), and store the fitted model (topic_model.save(...)) when others must see exactly the same topics.

The language of the embedding model. all-MiniLM-L6-v2 was trained on English. For Portuguese or mixed corpora, use a multilingual model such as paraphrase-multilingual-MiniLM-L12-v2 (or BERTopic(language='multilingual'), which picks one), and a Portuguese stop-word list in the vectorizer.

Next: the lab repeats the whole lesson on five other forums.