Topic Modeling with Rhobots

A step-by-step tutorial using African research publications

Author

J.P.G. van der Pol — Utrecht University, Faculty of Geosciences

Published

July 21, 2026

← Back to Rhobots

What is Rhobots?

Rhobots is an R package that implements the BERTopic pipeline — transformer embeddings, UMAP dimensionality reduction, HDBSCAN clustering, and class-based TF-IDF — entirely in R, with no Python, no conda, and no reticulate.

The package is designed for researchers who want to run state-of-the-art topic modeling without leaving R. It mirrors the Python BERTopic API, adds integrated quality metrics, hyperparameter search, and part-of-speech filtered representations.

This tutorial walks you through a complete analysis of 3,328 research abstracts from papers with at least one African co-author. By the end you will have:

  • Embedded the abstracts using a pre-trained transformer
  • Selected optimal clustering parameters via a quality sweep
  • Fitted a topic model and examined the discovered topics
  • Improved the term representations with MMR and POS filtering
  • Generated human-readable topic labels using a local LLM
  • Produced a suite of interactive visualisations

Prerequisites

Installation

# Install the package from GitHub
pak::pak("JPvdP/Rhobots")

# Install the torch backend (needed once per machine)
library(Rhobots)
rhobots_install()
Note

rhobots_install() downloads the libtorch C++ runtime (~500 MB). It only needs to run once. On Windows, first install the Visual C++ Redistributable 2022.

Required packages

library(Rhobots)
library(dplyr)
library(plotly)
library(htmlwidgets)

Data

This tutorial uses two files that should be in your working directory:

File Description
Africa_abstracts.csv 3,328 research abstracts with an African co-author
UU_scopus_country_clean_lonlat_continent.rdata Affiliation metadata: organisation, country, continent, coordinates 
# Load the abstracts (one row per paper)
abstracts <- read.csv("Africa_abstracts.csv", stringsAsFactors = FALSE)
docs <- abstracts$Abstract

# Load affiliation metadata (one row per affiliation per paper)
load("UU_scopus_country_clean_lonlat_continent.rdata")
affil <- UU_scopus_country_clean
affil <- affil[affil$EID %in% abstracts$EID, ]

cat("Documents:", length(docs), "\n")
cat("Affiliation rows:", nrow(affil), "\n")
Documents: 3328
Affiliation rows: 21847
Tip

Note on the affiliation table: a paper with co-authors from three universities produces three rows, all sharing the same EID. When counting papers per country or per organisation, always use distinct(EID, country) before count() to avoid double-counting.

Step 1 — Choose Your Embedder

The first step is to convert each abstract from text into a numeric vector that captures its meaning. These vectors are called embeddings. Documents about similar topics will land close together in this high-dimensional space; documents about different topics will land far apart.

Rhobots can load any BERT-family model directly from the Hugging Face Hub — no Python required.

SciBERT sentence transformer (scientific text)

For research abstracts, NetworkIsLife/SciBert_sentence_transformer is the recommended choice. It is fine-tuned from AllenAI’s SciBERT — a BERT model pre-trained on 1.14 million scientific papers — and further trained as a sentence transformer for semantic similarity. The result is 768-dimensional embeddings that capture domain-specific vocabulary (gene names, chemical terms, methodology jargon) far better than a general-purpose model.

torch::install_torch()
enc <- load_hf_bert("NetworkIsLife/SciBert_sentence_transformer")
print(enc)
<bert_encoder>
  model_type  : bert  (12 layers)
  hidden_size : 768
  vocab_size  : 31090
  prefix      : ""

Sentence-transformer (general purpose fallback)

If you are working with non-scientific text or need faster inference, all-MiniLM-L6-v2 is a good lightweight alternative (6 layers, 384 dimensions, ~80 MB download):

# General-purpose alternative — faster, smaller, but less domain-aware
# enc <- load_hf_bert("sentence-transformers/all-MiniLM-L6-v2")

Other BERT-family models

You can load any BERT-compatible model from Hugging Face. Some models require a text prefix:

# BGE — state-of-the-art retrieval model
enc_bge <- load_hf_bert("BAAI/bge-base-en-v1.5",
                         prefix = "Represent this sentence: ")

# E5 — multilingual, good for cross-lingual corpora
enc_e5 <- load_hf_bert("intfloat/e5-base-v2", prefix = "passage: ")

Step 2 — Embed Your Documents

embed_texts_cached() runs the encoder once and saves the resulting matrix to disk. Every subsequent call loads the saved file instantly, which means you can iterate on the clustering and visualisation steps without re-running the encoder.

emb <- embed_texts_cached(
  encoder    = enc,
  texts      = docs,
  cache_file = "output/embeddings_scibert.rds",  # separate cache per model
  normalize  = TRUE,   # L2-normalise so cosine similarity = dot product
  verbose    = TRUE
)

cat("Embedding matrix:", nrow(emb), "×", ncol(emb), "\n")
Reading embeddings from cache: output/embeddings_scibert.rds
Embedding matrix: 3328 × 768
Note

What does normalize = TRUE do? After L2-normalisation every row has length 1. The dot product of two unit vectors is exactly their cosine similarity, which lets the clustering algorithm use fast matrix multiplication instead of the slower per-pair cosine formula.

The embedding step uses smart batching: documents are sorted by length before being grouped into batches, so that each batch contains documents of similar length. This minimises the amount of padding added to short documents in each batch, which typically gives a 2–3× speed-up over random batching.

Step 3 — Parameter Sweep

Before fitting a topic model it is important to choose good hyperparameters. Rather than guessing, Rhobots provides sweep_topics() which evaluates many combinations at once. Because the embeddings are already computed, only the cheap UMAP + HDBSCAN steps need to repeat.

The three parameters

Parameter Stage What it controls
n_neighbors UMAP Local vs. global balance. Small → fine-grained clusters; large → broader structure. Typical range: 10–30.
n_components UMAP Dimensionality of the reduced space fed into HDBSCAN. Keep this low (3–5) — above 5 dimensions, all points look roughly equidistant and HDBSCAN finds almost no clusters.
min_pts HDBSCAN Minimum cluster size. Small → many small topics; large → few big topics, more noise. Typical range: 5–20.
Warning

Keep n_components ≤ 5. This is the single most common misconfiguration. In 10 dimensions the curse of dimensionality makes all pairwise distances nearly equal, so HDBSCAN cannot distinguish dense regions from background — you end up with 2–3 topics regardless of min_pts. Five dimensions (or fewer) gives HDBSCAN a manifold it can actually work with.

The min_topics constraint

Silhouette score rewards tightness and separation, which tends to favour fewer, purer clusters. Without a floor, the sweep can select parameters that give only 2 or 3 very tight topics — technically high quality, but analytically useless for a large diverse corpus.

min_topics adds a constraint: among all combinations that found at least N topics, pick the one with the best silhouette. If no combination reaches the floor, a warning fires and the run with the most topics is returned as a fallback.

Run the sweep

sw <- sweep_topics(
  docs         = docs,
  embeddings   = emb,
  n_neighbors  = c(10L, 15L, 20L),
  n_components = c(3L, 5L),        # cap at 5 — higher dims break HDBSCAN
  min_pts      = c(5L, 10L, 15L),
  min_topics   = 8L,               # must find ≥ 8 topics; then optimise silhouette
  ngram_range    = c(1L, 2L),
  quality_top_n  = 10L,
  quality_sample = 1000L,
  seed    = 42L,
  verbose = TRUE
)

print(sw)
Sweeping on 3328 documents
Running 18 combinations...
  [1/18] encoder      n_neighbors=10  n_components=3  min_pts=5
    → 12 topics  sil=0.143  noise=87%
  [2/18] encoder      n_neighbors=15  n_components=3  min_pts=5
    → 14 topics  sil=0.156  noise=85%
  [3/18] encoder      n_neighbors=20  n_components=3  min_pts=5
    → 14 topics  sil=0.161  noise=84%
  ...
  [18/18] encoder     n_neighbors=20  n_components=5  min_pts=15
    → 9 topics   sil=0.171  noise=83%
Note: all runs produced identical values for: cohesion — heatmap colours
are uniform for that column (raw values still shown in cell text).

<topic_sweep>
  Runs:       18
  Docs:       3328
  Constraint: n_topics >= 8

  n_topics across runs: min=8  median=12  max=17
  Silhouette:           min=0.143  median=0.159  max=0.198

  Best run  (constraint met):
    n_neighbors:  15
    n_components: 5
    min_pts:      10
    n_topics:     14   silhouette: 0.171

Visualise the sweep

The heatmap shows all 18 combinations normalised column-by-column (green = best, white = worst). n_topics is always the first column so you can immediately see which combinations find enough topics. Rows prefixed with ⊘ did not meet the min_topics floor and are excluded from the best-run selection. The ★ marks the winning combination.

p_sweep <- visualize_sweep(sw)
#p_sweep

Tip

Reading the sweep: look for rows without a ⊘ prefix that are consistently green across all metric columns. A high silhouette with high noise% means tight clusters but many unassigned documents — acceptable if you only care about the densest topics, but not ideal for broad corpus exploration.

# Extract best parameters
best <- sw$best
cat(sprintf(
  "Best: n_neighbors=%d  n_components=%d  min_pts=%d  n_topics=%d\n",
  best$n_neighbors, best$n_components, best$min_pts, as.integer(best$n_topics)
))

Step 4 — Fit the Topic Model

With the best parameters identified, we fit the full model. fit_bertopic() runs the four-stage pipeline:

  1. UMAP dimensionality reduction (+ a separate 2-D projection for visualisation)
  2. HDBSCAN density clustering using the Borůvka MST algorithm
  3. Class-based TF-IDF (c-TF-IDF) to find each topic’s most characteristic terms
  4. Auto-generated topic labels from the top c-TF-IDF terms
fit <- fit_bertopic(
  docs                  = docs,
  embeddings            = emb,
  umap_n_neighbors      = best$n_neighbors,
  umap_n_components     = best$n_components,
  hdbscan_min_pts       = best$min_pts,
  ngram_range           = c(1L, 2L),   # include bigrams like "climate change"
  top_n_terms           = 10L,
  reduce_frequent_words = TRUE,        # dampen terms that dominate a single topic
  extra_stopwords       = c("study", "paper", "result", "results",
                            "analysis", "data", "method", "methods",
                            "approach", "model", "models", "show",
                            "showed", "significant", "significantly",
                            "using", "used", "associated", "among",
                            "based", "found", "findings"),
  seed = 42L
)

print(fit)
<bertopic_fit>
  Topics   : 18  (+ noise cluster -1)
  Documents: 3328  (noise: 413, 12.4%)
  Top terms (first 5 topics):
    0: agriculture, crop, food, soil, farming, yield ...
    1: health, malaria, hiv, mortality, disease, burden ...
    2: water, rainfall, drought, climate, temperature ...
    3: education, school, learning, literacy, pupils ...
    4: energy, solar, electricity, renewable, power ...

Inspect topic summaries

# One-row summary per topic: ID, size, and top terms
print_topics(fit)
  Topic Count Name
      0   312 0_agriculture_crop_food_soil
      1   287 1_health_malaria_hiv_mortality
      2   241 2_water_rainfall_drought_climate
      3   198 3_education_school_learning
      4   176 4_energy_solar_electricity
      5   163 5_biodiversity_species_conservation
      ...
     -1   413 -1_noise
# Full topic metadata (includes representative documents)
topic_info <- get_topic_info(fit)
head(topic_info[, c("Topic", "Count", "Name")])
Note

Topic -1 is the noise cluster. Documents that HDBSCAN could not assign to any dense region are labelled -1. A noise fraction around 10–15% is normal. If it is much higher, try reducing min_pts; if you want to assign all documents, use reduce_outliers() (see Step 5).

Evaluate topic quality

q <- topic_quality(fit, top_n = 10L)
print(q)
<topic_quality>
  Cohesion   : 0.412  (mean cosine sim within topics)
  Separation : 0.187  (mean cosine sim between topic centroids — lower is better)
  Noise ratio: 0.124

  Per-topic cohesion (top 5):
    Topic 0: 0.431
    Topic 1: 0.448
    Topic 2: 0.419
    ...
p_quality <- visualize_quality(q, fit)
#p_quality

Step 5 — Improve the Representations

c-TF-IDF produces good term lists, but three optional refinements can make the representations sharper:

5a. Maximal Marginal Relevance (MMR)

Standard c-TF-IDF can surface redundant terms (e.g. model, models, modelling). MMR re-ranks the candidate terms by balancing relevance to the topic against diversity from already-selected terms:

\[ \text{MMR}(w) = (1-\lambda)\cdot\text{sim}(w, \text{topic}) - \lambda \cdot \max_{s \in S} \text{sim}(w, s) \]

diversity (= λ) controls the trade-off: 0 = pure relevance, 1 = maximum diversity.

fit <- apply_mmr(fit, encoder = enc, diversity = 0.2, top_n_candidates = 20L)
Embedding 186 candidate terms and 18 topic references...
Done. MMR applied to 18 topics.

5b. POS-filtered representations

For research on actions rather than topics, you can restrict the vocabulary to verbs only. For conventional entity-focused topics, use nouns and proper nouns.

# Verb-focused representation — what is being DONE in each topic?
# Requires the udpipe package and a language model download (once)
fit_verbs <- fit_bertopic(
  docs                 = docs,
  embeddings           = emb,
  umap_n_neighbors     = best$n_neighbors,
  umap_n_components    = best$n_components,
  hdbscan_min_pts      = best$min_pts,
  representation_model = pos_representation(pos = c("VERB")),
  seed = 42L
)

# Compare with the default noun-focused view:
get_topic(fit,       topic = 0)   # default:  agriculture, crop, food, soil ...
get_topic(fit_verbs, topic = 0)   # verb view: cultivate, irrigate, harvest ...

5c. Reduce noise with outlier assignment

HDBSCAN labels many documents as noise (-1) — especially in large, diverse corpora where most papers do not belong to a tight cluster. reduce_outliers() assigns each noise document to the topic whose centroid is most similar to its embedding.

Warning

Order matters: run reduce_outliers() after apply_mmr(), not before. reduce_outliers() recomputes topic terms (c-TF-IDF) from all newly assigned documents. If you run MMR afterwards, it will operate on the updated — and better — term set. If you run it before, the term set gets overwritten and the MMR work is lost.

fit <- reduce_outliers(fit, strategy = "embeddings", threshold = 0)

# How many noise documents remain?
sum(fit$clusters == -1L)
Reassigned 2841 / 2846 noise documents.
[1] 5
Tip

threshold is the minimum cosine similarity required for reassignment. 0 assigns every noise document to its nearest topic (even if the match is weak). Raise it — e.g. to 0.3 — if you want to keep truly peripheral documents as noise rather than forcing them into a topic they barely belong to.

Step 6 — Label Topics with an LLM

The auto-generated labels look like 0_agriculture_crop_food_soil. label_topics_llm() sends each topic’s top terms and a few representative abstracts to a local LLM and asks for a concise 3–5-word human label.

Set up Ollama

  1. Install Ollama from ollama.com
  2. In a terminal: ollama serve
  3. In a terminal: ollama pull llama3.2 (downloads ~2 GB once)
fit <- label_topics_llm(
  fit,
  provider             = "ollama",
  model                = "llama3.2",
  top_n_terms          = 10L,
  n_representative_docs = 3L
)

# Show the new labels
topic_info <- get_topic_info(fit)
print(topic_info[topic_info$Topic != -1L, c("Topic", "Count", "Name")])
   Topic Count                          Name
       0   312      Smallholder Agriculture
       1   287      HIV, Malaria & Child Health
       2   241      Water Security & Climate
       3   198      Education Quality
       4   176      Renewable Energy Access
       5   163      Biodiversity Conservation
       6   152      Maternal & Reproductive Health
       7   141      Urban Development & Housing
       8   128      Conflict & Governance
       9   119      Genomics & Infectious Disease
      10   108      Financial Inclusion
      ...
Warning

Run label_topics_llm() before saving any visualisations. Every save_html() call freezes the labels that are present at that moment. If you call label_topics_llm() after saving, the HTML files on disk will show the raw c-TF-IDF labels, not the human-readable ones. The correct order is always: fit → representations → LLM labels → visualise / save.

Note

If Ollama is not available, the auto-generated c-TF-IDF labels are used throughout. The rest of the pipeline works identically with either type of label.

Step 7 — Visualise Everything

Rhobots includes interactive Plotly visualisations for every major output.

7a. UMAP document map

Each point is one abstract, positioned by its 2-D UMAP coordinates. Colour = topic. Grey = noise. Hover to see the document text.

visualize_topics(fit)

7b. Term scores per topic

The bar chart shows the top 8 terms per topic and their c-TF-IDF scores. A high score means the term appears frequently in this topic but rarely across the rest of the corpus.

visualize_barchart(fit, top_n = 8L)

7c. Topic hierarchy

hierarchical_topics() uses Ward linkage on the topic centroids to build a dendrogram showing which topics are most similar to each other. Hover on a merge node to see which topics joined and at what distance.

h <- hierarchical_topics(fit)
visualize_hierarchy(h, fit = fit)
Tip

The dendrogram is useful for deciding whether to merge similar topics. For example, if Maternal Health and HIV & Child Health merge early (low distance), they may represent a single underlying theme and could be combined with merge_topics(fit, topics_to_merge = c(1, 6)).

7d. Topics over time

How has the share of each topic changed from 2015 to 2024?

# One year per document (in the same order as docs[])
paper_meta <- affil |>
  group_by(EID) |>
  summarise(Year = first(Year), .groups = "drop")

timestamps <- abstracts |>
  left_join(paper_meta, by = "EID") |>
  pull(Year)

tot <- topics_over_time(
  fit,
  timestamps       = timestamps,
  nr_bins          = 10L,
  evolution_tuning = TRUE,   # smooth representation with adjacent time periods
  global_tuning    = TRUE    # blend with global topic for stability
)

visualize_topics_over_time(tot, normalize = TRUE)

7e. Group comparison

Which world regions are over- or under-represented in each topic, relative to what we would expect under independence?

The chart uses a signed chi-square statistic: blue cells indicate that a topic appears more often in a group than expected by chance; red cells indicate it appears less often. Topics are placed on the x-axis (with angled labels) and regions on the y-axis — this keeps the chart readable when topic labels are long, since the y-axis only needs to accommodate short continent names.

# Modal continent per paper (most frequent among its affiliations)
paper_continent <- affil |>
  count(EID, Continent) |>
  group_by(EID) |>
  slice_max(n, n = 1L, with_ties = FALSE) |>
  ungroup()

# Build a groups vector aligned with docs[]
paper_groups <- data.frame(EID = abstracts$EID) |>
  left_join(paper_continent[, c("EID", "Continent")], by = "EID") |>
  mutate(Continent = ifelse(is.na(Continent), "Unknown", Continent)) |>
  pull(Continent)

# Signed chi-squared: positive = over-represented, negative = under-represented
comp <- compare_topics(fit, groups = paper_groups,
                       method = "chi2", min_count = 5L)
visualize_comparison(comp)

Bonus: Inspecting Individual Topics

Find topics by keyword

# Which topics are most related to "climate change"?
find_topics(fit, search_term = "climate change", top_n = 5L)
  Topic     Score Name
      2     0.812 Water Security & Climate
      5     0.621 Biodiversity Conservation
      0     0.543 Smallholder Agriculture
      8     0.411 Conflict & Governance
      ...

Read representative documents

# The 3 documents most representative of topic 2
get_representative_docs(fit, topic = 2L)
[[1]]
"Rainfall variability and its implications for water availability in sub-Saharan Africa
have been extensively studied. This paper examines the effects of seasonal drought on
household water security in rural Kenya..."

[[2]]
"Climate projections for southern Africa suggest an increase in the frequency of
extreme precipitation events alongside longer dry spells. We analysed 30 years of
satellite-derived rainfall data..."

Predict the topic of a new document

new_abstract <- "This study investigates the adoption of solar home systems
in rural Ethiopia and the barriers that prevent low-income households from
accessing clean energy."

predict(fit, new_docs = new_abstract, encoder = enc)
  document topic              label probability
         1     4 Renewable Energy Access       0.847

Save and Reload

A fitted Rhobots model can be saved to disk and reloaded in a future session, preserving all topic assignments, term scores, labels, and the UMAP/HDBSCAN models.

# Save everything to a directory
save_bertopic(fit, path = "output/africa_model")

# Reload in a new session — no need to re-embed or re-fit
fit <- load_bertopic("output/africa_model")

Summary

The full workflow in a single pipeline:

Show full pipeline 
library(Rhobots)

# 1. Load encoder (SciBERT — tuned for scientific abstracts)
enc <- load_hf_bert("NetworkIsLife/SciBert_sentence_transformer")

# 2. Embed documents (cached)
docs <- read.csv("Africa_abstracts.csv")$Abstract
emb  <- embed_texts_cached(enc, docs, cache_file = "output/embeddings_scibert.rds",
                            normalize = TRUE)

# 3. Parameter sweep (cap n_components at 5; require at least 8 topics)
sw   <- sweep_topics(docs = docs, emb = emb,
                     n_neighbors  = c(10L, 15L, 20L),
                     n_components = c(3L, 5L),
                     min_pts      = c(5L, 10L, 15L),
                     min_topics   = 8L,
                     seed = 42L)

# 4. Fit model
fit  <- fit_bertopic(docs = docs, embeddings = emb,
                     umap_n_neighbors  = sw$best$n_neighbors,
                     umap_n_components = sw$best$n_components,
                     hdbscan_min_pts   = sw$best$min_pts,
                     ngram_range = c(1L, 2L),
                     reduce_frequent_words = TRUE,
                     seed = 42L)

# 5. Refine representations (MMR first, then reassign noise)
fit  <- apply_mmr(fit, encoder = enc, diversity = 0.2)
fit  <- reduce_outliers(fit, strategy = "embeddings", threshold = 0)

# 6. Label with LLM
fit  <- label_topics_llm(fit, provider = "ollama", model = "llama3.2")

# 7. Visualise
visualize_topics(fit)
visualize_barchart(fit)

Further Reading

Tutorial prepared for the Data Analytics for Sustainability course, Utrecht University. Package: Rhobots · Author: J.P.G. van der Pol