Write the wiki: 14 pages plus a sidebar
Operational and scientific knowledge that does not belong in the README: architecture and why the pieces sit where they do, the pipeline and variant identity, the three execution backends, Argo and the GitOps loop, Terraform and cost, local development, data provenance, the ranking and its abstention rule, both benchmarks with their caveats, testing and CI, deployments including Railway's constraints, a Gotchas page, and a Roadmap that states the known gaps rather than hiding them. Every number came from a run in the repository and names the command that reproduces it.
+57
@@ -0,0 +1,57 @@
|
||||
# Architecture
|
||||
|
||||
One monorepo holds the scientific pipeline, the API, the interface, model serving and the
|
||||
infrastructure for two deployment tracks.
|
||||
|
||||
```
|
||||
scientist ──▶ SvelteKit UI ──▶ FastAPI ──▶ PostgreSQL
|
||||
(/api proxy) │ ▲
|
||||
│ │
|
||||
├──▶ MLflow registry (model)
|
||||
│
|
||||
└──▶ Nextflow pipeline ──▶ bcftools ──▶ VEP ──▶ loader
|
||||
```
|
||||
|
||||
The interface proxies `/api` from its own origin, so there is no CORS anywhere and one public
|
||||
address serves both.
|
||||
|
||||
## Why the pieces are where they are
|
||||
|
||||
**The pipeline is separate from the API on purpose.** Annotation takes minutes to hours and needs
|
||||
containers, reference data and a lot of memory; the API needs to answer in milliseconds. They
|
||||
communicate through a `jobs` row and nothing else, which is what lets the same pipeline run as a
|
||||
local subprocess, a Cloud Run job or an Argo Workflow without the API knowing.
|
||||
See [Execution backends](Execution-backends).
|
||||
|
||||
**Feature engineering travels with the model.** The registered MLflow artifact carries
|
||||
`rarelens_ml.features` as model code, so serving cannot drift from training. The API sends raw
|
||||
columns and gets a probability back; it does not know what the features are.
|
||||
|
||||
**The ontology work happens at load time, not query time.** HPO annotations are propagated up the
|
||||
ontology and each term's information content is computed once by `scripts/load-hpo.py`, so the
|
||||
ranking is a set lookup and a weighted sum. See [Ranking](Ranking).
|
||||
|
||||
## Data model
|
||||
|
||||
```
|
||||
cases ──┬── case_phenotypes the patient's HPO terms
|
||||
└── jobs ──── variants ──┬── predictions model score
|
||||
└── variant_decisions shortlist / dismiss + reason
|
||||
gene_phenotypes HPO gene→term, propagated (reference data)
|
||||
hpo_terms term → information content (reference data)
|
||||
```
|
||||
|
||||
`jobs` also records **what the annotation run produced** — `has_frequencies`,
|
||||
`has_effect_scores`. That is not bookkeeping: the ranking uses it to decide which components are
|
||||
allowed to score at all. See [Ranking](Ranking).
|
||||
|
||||
## Deployment tracks
|
||||
|
||||
| Track | What runs | When to use it |
|
||||
|---|---|---|
|
||||
| **Serverless** (default) | Cloud Run for api + web, Cloud Run job for the Nextflow driver, Google Batch for pipeline tasks | almost always; idles near £1/month |
|
||||
| **Kubernetes** (flagged off) | GKE Autopilot, Argo Workflows, Argo Events, ArgoCD | to demonstrate the GitOps path, then destroy |
|
||||
| **Railway** | Postgres + api + web, pre-seeded, no pipeline | putting it in front of people quickly |
|
||||
|
||||
Both cloud tracks run identical pipeline code; the executor is a Nextflow profile.
|
||||
See [Deployments](Deployments) and [Infrastructure](Infrastructure).
|
||||
+77
@@ -0,0 +1,77 @@
|
||||
# Benchmarks
|
||||
|
||||
Two separate measurements: how well the **phenotype ranking** retrieves the right gene, and how
|
||||
well the **model** classifies pathogenicity. They are not the same thing and should never be
|
||||
quoted as one.
|
||||
|
||||
## Phenotype ranking
|
||||
|
||||
```bash
|
||||
make benchmark # downloads phenopacket-store, exports gene_phenotypes, runs the ranking
|
||||
```
|
||||
|
||||
Given a real patient's reported HPO terms, where does the gene their authors diagnosed rank among
|
||||
all 5,269 HPO-annotated genes? Ties are reported as a range: optimistic counts a tie as a win,
|
||||
pessimistic counts every tied gene as ahead.
|
||||
|
||||
| cases | | top-1 | top-10 | MRR |
|
||||
|---|---|---|---|---|
|
||||
| all 10,178 | optimistic | 81.0% | 87.2% | 0.830 |
|
||||
| | pessimistic | 45.9% | 69.7% | 0.541 |
|
||||
| the 6,485 with ≥6 terms | optimistic | 77.1% | 85.4% | 0.797 |
|
||||
| | pessimistic | 59.5% | 81.0% | 0.670 |
|
||||
|
||||
Random guessing would be 0.02%.
|
||||
|
||||
### The benchmark is contaminated
|
||||
|
||||
HPO's gene-to-phenotype annotations are curated **from these same case reports**. The median
|
||||
causal gene already carries every one of its patient's terms. So this measures how well the
|
||||
ranking retrieves a gene HPO has already been told about — an upper bound. A prospective number,
|
||||
on a patient whose gene nobody has annotated yet, would be lower, and this corpus cannot say by
|
||||
how much.
|
||||
|
||||
Do not quote the headline without this paragraph.
|
||||
|
||||
### A change that did not work
|
||||
|
||||
Information-content weighting and ontology propagation both replaced plain term counting. Asked
|
||||
whether they helped — pessimistic figures, the 6,485 cases with six or more terms:
|
||||
|
||||
| scoring | top-1 | top-10 | MRR |
|
||||
|---|---|---|---|
|
||||
| count terms (the original) | 61.8% | 80.1% | 0.682 |
|
||||
| **+ information content** | **63.6%** | **83.5%** | **0.706** |
|
||||
| + propagation | 58.2% | 77.7% | 0.653 |
|
||||
| + both (shipped) | 59.5% | 81.0% | 0.670 |
|
||||
|
||||
Weighting earns its place. Propagation costs about what weighting gains. It was kept anyway,
|
||||
because this corpus **cannot show what propagation is for** — its term IDs were chosen by the same
|
||||
curators HPO records, so exact matching is flattered, while real users pick their own terms from a
|
||||
search box. That is an argument, not a measurement, and it is labelled as one.
|
||||
|
||||
## The model
|
||||
|
||||
```bash
|
||||
make training-set && make train
|
||||
```
|
||||
|
||||
Held out by gene, never by variant: a random split puts variants of the same gene on both sides and
|
||||
the model scores the gene instead of the variant (Grimm et al., *Hum Mutat* 2015).
|
||||
312,025 training and 74,239 test variants; 7,728 and 1,932 genes; no gene on both sides.
|
||||
|
||||
| | AUROC | AUPRC | missense AUROC | missense AUPRC |
|
||||
|---|---|---|---|---|
|
||||
| v2, with gnomAD allele frequency | 0.986 | 0.954 | 0.872 | 0.725 |
|
||||
| v3, allele frequency removed | 0.966 | 0.881 | **0.500** | 0.398 |
|
||||
|
||||
**0.500 is exactly random.** Strip frequency out and nothing is left but the consequence class, so
|
||||
every missense variant scores identically. The respectable-looking 0.872 was allele frequency, not
|
||||
variant-effect knowledge — and ACMG's BA1/BS1 criteria assign ClinVar's *benign* labels using
|
||||
frequency, so the feature had partly caused the label.
|
||||
|
||||
The overall 0.966 is the same trick one level up: ClinVar's pathogenic set is largely loss of
|
||||
function and its benign set largely is not.
|
||||
|
||||
This is why the model abstains unless it has CADD or AlphaMissense. See [Roadmap](Roadmap) for why
|
||||
that is not yet a complete fix.
|
||||
+55
@@ -0,0 +1,55 @@
|
||||
# Data sources
|
||||
|
||||
Everything is public, peer-reviewed and consented for open redistribution. **No patient data, no
|
||||
data access agreement, nothing that needs an application.** Full citations with DOIs are in
|
||||
`docs/data.md` in the repository — quote from there, not from here.
|
||||
|
||||
| Source | Used for | Licence |
|
||||
|---|---|---|
|
||||
| ClinVar | training labels, the ClinVar column | NCBI public domain |
|
||||
| Genome in a Bottle HG002 | the background genome | open, no restriction |
|
||||
| gnomAD v4 | allele frequency | free use |
|
||||
| Human Phenotype Ontology | what the phenotype half matches against | free with attribution |
|
||||
| Phenopacket Store (Monarch) | the published demonstration case, and the benchmark | BSD-3-Clause |
|
||||
| Ensembl VEP 113 | annotation | Apache 2.0 |
|
||||
| CADD, AlphaMissense | optional model features | see `docs/data.md` — CADD needs a commercial licence |
|
||||
|
||||
## Loading it
|
||||
|
||||
```bash
|
||||
make data # GIAB + ClinVar, chr22 only, ~100 MB
|
||||
make hpo # HPO annotations, propagated and IC-weighted
|
||||
make demo-case # a simulated proband: GIAB background + one ClinVar pathogenic variant
|
||||
make published-case # a real published patient, from a GA4GH phenopacket
|
||||
make training-set # ~370k labelled variants from ClinVar, 2-star and above
|
||||
```
|
||||
|
||||
## The demonstration cases
|
||||
|
||||
**The simulated proband** spikes one real ClinVar pathogenic variant into a real public genome.
|
||||
It has a right answer, so the ranking can be checked rather than admired.
|
||||
|
||||
**The published case** invents nothing. It reads a GA4GH phenopacket curated from a peer-reviewed
|
||||
case report and takes the patient's reported HPO terms and the authors' causal variant verbatim.
|
||||
The default is the *TGFBR2* proband from Loeys et al., *Nat Genet* 2005 — the paper that first
|
||||
defined Loeys–Dietz syndrome — with 30 reported terms and `NM_003242.6:c.1069G>T`.
|
||||
|
||||
The rest of that patient's genome is not public, and triage means nothing if the causal variant is
|
||||
the only variant in the file, so background variants come from GIAB HG002 around the locus, drawn
|
||||
from **coding exons** where possible. That last detail matters more than it sounds: of roughly
|
||||
4,000 HG002 variants in a 3 Mb window only 9 are coding, so a random sample is entirely intronic,
|
||||
the consequence filter discards all of it, and the causal variant is left as the only candidate —
|
||||
a funnel that proves nothing.
|
||||
|
||||
## Allele frequencies without the 25 GB cache
|
||||
|
||||
VEP's `--af_gnomade` is rejected outright with `--database`, and plain `--af` returns nothing even
|
||||
for common variants. But gnomAD's public bucket is tabix-indexed, so a range request works:
|
||||
|
||||
```bash
|
||||
bcftools view -r chr3:30672000-30673000 \
|
||||
https://storage.googleapis.com/gcp-public-data--gnomad/release/4.1/vcf/exomes/...chr3.vcf.bgz
|
||||
# 392 records, 5 KB -> feed to VEP with --custom alongside --database
|
||||
```
|
||||
|
||||
Verified end to end. Not yet wired into the pipeline — see [Roadmap](Roadmap).
|
||||
+64
@@ -0,0 +1,64 @@
|
||||
# Deployments
|
||||
|
||||
Three ways to run it. Pick on cost and on whether you need the pipeline.
|
||||
|
||||
| | Pipeline runs? | Idle cost | Effort |
|
||||
|---|---|---|---|
|
||||
| [GCP serverless](#gcp-serverless-default) | yes | ~£1/month | terraform apply |
|
||||
| [GCP Kubernetes](Kubernetes-and-GitOps) | yes | control plane + Cloud SQL | flag, then destroy |
|
||||
| [Railway](#railway) | **no** | $5/month flat | `railway up` |
|
||||
|
||||
## GCP serverless (default)
|
||||
|
||||
```bash
|
||||
cd infra/terraform
|
||||
terraform init -backend-config="bucket=<tfstate bucket>"
|
||||
export TF_VAR_database_url='postgresql+asyncpg://user:pass@host/db?sslmode=require'
|
||||
terraform apply -var project=<project id>
|
||||
cd ../.. && make serverless-deploy PROJECT=<project id> TAG=<commit sha>
|
||||
```
|
||||
|
||||
`terraform output web_url` is the address to share: it serves the UI and proxies `/api`, so there
|
||||
is one public address and no CORS. Upload the VEP cache to `gs://<project>-rarelens-data/refs/vep`
|
||||
before a real annotation, and set `-var model_uri=gs://.../models/pathogenicity/1` to score without
|
||||
an MLflow server.
|
||||
|
||||
**Set a billing budget first.** The demo has no authentication.
|
||||
|
||||
## Railway
|
||||
|
||||
Three services — Postgres, api, web — with the API on the private network only, so the UI's `/api`
|
||||
proxy is the single public entry point.
|
||||
|
||||
```bash
|
||||
railway link --project <id> --environment production --service api
|
||||
cd api && railway up --service api # the repo is on Gitea, so deploy from the working copy
|
||||
cd ../web && railway up --service web
|
||||
```
|
||||
|
||||
Variables: on the API, `DATABASE_URL` pointing at `postgres.railway.internal` and
|
||||
`DATABASE_IDLE_CONNECTIONS=false`; on the UI, `API_INTERNAL_URL=http://api.railway.internal:8000`,
|
||||
`PUBLIC_PIPELINE_ENABLED=false`, `BASIC_AUTH_USER` and `BASIC_AUTH_PASSWORD`.
|
||||
|
||||
Three things to know:
|
||||
|
||||
- **The pipeline cannot run there.** Nextflow shells out to `docker run` for VEP and bcftools, and
|
||||
a container platform gives you a container, not a Docker daemon. Cases are annotated elsewhere
|
||||
and copied up with `scripts/seed-remote.sh`; `PUBLIC_PIPELINE_ENABLED=false` hides the buttons
|
||||
that would otherwise be left to fail.
|
||||
- **`DATABASE_IDLE_CONNECTIONS=false` is what makes sleeping work.** Railway decides a service is
|
||||
idle from its *outbound* traffic, and a pooled database connection is outbound traffic, so the
|
||||
default pool keeps the service awake and billable for ever. It uses `NullPool` instead — a
|
||||
connection per request, which is nothing at demo traffic and the wrong trade under real load.
|
||||
- **Serverless must be enabled per service and only applies on the next deploy.** Leave it off for
|
||||
Postgres, which holds the volume. Expect the documented 502 on the first request to a slept
|
||||
service; it answers on retry.
|
||||
|
||||
## Seeding a remote database
|
||||
|
||||
```bash
|
||||
scripts/seed-remote.sh "postgresql://user:pass@host:port/db"
|
||||
```
|
||||
|
||||
Copies the analysed cases and the HPO reference data. `--data-only`: the target's schema comes from
|
||||
Alembic, so the two cannot disagree about it. Run `alembic upgrade head` against the target first.
|
||||
@@ -0,0 +1,34 @@
|
||||
# Execution backends
|
||||
|
||||
`POST /cases/{id}/annotate` writes a `jobs` row and hands off. `app/services/events.launch()`
|
||||
then picks a backend from configuration alone — the pipeline code is identical in all three cases.
|
||||
|
||||
| Setting present | Backend | Behaviour |
|
||||
|---|---|---|
|
||||
| `cloudrun_job` | Cloud Run job | executes the Nextflow driver with argument overrides via the Jobs API |
|
||||
| `pubsub_topic` | Pub/Sub | publishes the job; Argo Events triggers an Argo Workflow |
|
||||
| neither | local subprocess | runs Nextflow directly, streaming stdout into the job log |
|
||||
|
||||
They are checked in that order.
|
||||
|
||||
## Why three
|
||||
|
||||
**Cloud Run job** is the default because it scales to zero: nothing runs between pipelines. The
|
||||
service account holds `run.jobsExecutorWithOverrides` on exactly one job, not project-wide.
|
||||
|
||||
**Pub/Sub** is the decoupled path. Retries, ordering and back-pressure become the queue's problem
|
||||
instead of the API's, and the API can restart mid-pipeline without losing work. This is the right
|
||||
shape once more than one thing produces annotation requests.
|
||||
|
||||
**The local subprocess** is what a developer gets with nothing configured. `_watch` streams the
|
||||
last `PROGRESS_LINES` of Nextflow's output into the job row every `PROGRESS_INTERVAL_S` seconds so
|
||||
the UI can show live progress, and `_record_progress` only writes while the job is still `running`,
|
||||
so a late write cannot resurrect a finished job.
|
||||
|
||||
All three converge on the same `jobs` row. The interface polls one endpoint and does not know or
|
||||
care which ran.
|
||||
|
||||
## When the API has no Nextflow
|
||||
|
||||
If `shutil.which("nextflow")` finds nothing, the job is marked failed with the command to run by
|
||||
hand. That is the normal state in docker-compose, and on Railway — see [Deployments](Deployments).
|
||||
+70
@@ -0,0 +1,70 @@
|
||||
# Gotchas
|
||||
|
||||
Things that cost real debugging time here. Each one is a bug that was found and fixed; they are
|
||||
recorded because the symptom rarely points at the cause.
|
||||
|
||||
## Science
|
||||
|
||||
**A missing measurement is not a measurement of zero.** `rarity_score(None)` read "this run has no
|
||||
allele frequencies" as "absent from gnomAD, therefore maximally rare" and handed every variant a
|
||||
free 0.25. The fix is structural: jobs record what they looked up, and a component with no
|
||||
evidence abstains. If you add a component, decide what makes it abstain first. See
|
||||
[Ranking](Ranking).
|
||||
|
||||
**A DAG walked pre-order and read backwards loses ancestors.** Computing each HPO term's ancestors
|
||||
with a pre-order DFS reversed *looks* like a post-order, but on a DAG a term can be reached before
|
||||
one of its parents, and then inherits that parent alone instead of the parent's whole lineage. It
|
||||
silently dropped **399 terms** out of the phenotype branch — including Camptodactyly and Chiari
|
||||
malformation — so cases mentioning them scored lower for no visible reason. Now a true post-order,
|
||||
tested against a reference transitive closure.
|
||||
|
||||
**A feature can cause its own label.** The model took allele frequency as a feature, while ACMG's
|
||||
BA1/BS1 criteria assign ClinVar's *benign* labels using allele frequency. It was rediscovering the
|
||||
rule that generated its training data, and the impressive missense AUROC of 0.872 fell to exactly
|
||||
0.500 when the feature was removed. See [Benchmarks](Benchmarks).
|
||||
|
||||
**Split by gene, not by variant.** A random split puts variants of the same gene on both sides and
|
||||
the model learns the gene. Grimm et al. documented this inflation for exactly this class of tool.
|
||||
|
||||
## Pipeline
|
||||
|
||||
**VEP's coordinates are not a primary key.** `Location` and `Allele` trim indel alleles and shift
|
||||
positions; `CT>C` comes back as `->` elsewhere. Identity travels in the VCF ID instead. See
|
||||
[Pipeline](Pipeline).
|
||||
|
||||
**`--af_gnomade` is rejected with `--database`,** and plain `--af` silently returns nothing even
|
||||
for a variant at 15% global frequency. Frequencies genuinely need the cache — or a range request
|
||||
against gnomAD's public bucket, see [Data sources](Data-sources).
|
||||
|
||||
**Nextflow rejects two staged inputs sharing a filename.** The cache placeholder and the plugin
|
||||
placeholder must be different files (`assets/NO_CACHE` and `assets/NO_FILE`), or the run fails with
|
||||
an input collision that names neither.
|
||||
|
||||
**`head` plus `pipefail` kills a script silently.** `head` closes the pipe, the upstream process
|
||||
takes SIGPIPE, and with `pipefail` the whole script aborts with no message.
|
||||
|
||||
## Infrastructure
|
||||
|
||||
**Kustomize only substitutes hashed ConfigMap names in workloads it considers in scope.** Missing
|
||||
`namespace:` in an overlay means the hash is generated but never substituted, and pods mount a name
|
||||
that does not exist. The symptom is `CreateContainerConfigError` with correct-looking manifests.
|
||||
|
||||
**A connection pool stops a container ever sleeping.** Platforms that sleep idle containers measure
|
||||
*outbound* traffic; a pooled database connection is outbound traffic. See [Deployments](Deployments).
|
||||
|
||||
## Local environment
|
||||
|
||||
**macOS AirPlay Receiver owns port 5000** and answers 403, which looks exactly like an auth failure
|
||||
from MLflow. MLflow is on 5001 here.
|
||||
|
||||
**`*.rlwy.net` may resolve to 0.0.0.0** on a machine with DNS filtering, so Railway's TCP proxy
|
||||
looks refused when it is fine. Check with a public resolver and connect by IP.
|
||||
|
||||
## Security boundaries, not conveniences
|
||||
|
||||
**`vcf_uri` validation is a boundary.** A case's URI must be a `gs://` object or an absolute path
|
||||
beneath `LOCAL_DATA_ROOT`, ending in a VCF suffix. Without it, a crafted value becomes a Nextflow
|
||||
option or reads an arbitrary file. It has its own tests; do not relax it for convenience.
|
||||
|
||||
**The database URL never goes on a command line.** Environment or Nextflow secret only, so it stays
|
||||
out of `.command.sh` and the workflow logs.
|
||||
+46
-1
@@ -1,3 +1,48 @@
|
||||
# rarelens
|
||||
|
||||
Initialising the wiki.
|
||||
Rare-disease variant triage on public data: a proband's variants narrowed against their phenotype,
|
||||
each candidate carrying the evidence for its rank.
|
||||
|
||||
**This is a research and self-training platform. It is not a clinical tool and makes no diagnostic
|
||||
claim. No patient data is used, accepted, or possible to load.**
|
||||
|
||||
## Start here
|
||||
|
||||
| If you want to… | Read |
|
||||
|---|---|
|
||||
| understand how the pieces fit together | [Architecture](Architecture) |
|
||||
| run it on your machine | [Local development](Local-development) |
|
||||
| know what the ranking actually does | [Ranking](Ranking) |
|
||||
| know how well it works, and how well it doesn't | [Benchmarks](Benchmarks) |
|
||||
| change the pipeline | [Pipeline](Pipeline) |
|
||||
| deploy it | [Deployments](Deployments) |
|
||||
| avoid a trap someone already hit | [Gotchas](Gotchas) |
|
||||
|
||||
## What it does
|
||||
|
||||
1. A **case** is a proband: a VCF plus the patient's phenotype as HPO terms.
|
||||
2. A **Nextflow** pipeline normalises the VCF, annotates it with Ensembl VEP, and loads the
|
||||
variants into PostgreSQL.
|
||||
3. The **API** ranks each variant on four auditable components — phenotype fit, rarity,
|
||||
consequence severity, and a model score — and refuses to score any component the annotation run
|
||||
did not actually measure.
|
||||
4. The **interface** shows the narrowing as a funnel and the evidence behind each candidate. A
|
||||
reviewer shortlists or dismisses with a reason, and that trail becomes a case report.
|
||||
|
||||
## The shape of the repository
|
||||
|
||||
```
|
||||
api/ FastAPI, Pydantic v2, SQLAlchemy 2.0 async, Alembic (~1,380 lines)
|
||||
ml/ LightGBM model, HPO ontology handling, benchmark (~500 lines)
|
||||
pipeline/ Nextflow DSL2, bcftools + VEP + loader (~620 lines)
|
||||
web/ SvelteKit 5, adapter-node (~1,460 lines)
|
||||
infra/ Terraform (GCP), Kustomize, Argo Workflows, ArgoCD (~890 lines)
|
||||
scripts/ data preparation and seeding (~505 lines)
|
||||
docs/ architecture, cloud choice, data provenance, write-up
|
||||
```
|
||||
|
||||
## A note on how this wiki is written
|
||||
|
||||
Pages record what is true now, including what is broken. Where a number appears it came from a run
|
||||
in this repository and can be reproduced by the command next to it. If you change something that
|
||||
makes a number here wrong, change the number too — a wiki nobody trusts is worse than no wiki.
|
||||
|
||||
@@ -0,0 +1,54 @@
|
||||
# Infrastructure
|
||||
|
||||
`infra/terraform` provisions both deployment tracks from one root module, with the expensive half
|
||||
behind flags.
|
||||
|
||||
```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
|
||||
```
|
||||
|
||||
## What the default track creates
|
||||
|
||||
| Resource | Why |
|
||||
|---|---|
|
||||
| Cloud Run services (api, web) | `min-instances=0`, so nothing runs when nobody is looking |
|
||||
| Cloud Run **job** (nextflow driver) | exists only while a pipeline runs |
|
||||
| Google Batch | VEP tasks on Spot VMs |
|
||||
| GCS bucket | VCFs, the VEP cache, model artifacts |
|
||||
| Secret Manager | the database URL, and any Nextflow secrets |
|
||||
| Artifact Registry | images CI pushes |
|
||||
| Service accounts + IAM | one per component, least privilege |
|
||||
|
||||
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.
|
||||
|
||||
## Cost as a design constraint
|
||||
|
||||
A portfolio platform is idle more than 99% of the time, so idle cost is the only cost that matters.
|
||||
The reasoning, including why Google Cloud rather than AWS, is in `docs/cloud.md`. The short version:
|
||||
GKE's free tier covers one control plane where EKS charges about $73/month, and Cloud Run's
|
||||
scale-to-zero plus Batch on Spot puts the idle bill near £1/month — almost all of it the database.
|
||||
|
||||
The Kubernetes track is meant to be turned on, demonstrated, and destroyed:
|
||||
|
||||
```bash
|
||||
make serverless-destroy PROJECT=<project id>
|
||||
```
|
||||
|
||||
## Authentication from CI
|
||||
|
||||
Workload Identity Federation. **No service account key exists in the repository or in repository
|
||||
secrets.** CI exchanges its OIDC token for short-lived credentials.
|
||||
|
||||
## Configuration that matters
|
||||
|
||||
| Variable | Effect |
|
||||
|---|---|
|
||||
| `deploy_kubernetes` | GKE, Argo, the whole second track |
|
||||
| `deploy_cloud_sql` | managed Postgres instead of an external URL |
|
||||
| `database_url` | an external Postgres (a free tier, say) for the serverless track |
|
||||
| `model_uri` | score from an artifact, with no MLflow server running |
|
||||
|
||||
Terraform is `fmt -check`ed and validated on every pull request — see [Testing and CI](Testing-and-CI).
|
||||
@@ -0,0 +1,46 @@
|
||||
# Kubernetes and GitOps
|
||||
|
||||
Off by default. `terraform apply -var deploy_kubernetes=true` turns it on; it exists to
|
||||
demonstrate the GitOps path and is meant to be destroyed afterwards.
|
||||
|
||||
## Layout
|
||||
|
||||
```
|
||||
infra/k8s/base/ api, web, postgres, ingress
|
||||
infra/k8s/overlays/local/ kind: in-cluster Postgres, NodePort
|
||||
infra/k8s/overlays/gcp/ Cloud SQL, Workload Identity, image tags
|
||||
infra/argo-workflows/ annotate WorkflowTemplate, EventSource, Sensor, RBAC
|
||||
infra/argocd/app.yaml the Application that reconciles the cluster
|
||||
```
|
||||
|
||||
## Locally
|
||||
|
||||
```bash
|
||||
make kind # builds the images, loads them into a kind cluster, applies the local overlay
|
||||
```
|
||||
|
||||
## The GitOps loop
|
||||
|
||||
1. A pull request merges to `main`.
|
||||
2. CI builds and pushes images to Artifact Registry.
|
||||
3. CI bumps the image tags in `infra/k8s/overlays/gcp`.
|
||||
4. ArgoCD sees the commit and reconciles the cluster.
|
||||
|
||||
Deployment is therefore a commit, and rollback is a revert.
|
||||
|
||||
## Argo Workflows and Argo Events
|
||||
|
||||
`annotate.yaml` is a `WorkflowTemplate` running the same containers the Nextflow pipeline uses.
|
||||
`events.yaml` holds the Pub/Sub `EventSource` and the `Sensor` that triggers the template when the
|
||||
API publishes a job. `rbac.yaml` grants the sensor permission to create workflows in the
|
||||
`rarelens` namespace and nothing else.
|
||||
|
||||
## The ConfigMap trap
|
||||
|
||||
Kustomize appends a content hash to generated ConfigMap names so that a config change forces a
|
||||
rollout. That substitution only happens for workloads Kustomize considers in scope. Both overlays
|
||||
were missing `namespace: rarelens`, so the hashed name was generated but never substituted, and
|
||||
pods mounted a ConfigMap name that no longer existed.
|
||||
|
||||
The symptom is a pod stuck in `CreateContainerConfigError` while every manifest reads correctly.
|
||||
See [Gotchas](Gotchas).
|
||||
@@ -0,0 +1,51 @@
|
||||
# Local development
|
||||
|
||||
```bash
|
||||
make up # postgres + api + web + mlflow via docker-compose
|
||||
make migrate # alembic upgrade head
|
||||
make hpo # HPO annotations: what the ranking matches against
|
||||
make published-case # a real published patient
|
||||
make test # api, ml, loader and web tests
|
||||
```
|
||||
|
||||
Then open <http://localhost:5173>.
|
||||
|
||||
## Making the "Analyse case" button work
|
||||
|
||||
docker-compose has no Nextflow, so the button marks the job failed with the command to run instead.
|
||||
To get the real thing, run the API **on the host**, where Nextflow and Docker are:
|
||||
|
||||
```bash
|
||||
docker compose up -d db
|
||||
cd api && DATABASE_URL=postgresql+asyncpg://rarelens:rarelens@localhost:5432/rarelens \
|
||||
PIPELINE_DATABASE_URL=postgresql+asyncpg://rarelens:[email protected]:5432/rarelens \
|
||||
LOCAL_DATA_ROOT=$PWD/.. VEP_DATABASE=true \
|
||||
uv run --extra dev uvicorn app.main:app --port 8000
|
||||
```
|
||||
|
||||
- `PIPELINE_DATABASE_URL` is what the loader container gets. Inside it, the API's own `localhost`
|
||||
would be the container itself.
|
||||
- `LOCAL_DATA_ROOT` is the directory a case's `vcf_uri` must sit under. This is a security
|
||||
boundary, not a convenience — see [Gotchas](Gotchas).
|
||||
|
||||
## Ports
|
||||
|
||||
| Service | Port | Note |
|
||||
|---|---|---|
|
||||
| web | 5173 | |
|
||||
| api | 8000 | |
|
||||
| postgres | 5432 | |
|
||||
| MLflow | **5001** | not 5000: macOS AirPlay Receiver owns 5000 and answers 403 |
|
||||
|
||||
## Training a model
|
||||
|
||||
```bash
|
||||
make training-set # ClinVar-derived table, ~370k labelled variants
|
||||
make train # fits, reports held-out metrics by gene split, moves the production alias
|
||||
```
|
||||
|
||||
The API scores with `models:/rarelens-pathogenicity@production`. Setting `MODEL_URI` to an artifact
|
||||
path skips the registry entirely, which is how the serverless deployment scores with no MLflow
|
||||
server running.
|
||||
|
||||
What the metrics do and do not mean is in [Benchmarks](Benchmarks). Read that before quoting them.
|
||||
+69
@@ -0,0 +1,69 @@
|
||||
# Pipeline
|
||||
|
||||
`pipeline/main.nf` is a three-process Nextflow DSL2 workflow. Each process runs in a pinned
|
||||
container.
|
||||
|
||||
| Process | Container | Does |
|
||||
|---|---|---|
|
||||
| `NORMALISE` | bcftools 1.20 | left-align, split multi-allelics, set the VCF ID to `CHROM_POS_REF_ALT` |
|
||||
| `VEP` | ensembl-vep 113.0 | consequence, gene, HGVS, ClinVar; optionally CADD and AlphaMissense |
|
||||
| `LOAD_DB` | python + psycopg | parse VEP `--tab`, insert variants, mark the job succeeded |
|
||||
|
||||
## Running it
|
||||
|
||||
```bash
|
||||
# dry run: annotate without touching the database
|
||||
make pipeline VCF=data/example.vcf.gz
|
||||
|
||||
# annotate for a real job created in the UI
|
||||
make annotate JOB=<job id> VCF=data/example.vcf.gz
|
||||
|
||||
# no 25 GB VEP cache? query Ensembl's public database instead
|
||||
VEP_DATABASE=true make pipeline VCF=pipeline/tests/data/tiny.vcf
|
||||
```
|
||||
|
||||
Database mode takes roughly 25–35 seconds per variant and returns no allele frequencies and no
|
||||
plugin scores, so keep those runs to tens of variants. What that costs the ranking is described in
|
||||
[Ranking](Ranking).
|
||||
|
||||
## Variant identity
|
||||
|
||||
**VEP's `Location` and `Allele` columns cannot be used as a key.** They trim indel alleles and
|
||||
shift positions, so `CT>C` comes back as `->` at a different coordinate. `NORMALISE` therefore
|
||||
writes identity into the VCF ID, which VEP passes through untouched:
|
||||
|
||||
```bash
|
||||
bcftools annotate --set-id '%CHROM\_%POS\_%REF\_%FIRST_ALT'
|
||||
```
|
||||
|
||||
and `load_db.parse_variant_id` reads it back. It `rsplit`s on `_` from the right, because contig
|
||||
names can contain underscores (`chrUn_KI270742v1`) while positions and alleles cannot.
|
||||
|
||||
Verified on a real run: `22:42126611 CT>C` survives with its alleles intact.
|
||||
|
||||
## Executors
|
||||
|
||||
The workflow never names an executor. `pipeline/nextflow.config` maps profiles:
|
||||
|
||||
- **local / docker** — a laptop, or the API's subprocess backend
|
||||
- **gcp** — Google Batch, Spot VMs, needs `--project` and `--bucket`
|
||||
- **Argo** — the same containers as workflow steps, see [Kubernetes and GitOps](Kubernetes-and-GitOps)
|
||||
|
||||
## Stub runs
|
||||
|
||||
Every process has a `stub:` block, so CI can run the whole workflow with no containers, no cache
|
||||
and no database:
|
||||
|
||||
```bash
|
||||
cd pipeline && nextflow run main.nf -stub-run --vcf tests/data/tiny.vcf
|
||||
```
|
||||
|
||||
This runs on every pull request. It checks wiring and channel shapes — a renamed output, a process
|
||||
emitting the wrong cardinality — which is the class of breakage that otherwise only shows up in a
|
||||
two-hour annotation run.
|
||||
|
||||
## Secrets
|
||||
|
||||
The database URL is passed by environment or as a Nextflow secret, **never on a command line**, so
|
||||
it does not land in `.command.sh` or the workflow logs. `PIPELINE_DATABASE_URL` exists because the
|
||||
loader runs inside a container where the API's own `localhost` would be the container itself.
|
||||
+76
@@ -0,0 +1,76 @@
|
||||
# Ranking
|
||||
|
||||
`api/app/services/triage.py`. The rank is a **weighted mean of four components a reviewer can
|
||||
audit**, not a black box.
|
||||
|
||||
| Component | Weight | Scores when |
|
||||
|---|---|---|
|
||||
| phenotype fit | 0.35 | always |
|
||||
| rarity in gnomAD | 0.25 | the run looked up frequencies |
|
||||
| consequence severity | 0.20 | always |
|
||||
| model P(pathogenic) | 0.20 | the run has CADD or AlphaMissense |
|
||||
|
||||
**ClinVar is deliberately not an input.** It is displayed beside the result as independent
|
||||
confirmation, so nothing ranks highly merely because ClinVar already called it pathogenic.
|
||||
|
||||
Rarity and consequence **filter** — the usual first pass. Phenotype only **ranks**, because a real
|
||||
diagnosis can sit in a gene nobody has annotated yet, and filtering on phenotype would hide exactly
|
||||
that case.
|
||||
|
||||
## Evidence that was never looked up abstains
|
||||
|
||||
This is the most important rule in the file, and it replaced a real bug.
|
||||
|
||||
Run without a VEP cache there are no allele frequencies. `rarity_score(None)` read that as *absent
|
||||
from gnomAD, therefore maximally rare* and gave every variant a free 0.25. The model separately
|
||||
returned the same 0.887 for every variant, from features it had never been given. Two of four
|
||||
components were fiction and the total looked fully informed.
|
||||
|
||||
Now `jobs.has_frequencies` and `jobs.has_effect_scores` record what the run actually produced, a
|
||||
component with no evidence returns `None` rather than a number, and `combine()` renormalises the
|
||||
weights over whatever is left. The score stays on 0–1 and means the same thing; the UI prints
|
||||
"not looked up" instead of drawing a bar.
|
||||
|
||||
```python
|
||||
def combine(components):
|
||||
weight = sum(WEIGHTS[n] for n, v in components.items() if v is not None)
|
||||
return sum(WEIGHTS[n] * v for n, v in components.items() if v is not None) / weight
|
||||
```
|
||||
|
||||
**If you add a component, decide what makes it abstain before you decide its weight.**
|
||||
|
||||
## Phenotype fit
|
||||
|
||||
Information-content-weighted recall: the share of the *total specificity* of the patient's terms
|
||||
that this gene accounts for.
|
||||
|
||||
```
|
||||
score = Σ IC(matched terms) / Σ IC(all the patient's terms)
|
||||
```
|
||||
|
||||
`IC(term) = -ln(fraction of annotated genes carrying it)`, computed at load time. It matters:
|
||||
"Dilated left subclavian artery" scores 7.88, "Global developmental delay" 0.93 — an 8.5× gap that
|
||||
plain term counting threw away.
|
||||
|
||||
A term HPO has never annotated to any gene gets `DEFAULT_IC`, keeping it in the denominator so it
|
||||
depresses every gene equally. That is the neutral choice, not an oversight.
|
||||
|
||||
## The model
|
||||
|
||||
Trained on ClinVar, held out **by gene** rather than by variant. Allele frequency is deliberately
|
||||
not a feature — it double-counted against `rarity_score` and was circular, since ACMG assigns
|
||||
ClinVar's benign labels using frequency. Removing it dropped missense AUROC from 0.872 to **0.500,
|
||||
exactly random**, which is why the model now abstains without CADD or AlphaMissense.
|
||||
Full working in [Benchmarks](Benchmarks).
|
||||
|
||||
## Worked example
|
||||
|
||||
The published Loeys–Dietz case, run in VEP database mode:
|
||||
|
||||
| | score | phenotype (0.64) | rarity | consequence (0.36) | model |
|
||||
|---|---|---|---|---|---|
|
||||
| *TGFBR2* 3:30672252 missense | **0.855** | 1.00 (30/30) | not looked up | 0.60 | not looked up |
|
||||
| *OSBPL10* 3:31748090 missense | 0.218 | 0.00 | not looked up | 0.60 | not looked up |
|
||||
|
||||
Both are rare missense variants, identical on every piece of evidence this run holds except one.
|
||||
The phenotype is what separates a published diagnosis from an incidental variant.
|
||||
+50
@@ -0,0 +1,50 @@
|
||||
# Roadmap
|
||||
|
||||
Known gaps, most consequential first. Nothing here is hidden in an issue tracker.
|
||||
|
||||
## 1. 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** — `make training-set` builds from the ClinVar VCF, which carries neither.
|
||||
|
||||
Today the model abstains, so nothing is broken. But installing the VEP plugins flips
|
||||
`has_effect_scores` true, the model stops abstaining and takes 0.20 of every rank **while still
|
||||
knowing nothing**. That is the same class of bug as the one the abstention rule fixed, one step
|
||||
further out.
|
||||
|
||||
*Fix:* log at training time which features actually had values, and refuse to score at serving
|
||||
unless the model was trained on the features it is handed.
|
||||
|
||||
## 2. AlphaMissense in the training set
|
||||
|
||||
0.64 GB, tabix-indexed, tractable. It would close gap 1 and give the model something real to
|
||||
contribute — missense is exactly where it currently scores 0.500. CADD's whole-genome file is
|
||||
87.5 GB and is not worth it; its web API is an option for small volumes.
|
||||
|
||||
## 3. gnomAD frequencies without the cache
|
||||
|
||||
Proven to work (see [Data sources](Data-sources)) but not wired into the pipeline. It would make
|
||||
the rarity component real in the default demonstration instead of abstaining.
|
||||
|
||||
*Shape:* a pipeline step that streams the regions covering the case's variants into a small local
|
||||
slice, then `vep --custom` against it alongside `--database`. Note VEP's Perl HTS cannot open the
|
||||
remote index directly, and gnomAD's contigs need renaming to match.
|
||||
|
||||
## 4. Decontaminating the benchmark
|
||||
|
||||
HPO's `phenotype.hpoa` carries the source PMID for each annotation, so a **leave-one-publication-
|
||||
out** rebuild is possible: score each case against a gene-to-phenotype table with that case's own
|
||||
paper excluded. That would turn an upper bound into something much closer to a prospective number.
|
||||
|
||||
## 5. A baseline to beat
|
||||
|
||||
Scoring the same held-out rows with CADD and AlphaMissense would put the model's numbers next to
|
||||
something meaningful rather than leaving them unanchored.
|
||||
|
||||
## 6. Smaller things
|
||||
|
||||
- ClinVar review status is enforced at training (2-star and above) but is not shown in the UI,
|
||||
because VEP's `CLIN_SIG` does not carry `CLNREVSTAT`. A `--custom` annotation would fix it.
|
||||
- The deployed demo has one shared credential and no per-user accounts.
|
||||
- The funnel numbers come from tens of variants. A real exome has thousands of rare coding
|
||||
variants; nothing here has been run at that scale.
|
||||
@@ -0,0 +1,45 @@
|
||||
# Testing and CI
|
||||
|
||||
```bash
|
||||
make test # api (99 Python tests), ml, pipeline loader, web (49)
|
||||
make lint # ruff + mypy for the API, svelte-check for the front end
|
||||
```
|
||||
|
||||
## What is tested, and why there
|
||||
|
||||
Tests are concentrated where a mistake would be silent rather than loud.
|
||||
|
||||
- **`api/tests/test_triage.py`** — the ranking arithmetic, as pure logic. This is the scientific
|
||||
claim the app makes, so it is tested without a database in the way. Includes the regression for
|
||||
components scoring evidence nobody looked up.
|
||||
- **`ml/tests/test_hpo.py`** — ontology propagation and information content, including a property
|
||||
test against a reference transitive closure. That test exists because a subtly wrong graph walk
|
||||
silently deleted 399 terms; see [Gotchas](Gotchas).
|
||||
- **`pipeline/tests/test_load_db.py`** — the loader, against a real PostgreSQL. Idempotency is
|
||||
asserted directly: loading twice must not duplicate variants.
|
||||
- **`api/tests/test_security.py`** — the `vcf_uri` boundary.
|
||||
- **`api/tests/test_migrations.py`** — Alembic upgrade *and* downgrade against a live database.
|
||||
|
||||
Database tests use an embedded PostgreSQL (`pgserver`) locally and a service container in CI,
|
||||
because the schema is part of the behaviour and SQLite would not exercise it.
|
||||
|
||||
## CI
|
||||
|
||||
Runs on pull requests and pushes to `main`:
|
||||
|
||||
| Job | Does |
|
||||
|---|---|
|
||||
| api | ruff, mypy, pytest against a Postgres service |
|
||||
| ml | pytest |
|
||||
| web | `npm run check`, vitest, production build |
|
||||
| pipeline | loader tests against Postgres, then `nextflow run main.nf -stub-run` |
|
||||
| terraform | `fmt -check` and `validate` |
|
||||
| images | builds and pushes to Artifact Registry, `main` only |
|
||||
|
||||
CI authenticates to Google Cloud with Workload Identity Federation. **There is no service account
|
||||
key in the repository or in repository secrets.**
|
||||
|
||||
## A caution
|
||||
|
||||
CI does not run on branch pushes, only on pull requests and `main`. Work on a long-lived branch is
|
||||
therefore unverified by CI until the PR is opened — run `make test && make lint` locally.
|
||||
+23
@@ -0,0 +1,23 @@
|
||||
### rarelens
|
||||
|
||||
**Understanding it**
|
||||
- [Home](Home)
|
||||
- [Architecture](Architecture)
|
||||
- [Ranking](Ranking)
|
||||
- [Benchmarks](Benchmarks)
|
||||
- [Data sources](Data-sources)
|
||||
|
||||
**Working on it**
|
||||
- [Local development](Local-development)
|
||||
- [Pipeline](Pipeline)
|
||||
- [Execution backends](Execution-backends)
|
||||
- [Testing and CI](Testing-and-CI)
|
||||
|
||||
**Running it**
|
||||
- [Deployments](Deployments)
|
||||
- [Kubernetes and GitOps](Kubernetes-and-GitOps)
|
||||
- [Infrastructure](Infrastructure)
|
||||
|
||||
**When it goes wrong**
|
||||
- [Gotchas](Gotchas)
|
||||
- [Roadmap](Roadmap)
|
||||
Reference in New Issue
Block a user