ecc2331000637cbc2523653642d926935b5b83fc chmalee Sat Jun 13 00:42:03 2026 -0700 gnomAD v4.1.1 bigBed variant track for hg38, refs #37351 Co-Authored-By: Claude Opus 4.8 diff --git src/hg/makeDb/gnomad/gnomadVcfBedToBigBed src/hg/makeDb/gnomad/gnomadVcfBedToBigBed index 30026e5231d..40d0beb2374 100755 --- src/hg/makeDb/gnomad/gnomadVcfBedToBigBed +++ src/hg/makeDb/gnomad/gnomadVcfBedToBigBed @@ -12,42 +12,64 @@ Format: \ Allele|Consequence|IMPACT|SYMBOL|Gene|Feature_type|Feature|BIOTYPE|EXON|INTRON|HGVSc\ |HGVSp|cDNA _position|CDS_position|Protein_position|Amino_acids|Codons\ |Existing_variation|ALLELE_NUM|DISTANCE|STRAND|FLAGS|VARIANT_CLASS|MINIMISED\ |SYMBOL_SOURCE|HGNC_ID|CANONICAL|TSL|APPRIS|CCDS|ENSP|SWISSPROT|TREMBL|UNIPARC\ |GENE_PHENO|SIFT|PolyPhen|DOMAINS|HGVS_OFFSET|GMAF|AFR_MAF|AMR_MAF|EAS_MAF|EUR_MAF\ |SAS_MAF|AA_MAF|EA_MAF|ExAC_MAF|ExAC_Adj_MAF|ExAC_AFR_MAF|ExAC_AMR_MAF|ExAC_EAS_MAF\ |ExAC_FIN_MAF|ExAC_NFE_MAF|ExAC_OTH_MAF|ExAC_SAS_MAF|CLIN_SIG|SOMATIC|PHENO|PUBMED\ |MOTIF_NAME|MOTIF_POS|HIGH_INF_POS|MOTIF_SCORE_CHANGE|LoF|LoF_filter|LoF_flags|LoF_info"> """ import sys, argparse, hashlib,json from collections import defaultdict,namedtuple,OrderedDict # which version of gnomAD for parsing VEP string -versions = ["v2.1.1", "v3.1", "v3.1_chrM", "v3.1.1"] +versions = ["v2.1.1", "v3.1", "v3.1_chrM", "v3.1.1", "v4.1_genomes", "v4.1_exomes"] # the number of fields in the VEP string (depends on version): # how to count: # bcftools view -h in.vcf.gz | grep "^##INFO= 13 else ref displayName += "-" displayName += alt[:10]+"..." if len(alt) > 13 else alt name = hashlib.md5(str.encode(chrom+chromStart+chromEnd+ref+alt)).hexdigest() - if version == "v3.1.1" or version == "v3.1_chrM": + if version == "v3.1.1" or version == "v3.1_chrM" or version.startswith("v4.1"): outfh.write("\t".join([chrom, chromStart, chromEnd, name, score, strand, thickStart, thickEnd, color, ref, alt, filterTag, ac, an, af, faf95, nhomalt, rsId, ", ".join(genes.keys()), annot, ",".join(consList), str(unshiftedStart), displayName])) + if version.startswith("v4.1"): + # add the grp max information for the mouseovers + grpmax = popAbbr[fixedLine[18]] if fixedLine[18] in popAbbr else "N/A" + AC_grpmax = fixedLine[19] if fixedLine[19] else "N/A" + AN_grpmax = fixedLine[20] if fixedLine[20] else "N/A" + AF_grpmax = fixedLine[21] if fixedLine[21] else "N/A" + outfh.write("\t" + "\t".join([grpmax, AC_grpmax, AN_grpmax, AF_grpmax])) + # add the hemizygosity + nHomaltX = "N/A" + nHemi = "N/A" + if version == "v4.1_genomes": + nHomaltX = fixedLine[65] if fixedLine[65] else "N/A" + if chrom == "chrX" or chrom == "chrY": + nHemi = fixedLine[66] if fixedLine[66] else "N/A" + else: + nHomaltX = fixedLine[61] if fixedLine[61] else "N/A" + if chrom == "chrX" or chrom == "chrY": + nHemi = fixedLine[62] if fixedLine[62] else "N/A" + outfh.write("\t" + nHomaltX + "\t" + nHemi) popTable, otherFields = convertExtraFields(maybeExtraFields + [ac, an, af, nhomalt], version) if version == "v3.1_chrM": # otherFields is the haplotype table plus empties for the v3.1.1 fields: otherFields = ["" for i in range(5)] + [json.dumps(otherFields)] if extraFh: extraFh.write("\t".join([name, json.dumps(genes), json.dumps(popTable)] + otherFields)) extraFh.write("\n") else: outfh.write("\t" + "\t".join([name] + [json.dumps(genes), json.dumps(popTable)] + otherFields)) outfh.write("\n") elif version == "v3.1": hgvscList = [", ".join(list(genes[gene]["hgvsc"]))] hgvspList = [", ".join(list(genes[gene]["hgvsp"]))] pLoFList = [", ".join(list(genes[gene]["pLoF"]))] pLoFFlags = [", ".join(list(genes[gene]["Flag"]))] outfh.write("\t".join([chrom, chromStart, chromEnd, name, score, strand, thickStart, thickEnd, color, ref, alt, filterTag, ac, an, af, faf95, nhomalt, rsId, gene, annot, *consList, str(unshiftedStart), displayName])) if extraFh: extraFh.write("\t".join([name] + hgvscList + hgvspList + maybeExtraFields)) extraFh.write("\n") else: outfh.write("\t" + "\t".join([name] + hgvscList + hgvspList + maybeExtraFields)) outfh.write("\n") else: - pLoFCuration = getLofCuration(lofDict, version, chrom, + pLoFCuration = getLofCuration(lofDict, version, bed8Fields[0], str(unshiftedStart), ref, alt) - hgvscList = [x for x in genes[gene]["hgvsc"] for gene in genes] - hgvspList = [x for x in genes[gene]["hgvsp"] for gene in genes] - pLoFList = [x for x in genes[gene]["pLoF"] for gene in genes] - pLoFFlags = [x for x in genes[gene]["Flag"] for gene in genes] outfh.write("\t".join([chrom, chromStart, chromEnd, name, score, strand, thickStart, thickEnd, color, ref, alt, filterTag, ac, an, af, faf95, nhomalt, - rsId, ", ".join(genes.keys()), annot, ", ".join(consList)])) - outfh.write("\t" + str(unshiftedStart)) - outfh.write("\t" + displayName) - outfh.write("\t" + ", ".join(pLoFList)) - outfh.write("\t" + ", ".join(pLoFFlags)) - outfh.write("\t" + "\t".join(pLoFCuration)) + rsId, gene, annot, *consList, displayName])) + outfh.write("\t" + "\t".join(pLoFCuration, [str(unshiftedStart)])) if extraFh: - extraFh.write(", ".join(hgvscList) + "\t" + ", ".join(hgvspList) + "\t" + "\t".join(maybeExtraFields)) + extraFh.write("\t".join(name + hgvscList + hgvspList + maybeExtraFields)) extraFh.write("\n") else: - outfh.write("\t" + ", ".join(hgvscList)) - outfh.write("\t" + ", ".join(hgvspList)) - outfh.write("\t" + "\t".join(maybeExtraFields)) + outfh.write("\t" + "\t".join([name] + hgvscList + hgvspList + maybeExtraFields)) outfh.write("\n") def parseLofFile(fpath): """Make a struct of the different loss of function flags for a curated variant.""" gotHeader = False lofHeader = [] ret = {} with open(fpath) as fh: for line in fh: if not gotHeader: lofHeader = line.strip().split("\t") gotHeader = True else: lofDetails = line.strip().split("\t") ret[lofDetails[0]] = {lofHeader[x]: lofDetails[x] for x in range(len(lofHeader))} @@ -535,18 +610,21 @@ infh = sys.stdin else: infh = open(args.infile) if args.outfile == "stdout": outfh = sys.stdout else: outfh = open(args.outfile, "w") if args.extra_output_file: if args.extra_output_file == "stdout": extraFh = sys.stdout else: extraFh = open(args.extra_output_file, "w") gnomadVcfBedToBigBed(infh, outfh, extraFh, args.version, lofDict) infh.close() outfh.close() + # close the extra-output file too, else its final buffered line can be lost + if extraFh is not None and extraFh is not sys.stdout: + extraFh.close() if __name__ == "__main__": main()