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/cmpVCEPAnnotate.py src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPAnnotate.py index ff57aa546b2..b5c5f5a5a5c 100644 --- src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPAnnotate.py +++ src/hg/makeDb/scripts/cardiomyopathyVCEP/cmpVCEPAnnotate.py @@ -1,18 +1,18 @@ #!/usr/bin/env python3 """ -Phase 1 — hgVai annotation layer for the Cardiomyopathy VCEP hub (RM37446). +Phase 1 - hgVai annotation layer for the Cardiomyopathy VCEP hub (RM37446). Annotates the gnomAD-observed variant universe (the SAME universe used by the AF + Provisional tracks) with protein consequence via the internal UCSC tool hgVai (vai.pl), on the single MANE Select transcript per gene (ncbiRefSeqSelect). Reuses cmpVCEPAFfrequencies.fetch_gene_variants so the universe is identical. Output (consumed by Phase 2 cmpVCEPProvisionalClass): cmpVCEPAnnotate/cmpVCEPAnnotations.hg38.tsv one row per (chrom,pos,ref,alt) cmpVCEPAnnotate/phase1_unmapped.tsv universe variants with no VEP line cmpVCEPAnnotate/<gene>.vcf / <gene>.vep per-gene intermediates Annotation is computed on hg38 only; Phase 2 codes carry to hg19 via the existing liftOver of the universe (annotation is transcript-based, assembly-independent). Usage: @@ -68,31 +68,31 @@ def parse_extra(extra): """Pull HGVSP and EXON out of the VEP Extra column (key=value;key=value).""" hgvsp, exon_n, exon_total = '', '', '' for kv in extra.split(';'): if kv.startswith('HGVSP='): hgvsp = kv[len('HGVSP='):] elif kv.startswith('EXON='): m = re.match(r'(\d+)/(\d+)', kv[len('EXON='):]) if m: exon_n, exon_total = m.group(1), m.group(2) return hgvsp, exon_n, exon_total def parse_vep(vep_path, gene): - """Parse one gene's VEP output → dict key (chrom,pos,ref,alt) → annotation. + """Parse one gene's VEP output -> dict key (chrom,pos,ref,alt) -> annotation. Aggregates SO terms across lines; keeps the coding fields from the line(s) that carry them.""" ann = {} for line in open(vep_path): if line.startswith('#') or line.startswith('Uploaded'): continue f = line.rstrip('\n').split('\t') if len(f) < 14: continue parts = f[0].split(':') if len(parts) != 4: continue chrom, pos, ref, alt = parts[0], int(parts[1]), parts[2], parts[3] key = (chrom, pos, ref, alt) so = f[6] @@ -153,37 +153,37 @@ all_ann.update(ann) print(f' {gene} ({mane["chrom"]} {mane["strand"]}): ' f'{len(universe_keys)} universe, {len(ann)} annotated, ' f'{sum(1 for k in universe_keys if k not in ann_keys)} unmapped') tsv = os.path.join(work, 'cmpVCEPAnnotations.hg38.tsv') with open(tsv, 'w') as fh: fh.write('#chrom\tpos\tref\talt\tgene\tsoTerms\tproteinPos\taaRef\taaAlt\t' 'codonChange\texonNum\texonTotal\thgvsp\tcdnaPos\n') for (chrom, pos, ref, alt), r in sorted(all_ann.items(), key=lambda x: (x[0][0], x[0][1])): fh.write('\t'.join([ chrom, str(pos), ref, alt, r['gene'], ','.join(sorted(r['so'])), r['proteinPos'], r['aaRef'], r['aaAlt'], r['codon'], r['exonNum'], r['exonTotal'], r['hgvsp'], r['cdnaPos'], ]) + '\n') - print(f' wrote {len(all_ann)} annotations → {tsv}') + print(f' wrote {len(all_ann)} annotations -> {tsv}') unmapped_path = os.path.join(work, 'phase1_unmapped.tsv') with open(unmapped_path, 'w') as fh: fh.write('#gene\tchrom\tpos\tref\talt\n') for row in unmapped: fh.write('\t'.join(str(x) for x in row) + '\n') - print(f' {len(unmapped)} unmapped universe variants logged → {unmapped_path}') + print(f' {len(unmapped)} unmapped universe variants logged -> {unmapped_path}') # SO-term summary for sanity from collections import Counter so_counter = Counter() for r in all_ann.values(): for s in r['so']: so_counter[s] += 1 print(' SO-term distribution:') for s, c in so_counter.most_common(): print(f' {c:6d} {s}') if __name__ == '__main__': main()