Skip to content

Rewrite as a DVC pipeline over the complete FAERS corpus (2004-2026) - #1

Merged
neksa merged 10 commits into
masterfrom
feat/dvc-pipeline
Sep 14, 2026
Merged

neksa merged 10 commits into
masterfrom
feat/dvc-pipeline

Conversation

@neksa

@neksa neksa commented Aug 6, 2026 •

Copy link
Copy Markdown
Owner

Replaces the 2020 notebook with a DVC pipeline over the full FAERS corpus: 90 quarterly archives, 2004Q1–2026Q2, both the LAERS and FAERS eras, 20.3M deduplicated cases.

The notebook used the openFDA API, which caps pagination at 25k records and has no case-version fields — so its record-level results came from a 10k-report sample and it couldn't deduplicate. FDA's ASCII archives are 3.5 GB for the whole corpus (vs 113 GB of equivalent JSON) and carry caseid/caseversion.

Findings

  • 16.7% of raw records are redundant (4.0M superseded versions, 67,838 retracted)
  • Median 1 drug per report, mean 3.18 — the notebook couldn't compute this
  • Consumers now file more reports than physicians, 8.9M vs 4.5M
  • 0.147% of drug records have dechallenge + rechallenge, not the ~2% previously estimated
  • 2.00M drug–event pairs scored, 685,921 flagged at EB05 ≥ 2
  • 42% of reaction terms show MedDRA vocabulary drift, affecting 24% of mentions
  • Reaction profiles recover drug classes: statins, ACE inhibitors, TNF inhibitors, SSRIs and bisphosphonates are 4.8× more similar within class

Statistics

ROR, PRR+χ², IC/BCPNN, EBGM/MGPS and Fisher+FDR, all screened on lower bounds rather than point estimates. Validated against nine established label warnings — all nine flagged by all five rules. Mantel–Haenszel and Breslow–Day calibrated against simulated data with a known odds ratio.

Fixes the notebook's PRR variance (+ 1/(c+d) where it's a subtraction).

Running it

uv sync && dvc repro

No credentials or cloud storage — FDA hosts the inputs, dvc.lock pins every intermediate. ~25 min, mostly download.

Notes

  • Report renders to docs/, served by GitHub Pages on merge
  • 142 tests, CI runs tests + lint + pipeline contract + a secret scan
  • MIT code, CC BY 4.0 derived tables; the pt column is MedDRA, covered in DATASET_CARD.md
  • publish uploads to a HF dataset repo, opt-in and excluded from dvc repro

Limitations are in the report: no causality or incidence, 28% of drug labels are still product names (flagged), drift is flagged not corrected, independent duplicates estimated not removed.

neksa added 10 commits August 6, 2026 18:27
Replaces the 2020 exploratory notebook with a nine-stage DVC pipeline covering
all 90 quarterly FAERS/LAERS archives, 2004Q1-2026Q2.

The notebook worked through the openFDA API, which caps pagination at 25k
records and exposes no case-version fields. Its record-level results came from
a 10,000-report sample, and it could not deduplicate, compute reactions per
report, or adjust for the confounders it identified. FDA's ASCII archives are
3.5 GB for the whole corpus against 113 GB for equivalent openFDA JSON, and
carry caseid/caseversion.

Findings:
- 16.7% of raw records are redundant: 4.0M superseded versions, 67,838
  retracted. This is the artifact the notebook saw as "repetitive" clustergram
  rows and could not remove.
- Median 1 drug per report (mean 3.18), 2.98 reactions; tail to 171 drugs.
- Consumers now file more reports than physicians, 8.9M against 4.5M.
- 0.147% of drug records carry dechallenge and rechallenge validation, not the
  ~2% estimated from marginals treated as independent.
- 2.29M drug-event pairs scored; 835,291 flagged at EB05 >= 2.
- 98.5-100% of signals survive Mantel-Haenszel adjustment, but Breslow-Day
  rejects homogeneity for 50-76%.

Statistics screen on lower bounds (ROR025, PRR025, IC025, EB05), never point
estimates. Validated against nine established label warnings, all flagged by
all five rules; Mantel-Haenszel and Breslow-Day calibrated against simulated
data with a known common odds ratio.

Fixes the notebook's PRR variance, which used + 1/(c+d) where the term is a
subtraction.

Fixes a determinism bug found by diffing two pipeline runs: equal-support ties
in the drug-name dictionary were resolved by row arrival order, so identical
inputs produced different signal counts.

Removes the notebook, its HTML copy, environment.yml and the Binder/Colab
badges; the narrative moves into the Quarto report.

75 tests.
… report

Reviewing publication-readiness found that 37.3% of distinct "ingredients" in
the signal table were unnormalized product strings -- HUMIRA . PEN,
DURAGESIC-100, PARAGARD 380A, ADVAIR DISKUS 100/50. This splits marginals in
both directions at once: the substance is undercounted because reports filed
under the brand never reach it, and the brand fragment earns an inflated
disproportionality from its small marginal.

Extends the resolution cascade from three steps to five, using only free
sources:

- ndc_brand: FDA NDC directory (public domain), first-token matched. A token
  resolves only when one ingredient set dominates its records, so SODIUM and
  NEXIUM are declined rather than mapped to whichever formulation is most
  numerous. 7,958 of 10,660 tokens usable; 2,702 rejected as ambiguous.
- rxnav_brand: NLM RxNav for brands the NDC directory misses, notably
  biologics whose NDC entries carry no active ingredients. Answers cached,
  negatives included.

Also strips device and product designators before lookup (PEN, DISKUS, HFA,
SOLOSTAR, trailing product codes). The designator rule requires three or more
digits or a digit-letter suffix, because NONOXYNOL-9, OMEGA-3 and COQ-10 are
substances and stripping their suffix would invent new ones.

Results: rows resolved 97.55% -> 98.56%; unresolved rows 2.00M -> 1.18M;
distinct ingredients 460,755 -> 337,179; fallback-only share of signal-table
ingredients 37.3% -> 23.5%. All nine positive controls still flag on all five
rules, and several strengthen as their brands consolidate --
atorvastatin-rhabdomyolysis from n=2,397/EB05=15.5 to n=2,794/EB05=28.6.

Publishes the provenance: scored.parquet gains ingredient_curated, so a
consumer can distinguish a resolved substance from an unmapped product name.
That distinction existed upstream and was previously dropped at publication.

Corrects a false statement in the report, which claimed unresolved rows fall
below the minimum-support threshold and never enter the signal tables. They
do, and always did. The claim was inferred from the min_count rule rather than
measured. The report now says so explicitly rather than quietly editing it out,
and the figure is computed from the data so it cannot go stale.

Fixes a silent RxNav failure: urlencode escapes a literal "+" to %2B, so
tty="IN+MIN" asked for one nonexistent term type and every brand resolved to
None without raising. 250 lookups returned zero ingredients before the counter
showed it.

96 tests.
…tion

Closes the remaining gaps identified while reviewing publication-readiness.

Infrastructure:
- CI runs pytest, ruff, and a DVC graph + params contract check. Nothing ran
  automatically before.
- GitHub Pages was returning 404: it serves master root while the report
  rendered to report/_site. Now renders to docs/ with .nojekyll.

MedDRA, resolved: the pipeline never needed a licence. It reads preferred
terms out of FAERS as opaque strings with no dictionary lookup. PyMedTermino
does not help -- it ships no terminology, only a reader. A dictionary is
needed for version drift and hierarchy, and both have licence-free routes.

New stage `drift`: detects MedDRA vocabulary change empirically. A term that
goes from absent to established between adjacent quarters is behaving like a
revision, not epidemiology. 3,100 of 7,382 tested terms flagged, touching 24%
of reaction mentions -- a limitation this report previously only declared.
Tested against synthetic corpora: a smooth trend is not flagged, an abrupt
step is.

New stage `duplicates`: record linkage for the same incident reported under
different case ids, which version-based dedup cannot see. 558,132 redundant
records in 287,235 clusters, 9.1% of linkable records. Estimated and reported,
never removed -- dropping half a million records would move every published
number, so it is a params switch that defaults to off.

Getting the linkage tractable took three attempts. Demographic blocking alone
produced 128.8M candidate pairs and exhausted 739 GB of scratch. Folding a
shared drug into the blocking key, intersecting with a shared-reaction key,
and materializing the restricted attribute sets rather than re-scanning 82M
rows per reference brings it to 5.3M pairs.

Schema, following OFFSIDES/TWOSIDES so the output is joinable:
- a, b, c, d contingency cells published explicitly
- ingredient_rxcui on 96.7% of drug rows, up from 29.9% -- the brand-only
  lookup covered 560 of 337k ingredients, so ingredient names are now resolved
  directly against RxNorm

Fixes a metric that was misleading twice over. It reported 96.6% because it
counted all ingredient strings including the one-off tail; the replacement
proxy reported 86% because "appears in >= 3 reports" admits eight times as
many strings as the signal table. It now lives in the signals stage where that
population is defined, and reads 23.5%, matching the report's independent
computation.

Also fixes DuckDB reserved-word aliases, bounds scratch space so a runaway
query fails fast instead of filling the disk, and removes a report paragraph
claiming independent duplicates were not attempted.

Adds DATASET_CARD.md.

123 tests.
umap-learn pulled llvmlite 0.36, which refuses to build on Python 3.12 and
broke CI on its first run. None of the three is imported anywhere: they were
declared for a drug-similarity clustering stage that was planned but never
implemented, so the dependency outlived the feature.

The clustering analysis itself remains missing relative to the 2020 notebook,
which did UMAP + DBSCAN over drug ADR profiles. Tracked separately rather than
papered over.
CI could not run the pipeline-contract job: dvc was installed system-wide on
the development machine but never declared, so `uv run dvc` found nothing.
The kind of gap that only surfaces once something other than your own laptop
tries to build the project.
A non-commercial academic project, so the free MedDRA tier applies and the
derived tables can be released openly.

Adds MIT for the code and CC BY 4.0 for the derived tables, with a per-component
breakdown in both the README and the dataset card.

Corrects a false claim in the dataset card, which asserted that this project
holds a MedDRA subscription. It does not yet.

Also draws a distinction the earlier draft blurred: a licence set here covers
this project's own contribution -- the contingency counts, the measures, the
code. It cannot relicense MedDRA, whose preferred terms populate the pt column
and arrive through FDA's public release rather than from the MedDRA
distribution. Consumers are told which single column carries that constraint,
so anyone without a MedDRA position can drop it and use everything else.
Two things, found together: the similarity analysis was the one section of the
2020 notebook this rewrite had not reproduced, and verifying its output
surfaced a splitting bug in the drug column.

Salt forms were being counted as separate drugs. FAERS records ATORVASTATIN and
ATORVASTATIN CALCIUM, and treating them as two substances divided one drug
almost exactly in half -- 275,333 reports against 269,328 -- understating its
marginal and inflating the disproportionality of both halves. Warfarin,
amlodipine and metformin were split the same way.

Resolved through RxNorm's active moiety rather than by stripping suffixes,
because a regex gets the hard cases wrong: SODIUM CHLORIDE and CALCIUM
CARBONATE are drugs in their own right and would become "sodium" and "calcium".
RxNorm returns those unchanged and collapses only genuine salts. Combination
products with several base ingredients are left intact.

Effect: atorvastatin 539,618 reports where it was two entries of ~270k;
warfarin 221,287; amlodipine 496,013. Positive controls still 9/9 on all five
rules and better evidenced -- warfarin-haemorrhage n=2,778 -> 4,511,
atorvastatin-rhabdomyolysis 2,794 -> 5,183.

New stage `similarity`: drug-drug cosine similarity over information-component
reaction profiles, hierarchically clustered. Uses IC rather than raw counts, or
every drug resembles every other through nausea and fatigue. No UMAP: its
output depends on hyperparameters that were never reported, so clustering it
clusters an artefact, and it pins an llvmlite that will not build on current
Python.

Unlike the notebook, which described two clusters and left interpretation "out
of scope", this is validated. Statins, ACE inhibitors, TNF inhibitors, SSRIs
and bisphosphonates average 4.8x more similar within class than to the corpus,
with no class information given to the algorithm. Nearest neighbours are
recognisable pharmacology: adalimumab to infliximab and etanercept, sertraline
to fluoxetine and paroxetine, atorvastatin to other statins and ezetimibe.

The report states the caveat: this cannot separate shared mechanism from shared
indication or co-prescription, the same confound that puts indication bias at
the top of the disproportionality rankings.

Also notes why the fallback-only share rose from 23.5% to 28% while the count
fell from 2,040 to 2,008: merging salts shrank the curated denominator faster
than the numerator. The absolute problem got smaller; the ratio did not.

Drops the MedDRA subscription as a prerequisite. The pipeline performs no
terminology lookup and needs no licence to run.

140 tests.
The publish stage had never been executed. It has now run end to end against a
private Hugging Face dataset repo and the result reads back: 1,997,725 rows,
31 columns, values intact.

Two fixes it needed:
- The dataset card must be uploaded as README.md. Hugging Face renders the
  dataset page from that filename; under any other name the page stays blank.
- Honour a `private` flag, defaulting to true. Making a dataset public should
  be a separate deliberate act, not a side effect of the first upload.

Adds .env to .gitignore, with .env.example documenting the keys. The file was
untracked but not ignored, so any `git add -A` would have committed a live
token. It had not reached history -- confirmed by searching every commit -- but
that was timing rather than design.

HF_TOKEN is read from the environment or .env by a small local parser rather
than a dependency; existing environment variables win, so CI and one-off
overrides behave as expected.

Expands the published payload to cover the stages added since it was written:
similarity, drift and duplicates.
A live HF token sat in the working tree, untracked but not ignored, for several
days. It never reached history -- confirmed by searching every commit for the
literal value -- but only because no commit happened to land in that window.
That is luck, not a control.

Adds a CI job that fails if a .env file appears anywhere in history, if a
credential-shaped string is committed, or if .env stops being ignored. The
placeholder in .env.example is excluded by pattern length.

Also verified the published dataset carries no record-level data: all 22
published tables are aggregate, and results/duplicates/ -- which holds
record-id pairs from the linkage stage -- is not in the publish payload.
The README hardcodes headline numbers and had already gone stale: it quoted 287k
duplicate clusters after the stage re-ran and moved to 292k. The report does not
have this problem because it computes every figure at render time; the README
cannot, so scripts/check_docs.py verifies it against results/ instead. It skips
when results/ is absent, so CI runs it without a pipeline run.

It earned its place immediately, catching an eb05 drift of 685,128 -> 685,921
after the salt merge. That turned out to be the RxNav cache settling rather than
a defect -- two consecutive normalize runs are now byte-identical -- but the
number would otherwise have been published wrong. It then missed a mangled
figure in the adjacent sentence, so its coverage was extended to the pair count.

Clears both deprecation warnings rather than carrying them:
- explode(empty_as_null=...) is pinned, because the default flips in Polars 2.0.
  An empty fragment from a trailing separator must become null so the filter
  drops it; under the new default it would survive as an empty-string
  "ingredient".
- A class-scoped pytest fixture defined as an instance method, removed in
  pytest 10.

Pipeline re-run after the explode change: positive controls still 9/9 on all
five rules, outputs deterministic.

142 tests, no warnings.
@neksa
neksa merged commit 853c874 into master Sep 14, 2026
6 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant