aa3bab94fb50013ebdf152f761359925fc88b674 mspeir Tue Aug 4 14:53:24 2026 -0700 G2P otto: code-review fixes (logging, mixedCase names, itemRgb) - Log a per-value count when a confidence value is unrecognized and colored black. - Log a per-assembly count of G2P records with no HGNC coordinate match. - Unify function names to mixedCase (confidenceToColor, loadG2p, loadCoordinates, joinAndWrite). - Rename itemRGB to itemRgb in g2p.as for consistency with the rest of the tree. refs #36736 Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com> diff --git src/hg/utils/otto/g2p/doG2p.py src/hg/utils/otto/g2p/doG2p.py index a45cc683ec4..f23a5a82395 100755 --- src/hg/utils/otto/g2p/doG2p.py +++ src/hg/utils/otto/g2p/doG2p.py @@ -76,125 +76,148 @@ def validateColumns(csvFile): """Abort if any required column is missing from the CSV header.""" with open(EXPECTED_COLUMNS_FILE) as f: required = [line.strip() for line in f if line.strip()] with open(csvFile, newline="", encoding="utf-8") as f: header = next(csv.reader(f)) header = [h.strip() for h in header] missing = [c for c in required if c not in header] if missing: sys.exit("ERROR: G2P CSV is missing expected column(s): %s\n" "The source format may have changed; check %s" % (missing, csvFile)) -def confidence_to_color(confidence): - """Map a confidence string to an RGB color string for UCSC BED itemRgb.""" - color_map = { +# Confidence value -> itemRgb color. Unrecognized values fall back to DEFAULT_COLOR +# (black) and are counted/logged by joinAndWrite so a source change is visible. +CONFIDENCE_COLORS = { "definitive": "39,103,73", # dark green "strong": "56,161,105", # green "moderate": "104,211,145", # light green "limited": "252,129,129", # pink "disputed": "229,62,62", # red "refuted": "155,44,44", # dark red } - return color_map.get(confidence.lower(), "0,0,0") # default black +DEFAULT_COLOR = "0,0,0" # black, for unrecognized confidence values -def load_g2p(file_path): +def confidenceToColor(confidence): + """Return the itemRgb color for a confidence string, or None if unrecognized.""" + return CONFIDENCE_COLORS.get(confidence.lower().strip()) + + +def loadG2p(file_path): """Load G2P CSV into a dict keyed by HGNC ID (each value is a list of rows).""" g2p_map = {} numOfRows = 0 with open(file_path, newline="", encoding="utf-8") as csvfile: reader = csv.DictReader(csvfile) for row in reader: numOfRows += 1 hgnc_id = row["hgnc id"].strip() g2p_map.setdefault(hgnc_id, []).append(row) print("Number of rows in file: %s" % numOfRows) return g2p_map -def load_coordinates(db, hgnc_ids): +def loadCoordinates(db, hgnc_ids): """Build a dict of gene coordinates for the given HGNC IDs from the HGNC bigBed. One bigBedToBed pass over the whole track (~49k rows) instead of one bigBedNamedItems subprocess per HGNC ID. The bigBed name field is "HGNC:<id>"; the G2P CSV stores the bare numeric id, so we key on that. """ wanted = set(hgnc_ids) coord_map = {} hgncBB = "/gbdb/%s/hgnc/hgnc.bb" % db for line in bash("bigBedToBed %s stdout" % hgncBB).split("\n"): if not line.strip(): continue fields = line.split("\t")[:8] name = fields[3] # e.g. "HGNC:36036" id = name.split("HGNC:")[-1] if id in wanted: coord_map.setdefault(id, []).append(fields) return coord_map -def join_and_write(g2p_data, coords, output_file): - """Join G2P records and HGNC coordinates into BED 9+20 and write to output_file.""" +def joinAndWrite(g2p_data, coords, output_file): + """Join G2P records and HGNC coordinates into BED 9+20 and write to output_file. + + Returns a stats dict: + "unmatched" -> count of G2P records whose HGNC ID had no coordinate + match in this assembly's HGNC track (they are skipped). + "unknownConfidence" -> {confidence value: count} for values not in + CONFIDENCE_COLORS (colored black). + """ + unmatched = 0 + unknownConfidence = {} with open(output_file, "w", newline="", encoding="utf-8") as out: writer = csv.writer(out, delimiter="\t") for hgnc_id, rows in g2p_data.items(): + matches = coords.get(hgnc_id, []) + if not matches: + unmatched += len(rows) + continue for row in rows: - for coord in coords.get(hgnc_id, []): + for coord in matches: # BED 9 fields chrom = coord[0] chromStart = coord[1] chromEnd = coord[2] name = row["gene symbol"] score = coord[4] strand = coord[5] thickStart = coord[6] thickEnd = coord[7] - rgb = confidence_to_color(row["confidence"]) + rgb = confidenceToColor(row["confidence"]) + if rgb is None: + unknownConfidence[row["confidence"]] = \ + unknownConfidence.get(row["confidence"], 0) + 1 + rgb = DEFAULT_COLOR # G2P 20 fields g2p_id = row["g2p id"] gene_mim = row["gene mim"] hgnc_id_val = row["hgnc id"] prev_symbols = row["previous gene symbols"].replace(";", ",") disease_name = row["disease name"] disease_mim = row["disease mim"] disease_MONDO = row["disease MONDO"] allelic_req = row["allelic requirement"] cross_mod = row["cross cutting modifier"] confidence = row["confidence"] var_conseq = row["variant consequence"] var_types = row["variant types"] mol_mech = row["molecular mechanism"] mol_mech_cat = row["molecular mechanism categorisation"] mol_mech_ev = row["molecular mechanism evidence"] phenotypes = row["phenotypes"].replace(";", ",") publications = row["publications"].replace(";", ",") panel = row["panel"] comments = row["comments"] date_review = row["date of last review"] writer.writerow([ chrom, chromStart, chromEnd, name, score, strand, thickStart, thickEnd, rgb, g2p_id, gene_mim, hgnc_id_val, prev_symbols, disease_name, disease_mim, disease_MONDO, allelic_req, cross_mod, confidence, var_conseq, var_types, mol_mech, mol_mech_cat, mol_mech_ev, phenotypes, publications, panel, comments, date_review, ]) + return {"unmatched": unmatched, "unknownConfidence": unknownConfidence} def itemCount(bb): line = bash('bigBedInfo %s | grep "itemCount"' % bb) return int(line.rstrip().split("itemCount:")[1].replace(",", "").strip()) def checkItemCount(db, newBb): """Abort if the item count moved more than COUNT_TOLERANCE vs the live track.""" liveBb = GBDB_BB % db if not Path(liveBb).exists(): print("%s: no live bigBed yet, skipping item-count check" % db) return old = itemCount(liveBb) new = itemCount(newBb) @@ -217,45 +240,51 @@ print("Installed %s -> %s" % (liveBb, newBb)) def main(): if not updateNeeded(): # Silent no-op: nothing new from G2P this run. return validateColumns(NEW_CSV) date = str(datetime.now()).split(" ")[0] buildDir = "%s/%s" % (WORKDIR, date) bash("mkdir -p %s" % buildDir) bash("cp %s %s/AllG2P.csv" % (NEW_CSV, buildDir)) - g2p_data = load_g2p(NEW_CSV) + g2p_data = loadG2p(NEW_CSV) hgnc_ids = list(g2p_data.keys()) print("Number of HGNC IDs found: %s" % len(hgnc_ids)) - coordsByDb = {db: load_coordinates(db, hgnc_ids) for db in DBS} + coordsByDb = {db: loadCoordinates(db, hgnc_ids) for db in DBS} for db in DBS: print("Loaded %s %s HGNC IDs" % (len(coordsByDb[db]), db)) builtBb = {} for db in DBS: bedFile = "%s/%s_g2p_all.bed" % (buildDir, db) bbFile = "%s/%s_g2p.bb" % (buildDir, db) twoBit = "/gbdb/%s/%s.2bit" % (db, db) - join_and_write(g2p_data, coordsByDb[db], bedFile) + stats = joinAndWrite(g2p_data, coordsByDb[db], bedFile) print("Wrote %s" % bedFile) + if stats["unmatched"]: + print("%s: %d G2P record(s) had no HGNC coordinate match and were skipped" + % (db, stats["unmatched"])) + for conf, n in sorted(stats["unknownConfidence"].items()): + print("%s: unrecognized confidence value %r on %d record(s); colored black" + % (db, conf, n)) bash("bedToBigBed -type=bed9+20 -tab -sort " "-as=%s -sizesIs2Bit -extraIndex=name,g2p_id,gene_mim,hgnc_id %s %s %s" % (AS_FILE, bedFile, twoBit, bbFile)) print("Built %s" % bbFile) builtBb[db] = bbFile # Safety check before swapping anything live. for db in DBS: checkItemCount(db, builtBb[db]) for db in DBS: install(db, builtBb[db]) bash("mv %s %s" % (NEW_CSV, PREV_CSV)) print("G2P updated %s" % date)