feat(ingest): enrich mirrored GenBank metadata - #6981
Conversation
|
This PR may be related to: #2834 ("Include extra metadata from Genbank that is in the full Genbank files but not (yet) parsed and emitted by NCBI virus"), which requests exactly the kind of GenBank metadata (not present in NCBI Virus/Datasets) that this PR's Entrez enrichment sidecar adds. |
theosanderson-agent
left a comment
There was a problem hiding this comment.
Went over this with a focus on silent data corruption, since the output is a persistent cache that feeds a public mirror. One thing I'd treat as a blocker (the base-accession fallback), plus a cluster of robustness issues in the new retry path.
Good news first: the not_found retry at line 259 does close the cache-poisoning hole. I checked it against a fake Entrez that returns <eFetchResult><ERROR>ID list is empty!</ERROR></eFetchResult> — the miss is re-requested on the next run rather than being frozen into the uploaded JSONL forever. Worth keeping a test on that, because with the workflow doing s3cmd get / s3cmd put on the same key, a regression there is permanent in object storage rather than merely annoying.
The blocker is the base-accession fallback: a request for MN908947.3 can be answered with MN908947.4's metadata and written out labelled as .3, with no not_found and no warning. Details inline.
The rest, in rough priority order:
- Version mislabelling via the base-accession fallback (line 276) — silent, and it lands in the cache.
- Network errors that bypass the new retry logic entirely (line 206) — read-phase failures are the expected mode with multi-MB
gbwithpartsresponses. - All-or-nothing runs (line 284) — the output file is opened only after every batch succeeds, so one failure at batch 9,000 discards the previous 8,999. Combined with 2, a large taxon may never bootstrap its cache.
- A malformed
Retry-Afterheader raises inside the 429 handler (line 168). - A truncated cache file wedges all future runs (line 70).
All reproductions were run offline against a monkeypatched urlopen — no calls to NCBI.
| for requested in batch: | ||
| record = fetched.get(requested) | ||
| if record is None: | ||
| record = fetched.get(requested.split(".", 1)[0]) |
There was a problem hiding this comment.
This can attach a different version's metadata to the requested accession, silently.
When fetched.get("MN908947.3") misses, this falls back to the base accession, and line 280 then stamps requested_accession="MN908947.3" onto a record whose accession_version is something else entirely. Nothing compares the two.
This is the normal case, not an edge case: if the Datasets zip was cut before a version bump (or is simply a few hours stale), nuccore serves the current version, and EFetch resolves a versionless/superseded id to it.
Reproduced against this PR's head — input FASTA asks for MN908947.3, fake Entrez returns only MN908947.4:
requested_accession : MN908947.3
accession_version : MN908947.4
definition : THIS IS VERSION 4
strain : ['v4-strain']
not_found present : False
Exit 0, no warning. Any consumer joining on requested_accession — which is the documented key — now has v4's strain/host/collection dates attached to the v3 sequence in the mirror, and the mislabelled row is written straight into the JSONL that becomes next run's cache.
Suggest either dropping the fallback, or keeping it but only when the versions actually agree:
record = fetched.get(requested)
if record is None:
candidate = fetched.get(requested.split(".", 1)[0])
# Only accept a base-accession hit if it really is the requested version.
if candidate and candidate.get("accession_version") == requested:
record = candidateIf you'd rather keep the looser behaviour, it at least needs to record what was actually returned (e.g. a version_mismatch flag) so downstream can tell.
| records = {} | ||
| for sequence in root.findall("INSDSeq") + root.findall("GBSeq"): | ||
| record = parse_record(sequence) | ||
| for key in (record["accession_version"], record["accession_base"]): |
There was a problem hiding this comment.
Root cause of the fallback problem above: every record is indexed under both its version and its bare accession, so the base key is inherently ambiguous.
If one response contains two versions of the same accession, the base key holds whichever parsed last. Observed with MN908947.2 and MN908947.4 in a single response:
fetched keys: ['MN908947', 'MN908947.2', 'MN908947.4']
fetched['MN908947'] -> MN908947.4
So a request for MN908947.1 in that batch receives .4. Dropping the bare-accession key (or storing a list and refusing ambiguous hits) removes the whole class of problem.
| file=sys.stderr, | ||
| ) | ||
| time.sleep(delay) | ||
| except http.client.IncompleteRead as error: |
There was a problem hiding this comment.
The retry loop covers HTTPError 429 and IncompleteRead, but the failures most likely to occur here go straight through it and through the outer except (urllib.error.URLError, ET.ParseError) at line 270.
ET.parse(response) reads the body after urlopen returned, so read-phase errors are not wrapped in URLError. With timeout=90 against multi-MB gbwithparts responses, that's the expected failure mode. Tested each against this head:
TimeoutError *** ESCAPES uncaught ***
ConnectionResetError *** ESCAPES uncaught ***
SSLError *** ESCAPES uncaught ***
RemoteDisconnected *** ESCAPES uncaught ***
Each one is a raw traceback and a lost run rather than the intended RuntimeError message.
Separately, when IncompleteRead does exhaust its 10 retries it re-raises, and IncompleteRead is not a URLError either — so the exhaustion path also escapes line 270 uncaught. Widening the retry to except (http.client.IncompleteRead, TimeoutError, ConnectionError, ssl.SSLError) as error: and adding those to the outer handler would cover it.
Minor, same loop: the backoff is min(2 ** retry, 60) per attempt, so ten IncompleteRead retries is ~5 minutes of sleeping on a batch that may simply be too large — probably worth capping total retry time, or halving the batch size on truncation instead of retrying it unchanged.
| try: | ||
| return max(0.0, float(value)) | ||
| except ValueError: | ||
| retry_at = email.utils.parsedate_to_datetime(value) |
There was a problem hiding this comment.
On Python 3.10+, parsedate_to_datetime raises ValueError on an unparseable date rather than returning None — so a malformed Retry-After raises inside the except ValueError block that was meant to handle exactly that case, and the whole run dies at the moment we're already being rate-limited.
*** raises ValueError inside the 429 handler -> run aborts ***
(The if retry_at is None guard below is therefore also unreachable.) Wrapping the parse in its own try/except (ValueError, TypeError): return None restores the intended "fall back to exponential backoff" behaviour.
| records.append(record) | ||
| time.sleep(pause) | ||
|
|
||
| output = sys.stdout if args.output == "-" else open(args.output, "w", encoding="utf-8") |
There was a problem hiding this comment.
The output file is opened only after every batch has succeeded, and fetch_batch has no retry for anything except 429/truncation — so a single unrecoverable failure discards every batch already fetched.
Reproduced: 250 accessions, batch size 100, HTTP 429 injected past the retry budget on the 3rd batch → 2 batches (200 records) fetched, output file exists: False. Nothing to upload, and the next run starts from zero.
That matters more than usual given the workflow only uploads on success: for a large taxon a cold run has ~10⁴ batches, and at ≥0.35 s each plus fetch time it's already brushing the 6-hour job limit. If it can't complete in one clean pass, the cache never bootstraps at all.
Streaming each batch to the output as it completes (and uploading whatever was written) would make partial progress durable, and would turn the next run's cache lookup into the natural resume point.
| for line_number, line in enumerate(handle, 1): | ||
| if not line.strip(): | ||
| continue | ||
| record = json.loads(line) |
There was a problem hiding this comment.
json.loads per line with no error handling, against a file fetched with s3cmd get. Combined with the plain open(..., "w") at line 284 (no temp-file-and-rename), a run interrupted mid-write — or a partial S3 download — leaves a half-written final line:
load_previous_enrichment -> JSONDecodeError: Unterminated string starting at: line 1 column 39
Every subsequent run then hard-fails until someone manually deletes the object, and the thousands of valid records earlier in the file aren't salvaged. Skipping (and warning about) an unparseable trailing line, plus writing via .tmp + os.replace, makes this self-healing.
|
|
||
|
|
||
| EUTILS_URL = "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi" | ||
| ACCESSION_RE = re.compile(r"\b([A-Za-z]{1,6}_?\d+(?:\.\d+)?)\b") |
There was a problem hiding this comment.
Worth a comment noting this regex is only safe for the datasets download virus genome header format.
It can't match a two-part RefSeq prefix — after NZ_ come letters, not digits, and _ is a word character so \b blocks matching CP007592.1 mid-token. search() then either returns None (header skipped entirely — no row at all, not even not_found) or matches a token out of the free-text description:
None <- >NZ_CP007592.1 Escherichia coli chromosome, complete genome
'contig_1' <- >NZ_JAAABC010000001.1 Vibrio cholerae strain X contig_1, whole genome shotgun sequence
'NC_045512' <- >lcl|NC_045512.2_cds_YP_009724389.1_1 [gene=ORF1ab]
I checked the forms this workflow actually produces (NC_045512.2, MN908947.3, AB012345.1, JAAABC010000001.1, AAAA01000001.1) and they all match correctly, so this is not a live bug for datasets-mirror.yml today. It does bite the --assembly-source refseq output (NZ_*) and the FASTA-file input mode the CLI advertises, and it fails open — silent data loss, no error. Given the .zip/FASTA/accession-list input modes are all documented as supported, either tightening the regex or erroring on an unmatched > line would be safer than the current silent skip.
| query = { | ||
| "db": "nuccore", | ||
| "id": ",".join(accessions), | ||
| "rettype": "gbwithparts", |
There was a problem hiding this comment.
gbwithparts returns the full sequence, which the parser never reads — and for CON/WGS records it forces NCBI to assemble the sequence server-side, which is the slowest EFetch path and the one most likely to hit the 90 s timeout.
Measured on a synthetic 100-record SARS-CoV-2-sized batch: 3.0 MB response, 15.9 kB of fields actually retained (~190×). rettype=gb gives the same annotations without the sequence body. (I checked the memory angle too — ET.parse peaks at only ~1.1× the response size, so this is bandwidth and wall-clock, not OOM.)
| if record is None: | ||
| record = fetched.get(requested.split(".", 1)[0]) | ||
| if record is None: | ||
| record = {"requested_accession": requested, "not_found": True} |
There was a problem hiding this comment.
Minor schema wrinkle: a not_found row has 2 keys where every other row has 15, and both land in the same JSONL. A consumer doing record["division"] or record["organism"] gets a KeyError on these rather than a null. Emitting the full key set with null values (plus not_found: true) would make the file uniformly parseable.
corneliusroemer
left a comment
There was a problem hiding this comment.
Nice, some thoughts:
- Ingest would ideally see the consolidated (NCBI Virus + Entrez) JSONL to abstract away the details
- Entrez enrichment would always run on new sequences before they become available to ingest - only entrez enriched sequences would be added to JSONL used by ingest
- Entrez enrichment results cache can be keyed on a hash of NCBI virus JSON output - when the virus output changes then entrez would be able to refetch
- In addition we could run hard refreshes of cache periodically - potentially on rolling basis so as not to overwhelm Github actions.
Summary
This adds a GenBank/Entrez metadata sidecar to the NCBI Datasets mirror. The existing NCBI Datasets ZIP and compressed archive remain the source of sequence data, while
mirroring/enrich_genbank.pyfetches the corresponding nuccore records and writes a taxon-specific.entrez.jsonlfile containing metadata that may be absent from NCBI Virus/Datasets.The enrichment parser supports both INSDSeq and GBSeq XML responses and records accession identifiers, definition, organism, taxonomy, division, dates, source qualifiers, and authors. Existing enrichment output is downloaded from object storage and reused by accession version, while prior misses are retried so transient failures and parser improvements can fill them later.
Reliability
Entrez requests now:
Retry-Aftersupport;IncompleteRead);The workflow uploads the base dataset ZIP and tar archive before starting enrichment. This ensures the ordinary mirror is updated even if enrichment later fails or exceeds the Actions time limit. The enrichment JSONL is uploaded separately after it completes.
Validation
python -m unittest discover -s mirroring -p 'test_*.py' -v— 4 tests pass.3048148produced 15,800 records in 8m37s and recovered from a deliberately observed truncated upstream response.3048148completed successfully end-to-end.197911) mirror completed successfully with batch size 100, where batch size 20 had previously exceeded the six-hour job limit.Impact
Consumers can join the new JSONL sidecar to mirrored accessions for richer INSDC/GenBank metadata. Existing ZIP and tar artifact names and locations are unchanged, and their upload no longer depends on enrichment succeeding.
🚀 Preview: Add
previewlabel to enable