kemal/rarelens
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

docs: a technical introduction to the platform, as a LaTeX article

Browse Source 
docs/blog/rarelens.tex builds a 10-page write-up covering the Nextflow DSL2
pipeline, the three execution backends behind one API call, the Argo and
ArgoCD track, Terraform, the external systems integrated, the ranking, the
benchmark and the model.

It is deliberately as much an account of what was wrong as of what works: the
components that scored evidence nobody had looked up, the missense AUROC that
fell from 0.872 to 0.500 once allele frequency was removed, the ontology walk
that silently dropped 399 terms, and the propagation change that measured
slightly worse than what it replaced and was kept anyway with the numbers
published.

Figures come from one headless-browser script. Screenshots are of the running
application; the six diagrams are hand-written HTML styled from the same
palette, rendered by the same script. Keeping both in one place is what stops
the article drifting from the system, and every number in a diagram has to be
changed deliberately.
This commit is contained in:
Kemal Yaylali
2026-09-12 13:10:47 +01:00
parent c25fb53666
commit 1e46fff2ff
25 changed files with 729 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files
docs/blog/rarelens.tex
+382
View File
@@ -0,0 +1,382 @@
% rarelens - a technical introduction
% Build: make -C docs/blog (or: pdflatex rarelens.tex, twice)
\documentclass[11pt,a4paper]{article}
\usepackage[T1]{fontenc}
\usepackage[utf8]{inputenc}
\usepackage{lmodern}
\usepackage{microtype}
\usepackage[margin=2.4cm,top=2.6cm,bottom=2.6cm]{geometry}
\usepackage{graphicx}
\usepackage{booktabs}
\usepackage{xcolor}
\usepackage{caption}
\usepackage{enumitem}
\usepackage{fancyhdr}
\usepackage{titlesec}
\usepackage{listings}
\usepackage{float}
\usepackage[hidelinks]{hyperref}
\definecolor{ink}{HTML}{16243A}
\definecolor{soft}{HTML}{52627A}
\definecolor{plum}{HTML}{7A1F5C}
\definecolor{line}{HTML}{CFD6DF}
\definecolor{paper}{HTML}{F3F5F7}
\hypersetup{colorlinks=true,linkcolor=plum,urlcolor=plum,citecolor=plum}
\color{ink}
\titleformat{\section}{\normalfont\Large\bfseries\color{ink}}{\thesection}{0.7em}{}
\titleformat{\subsection}{\normalfont\large\bfseries\color{ink}}{\thesubsection}{0.6em}{}
\titlespacing*{\section}{0pt}{1.6em}{0.6em}
\captionsetup{font=small,labelfont={bf,color=soft},textfont={color=soft},width=0.94\textwidth}
\setlist[itemize]{leftmargin=1.2em,itemsep=0.25em,topsep=0.35em}
\setlength{\parskip}{0.55em}
\setlength{\parindent}{0pt}
\pagestyle{fancy}
\fancyhf{}
\fancyhead[L]{\small\color{soft}rarelens --- a technical introduction}
\fancyfoot[C]{\small\color{soft}\thepage}
\renewcommand{\headrulewidth}{0.4pt}
\lstdefinestyle{code}{
basicstyle=\ttfamily\footnotesize\color{ink},
backgroundcolor=\color{paper},
frame=single, rulecolor=\color{line}, framesep=5pt,
breaklines=true, showstringspaces=false,
keywordstyle=\color{plum}, commentstyle=\color{soft}\itshape,
xleftmargin=0pt, aboveskip=0.9em, belowskip=0.9em
}
\lstset{style=code}
% LaTeX is conservative about letting text share a page with a large float, which on a document
% this figure-heavy strands them on near-empty pages of their own. Loosen the thresholds.
\renewcommand{\topfraction}{0.92}
\renewcommand{\bottomfraction}{0.80}
\renewcommand{\textfraction}{0.06}
\renewcommand{\floatpagefraction}{0.82}
\setcounter{topnumber}{2}
\setcounter{totalnumber}{3}
% htbp, not H: a full-width figure that will not fit pinned in place leaves half a page blank.
\newcommand{\fig}[3]{%
\begin{figure}[htbp]\centering
\includegraphics[width=#2\textwidth]{#1}
\caption{#3}\end{figure}}
\begin{document}
\begin{center}
{\LARGE\bfseries rarelens}\\[0.35em]
{\large\color{soft} An end-to-end rare-disease variant triage platform}\\[1.1em]
{\color{soft}\small A technical introduction --- pipeline, infrastructure, model and the evidence behind the ranking}
\end{center}
\vspace{0.6em}
\hrule
\vspace{1.2em}
\section*{In one paragraph}
\textbf{rarelens takes a patient's genome and their clinical phenotype and narrows thousands of
variants to a handful a scientist can actually review, showing the evidence behind every rank.}
It is a working monorepo: a Nextflow DSL2 pipeline that annotates variants with Ensembl VEP, a
FastAPI service that ranks them against Human Phenotype Ontology annotations, a LightGBM model
served from an MLflow registry, a SvelteKit interface built around triage decisions rather than
table filtering, and Terraform for two deployment tracks on Google Cloud --- a serverless one that
idles at roughly \pounds1 a month, and a Kubernetes one with Argo Workflows, Argo Events and
ArgoCD behind a feature flag. Everything runs on public, openly licensed data. No patient data is
used, accepted, or possible to load.
This document is a technical tour. It is also, deliberately, an account of the things that turned
out to be wrong --- including a model metric that collapsed from 0.872 to 0.500 the moment I took
away the feature that was doing all the work, and an ontology bug that silently deleted 399 terms.
\section{The problem}
A rare-disease proband's exome contains something like twenty thousand coding variants. Perhaps a
few hundred are rare and protein-altering. Exactly one, usually, explains why the patient is ill.
Finding it is not a filtering problem so much as an evidence problem: the same variant is
uninteresting in one patient and diagnostic in another, and what changes between them is the
clinical phenotype.
That observation drives the whole design. The interface is not a variant table with filters --- it
is a ranked shortlist where each candidate carries the four pieces of evidence that put it there,
and where a reviewer shortlists or dismisses with a reason that ends up in a case report.
\fig{figures/arch.png}{1.0}{The system. One repository holds the pipeline, the API, the interface,
model serving and the infrastructure for two deployment tracks. The scientist never sees any of it
except the middle column.}
\section{The Nextflow pipeline}
Annotation is a three-process DSL2 workflow, each process a pinned container: \texttt{bcftools}
normalises, Ensembl VEP 113 annotates, and a loader writes results into PostgreSQL and marks the
job succeeded.
\fig{figures/pipeline.png}{1.0}{The pipeline. The same workflow file runs on a laptop, on Google
Batch and under Argo Workflows; the executor is a profile, not a rewrite.}
Two details are worth pulling out, because both cost real debugging time.
\subsection{Variant identity survives annotation}
VEP's \texttt{Location} and \texttt{Allele} output columns trim indel alleles and shift positions
by one, which is correct for VEP's own purposes and wrong as a primary key. A deletion written by
the caller as \texttt{CT>C} comes back as \texttt{->} at a different coordinate, and the loader
then stores a variant that does not exist in the input file.
The fix is to carry identity through a field VEP does not touch. \texttt{NORMALISE} sets the VCF ID
to \texttt{CHROM\_POS\_REF\_ALT}, and the loader parses identity back out of it:
\begin{lstlisting}[language=bash]
bcftools norm -m -any -f $genome | \
bcftools annotate --set-id '%CHROM\_%POS\_%REF\_%FIRST_ALT'
\end{lstlisting}
Verified on a real run: \texttt{22:42126611 CT>C} round-trips with its alleles intact.
\subsection{Stub blocks make continuous integration possible at all}
A VEP cache is 25\,GB. No pull request is going to download one. Every process therefore carries a
\texttt{stub:} block, and CI runs the workflow with \texttt{-stub-run} against a three-record
fixture, checking channel wiring and process contracts in seconds without a container or a cache in
sight. It catches the failure that actually happens in practice --- a renamed output, a channel
that emits the wrong cardinality --- while leaving the scientific correctness to the tests that can
afford to be slow.
\section{Event-driven execution, three ways}
\texttt{POST /cases/\{id\}/annotate} writes a job row and hands off. What happens next is a
configuration decision, not a code path the caller knows about.
\fig{figures/events.png}{1.0}{One entry point, three backends, chosen from settings alone.}
\begin{itemize}
\item \textbf{Cloud Run job.} The Nextflow driver runs as a Cloud Run job started through the
Jobs API with argument overrides. It scales to zero between runs, and its service account
holds \texttt{run.jobsExecutorWithOverrides} on exactly one job --- not project-wide.
\item \textbf{Pub/Sub.} The job is published as an event. An Argo Events sensor subscribes and
triggers an Argo Workflow on GKE. This is the decoupled path: retries, ordering and
back-pressure become the queue's problem rather than the API's, and the API can be
restarted mid-pipeline without losing work.
\item \textbf{Local subprocess.} Nextflow runs directly and its stdout is streamed into the job
log, so the interface shows live progress. This is what a developer gets with nothing
configured, and it is the same code path the other two wrap.
\end{itemize}
All three converge on the same \texttt{jobs} row, so the interface polls one endpoint regardless.
The database URL is passed by environment or as a Nextflow secret and never appears on a command
line, keeping it out of \texttt{.command.sh} and the workflow logs.
\section{Kubernetes, Argo and GitOps}
The Kubernetes track is Kustomize bases with a local overlay (kind, an in-cluster Postgres) and a
GCP overlay (Cloud SQL, Workload Identity). Argo Workflows runs the annotation
\texttt{WorkflowTemplate}; Argo Events holds the Pub/Sub EventSource and the sensor that triggers
it; ArgoCD reconciles the cluster from the repository, and a green CI run on \texttt{main} bumps
image tags in the GCP overlay so that deployment is a commit rather than a command.
One bug from this area is worth recording because it is invisible until it bites: Kustomize
generates hashed ConfigMap names so that a configuration change forces a rollout, but the hashed
name is only substituted into workloads Kustomize believes are in scope. The overlays were missing
\texttt{namespace: rarelens}, so the substitution silently did not happen and pods mounted a
ConfigMap name that no longer existed. The symptom was a pod stuck in \texttt{CreateContainerConfigError}
with nothing wrong in the manifests as written.
\section{Infrastructure as code}
Terraform provisions both tracks from one root module, with the expensive half behind flags:
\begin{lstlisting}[language=bash]
terraform apply -var project=<id> # serverless: Cloud Run + Batch
terraform apply -var project=<id> -var deploy_kubernetes=true \
-var deploy_cloud_sql=true # adds GKE, Argo, Cloud SQL
\end{lstlisting}
The default track provisions Cloud Run services for the API and interface, a Cloud Run job for the
Nextflow driver, Google Batch for pipeline tasks on Spot VMs, a GCS bucket, Secret Manager entries
and Artifact Registry. GKE Autopilot and Cloud SQL are opt-in, because a Kubernetes control plane
and a managed database are most of what a demonstration estate costs.
CI authenticates to Google Cloud through Workload Identity Federation, so there is no service
account key anywhere in the repository or in GitHub secrets. Terraform is validated and
format-checked on every pull request.
\subsection{Cost as a design constraint}
A portfolio platform is idle more than 99\% of the time, which makes idle cost the only cost that
matters. The serverless track is built around that: Cloud Run at \texttt{min-instances=0}, a driver
that exists only while a pipeline runs, and Batch on Spot. Idle cost lands near \pounds1 a month,
almost all of it the database. The Kubernetes track exists to demonstrate the GitOps path and is
meant to be destroyed afterwards.
\section{Integrating external systems}
Almost nothing here is self-contained, and integrating public biological infrastructure is most of
the work:
\begin{itemize}
\item \textbf{Ensembl VEP} in cache mode, or against Ensembl's \emph{public database server}
with \texttt{--database} when 25\,GB is not available --- slower per variant, no plugin
scores, but no download.
\item \textbf{gnomAD v4.1} allele frequencies streamed by genomic region straight out of the
public Google Cloud bucket. The files are tabix-indexed, so a range request returns a few
hundred records without fetching the file. I had previously written off frequencies as
impossible without the full cache; testing that assumption disproved it.
\item \textbf{Human Phenotype Ontology} gene-to-phenotype annotations plus the ontology itself,
propagated and weighted at load time.
\item \textbf{ClinVar} for model training labels, filtered to two-star review status and above.
\item \textbf{Phenopacket Store} (Monarch Initiative), which curates published case reports into
GA4GH phenopackets --- the source of both the demonstration case and the benchmark.
\item \textbf{MLflow} as a model registry, resolved by alias, with a registry-free path for
deployments that should not run a tracking server.
\end{itemize}
\begin{lstlisting}[language=bash]
# real gnomAD frequencies with no bulk download, verified end to end
bcftools view -r chr3:30672000-30673000 \
https://storage.googleapis.com/gcp-public-data--gnomad/release/4.1/...chr3.vcf.bgz
# -> 392 records, 5 KB, then fed to VEP with --custom alongside --database
\end{lstlisting}
\section{The ranking, and what it refuses to claim}
The rank is a weighted mean of four components a reviewer can audit. ClinVar is deliberately
\emph{not} one of them --- it sits beside the result as independent confirmation, so nothing ranks
highly merely because ClinVar already called it pathogenic.
\fig{figures/evidence.png}{1.0}{Rarity and consequence filter; phenotype only ranks, because a real
diagnosis can sit in a gene nobody has annotated yet. A component with no evidence behind it
abstains, and the remaining weights renormalise.}
The abstention rule replaced a genuine bug. Run without a VEP cache, there are no allele
frequencies --- and the code read a missing frequency as \emph{absent from gnomAD, therefore
maximally rare}, handing every variant a free quarter of its score. The model, separately, was
returning 0.887 for every variant on features it had never been given. Two of four components were
fiction, and the total looked fully informed.
The fix is structural rather than cosmetic: the loader records what the annotation run actually
produced (\texttt{has\_frequencies}, \texttt{has\_effect\_scores}), components without evidence
return null instead of a number, and the weights renormalise over whatever is left. The interface
prints ``not looked up'' where it would otherwise have drawn a bar.
\section{Measuring the ranking}
One demonstration case ranking correctly is an anecdote. The benchmark asks the only question that
matters for a phenotype-driven tool, across every usable case in Phenopacket Store.
\fig{figures/benchmark.png}{1.0}{Retrieval against 5,269 candidate genes, and an honest report of
what the measurement cannot tell you.}
Two things in that figure matter more than the headline. The first is the contamination: HPO's gene
annotations are curated from these same case reports, so the median causal gene already carries
every one of its patient's terms. The number is an upper bound and is labelled as one. The second
is that when I measured my own improvements, one of them did not work --- information-content
weighting helped, ontology propagation cost about as much as weighting gained. That result is in
the documentation with the table, rather than quietly dropped.
\section{The model, and the number that changed it}
The pathogenicity model is LightGBM trained on ClinVar labels, held out \emph{by gene} rather than
by variant. That distinction is not pedantry: a random split puts variants of the same gene on both
sides, and the model then scores the gene instead of the variant, which is exactly the inflation
Grimm \emph{et al.} documented for this class of tool in 2015.
\fig{figures/model.png}{1.0}{Removing one feature moved missense AUROC from 0.872 to 0.500.}
The model originally took allele frequency as a feature and looked respectable. Probing it showed
frequency dominating everything --- the same missense variant scored 0.887 at frequency zero and
0.0003 at one per cent. That is two separate problems. It double-counted, because the ranking
already scores frequency explicitly, putting roughly 45\% of every rank on one measurement. And it
was circular, because ACMG's BA1/BS1 criteria assign ClinVar's benign labels \emph{using} allele
frequency, so the model was rediscovering the rule that had generated its own labels.
Retraining without it returned missense AUROC 0.500 --- exactly random. With frequency gone and no
CADD or AlphaMissense scores in the training table, nothing is left but the consequence class, so
every missense variant scores identically. The conclusion is unambiguous and slightly
uncomfortable: the model never had variant-effect knowledge. It now abstains from the ranking
unless it has a predictor the other components do not already provide, and 0.500 is the measurement
that justifies the abstention.
\section{The application}
The interface is built around the narrowing and the decision, not around the table.
\fig{figures/02-case.png}{1.0}{A published case. The funnel across the top is the story --- variants
called, rarity, coding candidates, phenotype-matched --- and each candidate carries its evidence as
chips. The note under the header states plainly which components did not score and why.}
\fig{figures/04-panel-only.png}{0.62}{Selecting a candidate opens the arithmetic. Every component is
shown with its weight and its contribution, and the two that had nothing to go on say so rather
than displaying a zero.}
\fig{figures/06-report.png}{0.92}{The case report: shortlisted and dismissed variants with the
reviewer's reasons, the funnel counts, and provenance --- VEP version, model version, run time ---
so a result can be reproduced or challenged later.}
\section{Engineering practice}
\begin{itemize}
\item \textbf{Tests where the risk is.} 99 Python tests plus 49 in the front end, concentrated on
the ranking arithmetic, the ontology handling, the loader's idempotency and the security
boundary. Database tests run against a real PostgreSQL --- an embedded server locally, a
service container in CI --- because the schema is part of the behaviour.
\item \textbf{Migrations are tested, both ways.} Alembic upgrade and downgrade are exercised
against a live database, which is how a downgrade that left a stale enum type behind was
caught.
\item \textbf{CI runs the whole estate:} lint and types for the API, unit tests for the model
code, a real Postgres for the loader, a Nextflow stub run, and \texttt{terraform validate}.
\item \textbf{Input validation as a security boundary.} A case's \texttt{vcf\_uri} must be a
\texttt{gs://} object or an absolute path beneath a configured data root, with a VCF
suffix. That is what stops a crafted path becoming a Nextflow option or reading an
arbitrary file, and it has its own test file.
\end{itemize}
\section{Deployment}
The demonstration runs on Railway --- three services, with the API reachable only on the private
network so the interface's \texttt{/api} proxy is the single public entry point. One shared
credential sits in front of it.
Railway cannot run the pipeline: Nextflow shells out to \texttt{docker run} for VEP and bcftools,
and a container platform gives you a container, not a Docker daemon. Rather than leave a button
that always fails, the cases are annotated where the pipeline works and copied up, and the
pipeline actions are hidden. Making it idle correctly needed one real change --- the platform
decides a service is idle from its \emph{outbound} traffic, and a pooled database connection is
outbound traffic, so the default connection pool would have kept the service awake and billable
for ever.
\section{What is still wrong}
An honest introduction should end with the open problems, not the achievements.
\begin{itemize}
\item \textbf{The model has no features in its training data.} Its only two remaining inputs,
CADD and AlphaMissense, are absent from all 688,362 training rows. It abstains today, so
nothing is broken --- but installing the VEP plugins would flip it on while it still knows
nothing. Serving should refuse unless the model was trained on the features it is handed.
\item \textbf{The benchmark is contaminated} and no amount of code fixes it. A leave-one-publication-out
rebuild is possible --- HPO's annotation file carries the source PMID ---
and is the honest next step.
\item \textbf{No baseline comparison.} Scoring the same held-out rows with CADD and AlphaMissense
would give the model something to beat. AlphaMissense is a 0.64\,GB download and tractable;
CADD's whole-genome file is 87.5\,GB and is not.
\item \textbf{Rarity contributes nothing without the cache} in the default demonstration, though
the gnomAD streaming route above now makes that solvable.
\end{itemize}
\vspace{1.4em}
\hrule
\vspace{0.9em}
{\small\color{soft}
rarelens is a self-training project built in the open on public data, licensed AGPL-3.0. It is not
a clinical tool and makes no diagnostic claim: ACMG/AMP treats computational predictions as
supporting evidence only, never sufficient alone. Data sources, licences, citations and the
evaluation caveats are documented in \texttt{docs/data.md}; the architecture decisions, including
why Google Cloud rather than AWS, are in \texttt{docs/architecture.md} and \texttt{docs/cloud.md}.}
\end{document}