fcf788d7357f6f207f2bb0b61acde9dc0507d89f
lrnassar
  Mon Sep 21 04:04:10 2026 -0700
Fix two stale figures and guard the accession interpolation in the MaveMD build, per CR. refs #37800

The makeDoc script listing still said the variant converter writes bed12+31. It writes
bed12+34, as the autoSql, mavemd.ra, runBuild.sh and the makeDoc's own build-results
line all already said.

mavemdLib.py's docstring claimed the codon projection is cross-checked against "~154k
variants that carry both a genomic and a protein term". 154,895 is the number placed by
the genomic route; the cross-check set is the smaller number that also has a resolvable
protein term. The docstring now describes the set rather than quoting a figure that
drifts with every build.

Also adds checkAccession() and calls it before either query that interpolates an
accession into SQL. Nothing can currently reach those queries with a quote in it, since
the accessions come from PROTEIN_TERM whose character class excludes one, but the regex
is a hundred lines from the query and a later edit to it should not be able to open this
up silently. Output is byte-identical to the previous build.

diff --git src/hg/makeDb/doc/hg38/mavemd.txt src/hg/makeDb/doc/hg38/mavemd.txt
index 3201e1d233a..c2b3da58f48 100644
--- src/hg/makeDb/doc/hg38/mavemd.txt
+++ src/hg/makeDb/doc/hg38/mavemd.txt
@@ -1,228 +1,228 @@
 # [Claude/lrnassar] MaveMD - clinically calibrated multiplexed variant effect measurements (2026-09-17)
 # refs #37800
 
 # MaveMD ("MAVEs for MeDicine") is not a separate database. It is a curated collection inside
 # MaveDB holding the score sets that carry clinical evidence calibrations, so everything comes
 # from the ordinary MaveDB API. The collection URN is
 #   urn:mavedb:collection-603dafbf-4a3f-4d70-ab8c-aafb226fbff4
 # which is not documented; it was read out of the MaveDB front end, which fetches it to render
 # the mavedb.org/mavemd page.
 # Reference: McEwen et al. (2025) medRxiv, PMID 41332838. A local copy of the paper is in
 # paper/ alongside this build:
 #   mkdir -p paper && cd paper
 #   curl -s "https://www.ebi.ac.uk/europepmc/webservices/rest/PMC12668102/fullTextXML" -o paper.xml
 #   # paper.md is that XML flattened to markdown (title, abstract, sections) for grepping
 #
 # Two things from the paper shaped the build:
 #  1. The paper counts a "dataset" at the MaveDB experiment level, not the score set level.
 #     The collection is 84 score sets drawn from 81 experiments; three experiments contributed
 #     two score sets each. Say "score set" in the docs, not "dataset".
 #  2. Counts are measurement-level. A variant measured by several assays contributes several
 #     measurements. Our 460,490 items are 275,875 distinct ClinGen alleles, which is why a
 #     codon studied by several groups stacks over a hundred items deep.
 
 # Request pacing is what MaveDB asked us for directly (Benjamin Capodanno, 21 Aug 2026 email):
 # metadata at no more than two concurrent requests, score and variant endpoints sequential,
 # page at 100k variants, and back off rather than retry immediately on a 504, because the
 # gateway gives up while their query keeps running. fetchMaveMd.py implements all of this.
 
 mkdir -p /hive/data/outside/mavemd/2026-09-17
 mkdir -p /hive/data/genomes/hg38/bed/mavemd/2026-09-17
 
 # Download: the collection, one metadata record per score set, and one flat CSV per score set
 # from score-sets/{urn}/variants/data. The CSV joins the assay score to the post-mapped HGVS
 # terms, ClinGen allele IDs, gnomAD frequencies, the most recent ClinVar release, and every
 # calibration's functional class and ACMG call. All calibrations are requested, research-use-only
 # ones included; only the newest of the twelve available ClinVar releases is requested.
 ~/kent/src/hg/makeDb/scripts/mavemd/fetchMaveMd.py /hive/data/outside/mavemd/2026-09-17
 
 #   MaveDB API 2026.2.7.2
 #   Collection MaveMD: 84 score sets, modified 2026-04-30
 #   Downloaded 242.1 MB in 257 seconds
 
 # Build both tracks: the per-variant bigBed and the per-score-set heatmap.
 ~/kent/src/hg/makeDb/scripts/mavemd/runBuild.sh \
     /hive/data/outside/mavemd/2026-09-17 \
     /hive/data/genomes/hg38/bed/mavemd/2026-09-17
 
 # Scripts (all in ~/kent/src/hg/makeDb/scripts/mavemd/):
 #   fetchMaveMd.py         download the collection from the API
 #
 # GENCODE VERSION: mavemdLib.py hardcodes GENCODE_ATTRS/GENCODE_GENEPRED to the V50 tables,
 # used only for the six ENSP-stated score sets. Bump both when hg38 moves to a newer GENCODE;
 # the RefSeq route (ncbiRefSeqLink/ncbiRefSeqCurated) is unversioned and needs no attention.
 #   mavemdLib.py           shared: codon projection, color tables, BED field hygiene
-#   makeMaveMdVariants.py  per-variant bed12+31 -> mavemdVar.bb
+#   makeMaveMdVariants.py  per-variant bed12+34 -> mavemdVar.bb
 #   makeMaveMdHeatmap.py   per-score-set bed12+20 -> mavemdMap.bb
 #   mavemdVariants.as, mavemdHeatmap.as
 #   runBuild.sh            driver for the two converters plus bedToBigBed
 
 # COORDINATES. MaveDB's VRS mapper resolves about a third of the collection to the genome and
 # the rest only to a protein sequence, so placement takes three routes:
 #   154,895  MaveDB's own genomic HGVS term, run through kent's hgvsToVcf
 #   296,516  projected from the protein term onto its codon
 #     9,079  the submitter's original transcript term, run through hgvsToVcf
 #   460,490  placed, of 481,717 measurements (95.6%)
 #
 # The codon projection resolves the protein accession to its transcript and walks that CDS to
 # the requested codon. MaveDB states protein terms against RefSeq for most score sets and
 # Ensembl for six, so both routes exist: NP_ through ncbiRefSeqLink + ncbiRefSeqCurated, ENSP
 # through wgEncodeGencodeAttrsV50 + wgEncodeGencodeCompV50. All 40 protein accessions resolve.
 #
 # 2,986 items sit in a codon split across an exon junction, so its three bases are not a 3 bp
 # run. Those are written as BED12 with the two real blocks, so the intron draws as a thin
 # connector rather than as a solid item up to 95 kb wide. The heatmap column is drawn on the
 # longest contiguous run of the codon's bases (CodonMap.codonBlock), keeping the block off
 # the intron.
 #
 # CROSS-CHECK. 130,955 measurements carry both a genomic and a protein term, so the projection
 # is run on those too and must land on the genomic coordinate. 41 disagreed (0.0313%), all
 # nonsense terms (p.XxxNNNTer) where MaveDB's protein term names a codon one position off from
 # the nucleotide change; mostly the ENSP-mapped score sets 00000662-0-1 and 00000673-0-1.
 # makeMaveMdVariants.py aborts above 0.5%.
 
 # NOT PLACED (21,227 of 481,717):
 #   18,963  haplotypes: several substitutions measured as one unit, no single position.
 #           All from urn:mavedb:00000665-a/b/c/d-1. The earlier MaveDB track excluded
 #           haplotype datasets for the same reason.
 #    1,961  submitted term rejected by hgvsToVcf. These were qualified with a transcript and
 #           submitted; they are mostly PTEN NM_000314.8:c.79_80delins... terms spanning an
 #           exon junction, which hgvsToVcf expands to a ~2 kb REF and then rejects on a
 #           reference mismatch. Dropping them is correct.
 #      114  genomic term rejected by hgvsToVcf (terms like NC_000003.12:g.=)
 #      126  no usable HGVS term of any kind (includes rows named "_wt")
 #       63  protein term is not a single substitution
 
 # Eleven score sets state their variants as bare terms with no accession (c.34_36delinsCTT).
 # MaveDB leaves those unmapped, but their neighbours in the same score set carry a post-mapped
 # cDNA term naming the transcript, so the bare terms are qualified from it and then placed.
 # This is what recovered most of the 9,079 hgvsToVcf placements.
 
 # HEATMAP. All 84 score sets produce a map. Rows are the 20 standard amino acids ordered by
 # class (A V L I M F Y W R H K D E S T N Q G C P), matching the MaveDB and popEVE heatmap
 # tracks, plus a final row for nonsense. A synonymous measurement fills the reference
 # residue's cell.
 #   427,471  measurements with a single amino acid substitution
 #   385,732  cells drawn. 41,739 measurements are a later nucleotide change encoding a
 #            substitution already measured at that position. The strongest evidence wins the
 #            cell (CSV order is not a defined rule and flipped 611 cells between pathogenic
 #            and benign); the cell is marked '+', or '!' where the measurements disagree in
 #            direction. 1,865 cells disagree.
 #
 # One calibration colors a whole map, or its cells are not comparable to each other. 16 score
 # sets have no MaveDB primary and their rows fall to different calibrations, so the one that
 # classifies the most variants wins and every cell is recolored to it.
 #
 # The renderer splits the score array with chopCommas (keeps a trailing empty field) but the
 # label array with chopByCharRespectDoubleQuotesKeepEmpty (drops one). If the last cell of the
 # last row is empty the counts disagree and the track aborts, so the converter puts a
 # placeholder label on that one cell. 40 score sets needed it. Same trap as popEVE.
 
 # ASSAY METADATA. The experiment record's controlled keywords carry the fields McEwen et al.
 # argue a clinician judges an assay by, so four of them ride on every item: Phenotypic Assay
 # Method, Phenotypic Assay Model System, Molecular Mechanism Assessed, Variant Library Creation
 # Method. The last matters most: an in vitro construct library introduces a synthetic copy of
 # the target, so it cannot see an effect on splicing or NMD and can read falsely normal for
 # one; an endogenous locus library can. 61 score sets use construct libraries, 16 endogenous,
 # 7 unstated. The keyword value is in keyword.label, not keyword.value.
 #
 # Not exposed by the API: the paper's per-assay "can this assay detect splicing / NMD
 # variants" curation, and the rest of the ~180 curated fields. Only the controlled keywords
 # above come through, so the library route is the closest available proxy.
 #
 # NON-ASCII. MaveDB free text carries en dashes and accented author names (Gronbaek-Thygesen
 # and Olvera-Leon are spelled with o-slash, ae and o-acute upstream). The browser does not
 # transcode UTF-8, so these reach the details page as
 # mojibake. mavemdLib.bedField() converts every non-ASCII character to a numeric HTML entity;
 # 35,694 variant rows and 6 heatmap rows needed it. Re-check after any rebuild with:
 #   LC_ALL=C grep -cP '[^\x00-\x7F]' mavemdVar.bed mavemdMap.bed     # must be 0
 
 # SCORE DIRECTION. Functional scores have no common scale and no common direction across score
 # sets; the paper calls out MSH2, where LOWER scores mean normal activity. Nothing in the build
 # compares raw scores across score sets, and the colors come from the functional class and the
 # ACMG code rather than from the score, so this is a documentation problem only. It is called
 # out in a warn-note on mavemdVar.html.
 
 # COLORS. Cells and items are colored by the ACMG evidence code on an 11-class red-to-blue
 # diverging ramp, red toward pathogenic and blue toward benign, darkening with strength. The
 # strength ladder is MaveDB's own StrengthOfEvidenceProvided enum, taken from their OpenAPI
 # spec, not assumed: VERY_STRONG, STRONG, MODERATE_PLUS, MODERATE, SUPPORTING. Score sets that
 # classify variants without assigning an ACMG code use a separate purple-green palette, so a
 # measurement is not mistaken for a clinical claim.
 
 # WHICH CALIBRATION IS SHOWN. A score set can carry up to six calibrations. MaveDB curates a
 # `primary` flag with its own promote and demote API endpoints, so where that primary actually
 # classifies the variant it is their editorial choice and is used as-is: 170,197 items.
 #
 # Two traps. 21 score sets have no primary at all. A further 33 designate one that carries no
 # functionalClassifications, so it says nothing about any variant and arrives as NA on every
 # CSV row - it exists but cannot be displayed. Those are different situations and the item's
 # "Calibration chosen by" field distinguishes them, because telling a user "no primary
 # designated" about a score set that has one is simply false. 171,760 items are in that second
 # case; an earlier build mislabelled all of them.
 #
 # Every calibration's call is listed on the details page regardless of which one is displayed.
 
 # Build results:
 #   mavemdVar.bb   460,490 items, 46 fields (bed12+34)
 #
 # PROVENANCE OF THE ANNOTATIONS. clinvarRelease and gnomadVersion ride on each item that has a
 # ClinVar or gnomAD annotation (60,737 and 3,979 respectively), so a stale annotation is visible
 # rather than silent. This build used ClinVar 2026-01 and gnomAD v4.1, whichever the API served;
 # fetchMaveMd.py always asks for the newest ClinVar namespace the score set offers, so the
 # release moves on its own between builds and is not pinned anywhere.
 
 # LINKOUTS. The urls setting linkifies publicationUrl, scoreSet (to its MaveDB page) and
 # clinGenId (to the ClinGen Allele Registry, which serves both CA genomic and PA protein
 # alleles). hgc percent-encodes the substituted value, so the score set URN arrives as
 # urn%3Amavedb%3A...%2D0%2D1; mavedb.org resolves that fine. There is no separate mavedbUrl
 # field: it used to store 460,490 copies of 84 distinct strings that are just a prefix plus
 # scoreSet, so linkifying scoreSet replaces it.
 #
 # DEFAULT DISPLAY. mavemdMap pack at priority 1, mavemdVar dense at priority 2. The variant
 # track in pack is ~5,900 px tall over a well studied gene, so it is not a reasonable thing to
 # open onto; dense gives a compact colored strip and the user switches to pack once zoomed to
 # a few codons or filtered to one gene.
 #   mavemdMap.bb        84 items, 32 fields (bed12+20)
 #
 # Install: the bigDataUrl and dataVersion settings point at /gbdb, so the links have to be
 # repointed at each new dated build directory. runBuild.sh does not do this.
 mkdir -p /gbdb/hg38/mavemd
 cd /gbdb/hg38/mavemd
 ln -sf /hive/data/genomes/hg38/bed/mavemd/2026-09-17/mavemdVar.bb mavemdVar.bb
 ln -sf /hive/data/genomes/hg38/bed/mavemd/2026-09-17/mavemdMap.bb mavemdMap.bb
 ln -sf /hive/data/genomes/hg38/bed/mavemd/2026-09-17/version.txt  version.txt
 
 # Load trackDb. "make DBS=hg38" writes the per-user sandbox tables (trackDb_$USER); the track
 # is gated "include mavemd.ra alpha" so "make alpha DBS=hg38" is what puts it on dev.
 cd ~/kent/src/hg/makeDb/trackDb
 make DBS=hg38
 #
 # STRAND comes from the transcript of the protein term, so it is set on every item that has a
 # usable one (427,471), not only on the ones the codon projection ultimately placed: 186,869 +,
 # 240,602 -. The remaining 33,019 are '.' - 23,940 placed from a genomic HGVS term with no
 # protein term, and 9,079 placed by hgvsToVcf from the submitted term.
 
 # Version file for the dataVersion trackDb setting.
 cd /hive/data/genomes/hg38/bed/mavemd/2026-09-17
 python3 -c "
 import json
 c=json.load(open('/hive/data/outside/mavemd/2026-09-17/collection.json'))
 api=open('/hive/data/outside/mavemd/2026-09-17/apiVersion.txt').read().split()[-1]
 open('version.txt','w').write('MaveDB collection modified %s (MaveDB API %s)\n'%(c['modificationDate'],api))
 "
 
 # trackDb: human/hg38/mavemd.ra, included from trackDb.ra as "include mavemd.ra alpha".
 # The filterValues blocks in that file are generated by makeMaveMdVariants.py --raFragment;
 # paste mavemdFilters.ra in rather than maintaining the value lists by hand. Paste it INSIDE
 # the stanza with no blank lines: a blank line ends a .ra stanza, so blank separators between
 # the filter groups silently orphan every setting after the first one. hgTrackDb loads without
 # complaint and the filters just never appear in hgTrackUi. Check after a make with:
 #   hgsql hg38 -Ne "select settings from trackDb_$USER where tableName='mavemdVar'" \
 #     | tr '\\' '\n' | grep -c filterValues ClinVar terms
 # containing a comma are escaped by doubling it, which is what slNameListFromCommaEscaped wants.
 # No filterType is set on those fields on purpose: the list types run COMPARE_HASH_LIST_OR in
 # hgTracks/bigBedTrack.c, which splits the field VALUE on commas before looking it up, so a
 # ClinVar term with a comma could never match its own menu entry. The default whole-field
 # match is correct here and still multi-select.
 # relatedTracks.ra carries reciprocal links between mavedb and mavemd.