6577d5ee1319bbea85988c1c89179436b4a94edf lrnassar Tue Jul 14 11:27:59 2026 -0700 Address code-review feedback on the Cardiomyopathy VCEP build scripts. refs #37446 - cmpVCEPCardioBoost.py: add the standard --db/--output-dir CLI. It previously hardcoded the working directory for both its input TSV and its output (unlike the 11 sibling scripts, and contrary to the makedoc's documented interface); the build loop's flags were silently ignored. Output is unchanged (31,236 variants per assembly). - Decode leftover HTML entities (arrows, >=, <=, +/-, x) in print/stderr diagnostics, comments, and docstrings across all scripts so build logs read cleanly. The mouseOver / bigBed display strings intentionally keep their entities. - cmpVCEPWalsh2019.py: fix the stale docstring that described the ClinVar-unmatched entries as "deferred" (they are mapped via the hgvsToVcf fallback, item L) and drop the unverified "163 rows" count. Per code-review feedback on commit aa5669fe64. No track data changed. diff --git src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPPVS1.py src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPPVS1.py index e8cbef504da..52c62423626 100644 --- src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPPVS1.py +++ src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPPVS1.py @@ -1,55 +1,55 @@ #!/usr/bin/env python3 """ -B.2 — MYBPC3 PVS1 Evidence track builder (MYBPC3-only, the lone PVS1-applicable gene). +B.2 - MYBPC3 PVS1 Evidence track builder (MYBPC3-only, the lone PVS1-applicable gene). Per GN095 (verified at A.1): - NMD-escape region: codons 1254+ (50 nt upstream of last exon-exon junction at exon 33:34) PVS1 may not apply to LoF variants in this region (they escape NMD) - Micro-exons 10, 11, 14: in silico splice predictions less reliable - In-frame exons 2-4, 8-11, 14, 20, 22, 24-27: in-frame deletions may not be true LoF (these encode functionally important domains but consequences vary) Outputs (per-assembly bigBed): cmpVCEPPVS1/cmpVCEPPVS1.as cmpVCEPPVS1/cmpVCEPPVS1Hg{38,19}.bed + .bb Total features: ~15 (1 NMD-escape region + 14 in-frame/micro exon caveats) -The 3 micro-exons (10, 11, 14) are tagged as in-frame too — single feature per exon +The 3 micro-exons (10, 11, 14) are tagged as in-frame too - single feature per exon with both caveats listed in the mouseover. """ import argparse, os, subprocess, sys # Re-use B.1's MANE parsing helpers sys.path.insert(0, os.path.dirname(os.path.abspath(__file__))) from cmpVCEPClinDomains import parse_mane_record, cds_exons, aa_to_genomic_segments -# MYBPC3 GN095 PVS1 caveats — hand-transcribed (will be re-verified at D.0) +# MYBPC3 GN095 PVS1 caveats - hand-transcribed (will be re-verified at D.0) # In-frame exons (1-based exon numbering per CDS) IN_FRAME_EXONS = [2, 3, 4, 8, 9, 10, 11, 14, 20, 22, 24, 25, 26, 27] # Micro-exons (subset of in-frame; splice prediction unreliable) MICRO_EXONS = {10, 11, 14} -# NMD-escape: codons 1254+ — PVS1 should not be applied here +# NMD-escape: codons 1254+ - PVS1 should not be applied here NMD_ESCAPE_CODON_START = 1254 -PVS1_NMD_ESCAPE_COLOR = '136,136,136' # gray — "PVS1 not applicable" -PVS1_MICRO_EXON_COLOR = '255,140,40' # orange — splice unreliable -PVS1_INFRAME_EXON_COLOR = '210,80,40' # red-orange — caveat for in-frame del +PVS1_NMD_ESCAPE_COLOR = '136,136,136' # gray - "PVS1 not applicable" +PVS1_MICRO_EXON_COLOR = '255,140,40' # orange - splice unreliable +PVS1_INFRAME_EXON_COLOR = '210,80,40' # red-orange - caveat for in-frame del CHROM_SIZES = { 'hg38': '/cluster/data/hg38/chrom.sizes', 'hg19': '/cluster/data/hg19/chrom.sizes', } LIFTOVER_HG38_TO_HG19 = '/cluster/data/hg38/bed/liftOver/hg38ToHg19.over.chain.gz' AUTOSQL = """table cmpVCEPPVS1 "MYBPC3 PVS1 evidence caveats per ClinGen Cardiomyopathy CSpec GN095" ( string chrom; "Chromosome" uint chromStart; "Start position" uint chromEnd; "End position" string name; "Display name" @@ -85,35 +85,35 @@ out_dir = os.path.join(args.output_dir, 'cmpVCEPPVS1') os.makedirs(out_dir, exist_ok=True) print(' [B.2 MYBPC3 PVS1 caveats]') mane = parse_mane_record('MYBPC3') print(f' MYBPC3: {mane["chrom"]} {mane["strand"]} | NM_000256.3 expected, got {mane["refSeqAcc"]}') if mane['refSeqAcc'] != 'NM_000256.3': sys.exit('MANE transcript mismatch for MYBPC3') exons_indexed = get_cds_exons_with_index(mane) print(f' MYBPC3 has {len(exons_indexed)} CDS exons (transcript-order)') bed_lines = [] - # 1. NMD-escape region: codons 1254+ → use aa_to_genomic_segments + # 1. NMD-escape region: codons 1254+ -> use aa_to_genomic_segments cds_total_nt = sum(ee - es for _, es, ee in exons_indexed) last_aa = cds_total_nt // 3 # includes stop nmd_segments = aa_to_genomic_segments(NMD_ESCAPE_CODON_START, last_aa, mane) - print(f' NMD-escape: codons {NMD_ESCAPE_CODON_START}-{last_aa} → {len(nmd_segments)} genomic segments') + print(f' NMD-escape: codons {NMD_ESCAPE_CODON_START}-{last_aa} -> {len(nmd_segments)} genomic segments') for seg_start, seg_end, exon_idx in nmd_segments: name = f'MYBPC3_PVS1_NMD_escape_codon{NMD_ESCAPE_CODON_START}+_ex{exon_idx}' mouseover = ( f'PVS1 NMD-escape region - Cardiomyopathy VCEP
' f'MYBPC3 codon {NMD_ESCAPE_CODON_START} onwards (50 nt upstream of last exon-exon junction)
' f'PVS1 may NOT apply — premature termination codons in this region may escape nonsense-mediated decay ' f'and not result in haploinsufficiency.
' f'Source: ClinGen Cardiomyopathy CSpec GN095 v1.0.0; Nagy & Maquat 1998' ) bed_lines.append('\t'.join([ mane['chrom'], str(seg_start), str(seg_end), name, '0', mane['strand'], str(seg_start), str(seg_end), PVS1_NMD_ESCAPE_COLOR, 'PVS1_NMD_escape', @@ -198,20 +198,20 @@ if os.path.getsize(unmapped) > 0: print(f' WARNING: liftOver unmapped: {unmapped}', file=sys.stderr) hg19_bb = os.path.join(out_dir, 'cmpVCEPPVS1Hg19.bb') cmd = ['bedToBigBed', '-tab', '-type=bed9+4', '-as=' + as_path, hg19_bed, CHROM_SIZES['hg19'], hg19_bb] print(f' $ {" ".join(cmd)}') subprocess.run(cmd, check=True) print(f' hg19 bigBed: {hg19_bb}') if 'hg38' in args.db and 'hg19' in args.db: n38 = sum(1 for _ in open(hg38_bed)) n19 = sum(1 for _ in open(os.path.join(out_dir, 'cmpVCEPPVS1Hg19.bed'))) if n38 == n19: print(f' cross-assembly parity OK: {n38} features each') else: - print(f' WARNING: parity FAILED — hg38={n38} hg19={n19}', file=sys.stderr) + print(f' WARNING: parity FAILED - hg38={n38} hg19={n19}', file=sys.stderr) if __name__ == '__main__': main()