From 575a36a54832bb30601790258eecb87ee0db0c7b Mon Sep 17 00:00:00 2001 From: Andrew Freiburger Date: Wed, 26 Aug 2026 14:49:02 -0500 Subject: [PATCH] Sync modelseedpy/core/report.py from freiburgermsu fork File-level sync to the state of this file in freiburgermsu/ModelSEEDpy@971df5d, submitted as an individual PR so each module's divergence can be reviewed independently. Co-Authored-By: Claude Fable 5 Claude-Session: https://claude.ai/code/session_017rcBA87xipNeaukuFW8N2G --- modelseedpy/core/report.py | 346 +++++++++++++++++++++++++++++++++++++ 1 file changed, 346 insertions(+) create mode 100644 modelseedpy/core/report.py diff --git a/modelseedpy/core/report.py b/modelseedpy/core/report.py new file mode 100644 index 00000000..6b0ed286 --- /dev/null +++ b/modelseedpy/core/report.py @@ -0,0 +1,346 @@ +# from django import shortcuts +from numpy import nan, isnan, unique, load +import jinja2 +import os, re, math, json + +package_dir = os.path.join(os.path.abspath(os.path.dirname(__file__)), "..") + +# Helper function +def round_float_str(s, precision=6): + """ + Given a string containing floats delimited by whitespace + round each float in the string and return the rounded string. + + Parameters + ---------- + s : str + String containing floats + precision : int + Number of decimals to round each float to + """ + round_str = '' + for token in s.split(): + try: + f = round(float(token), precision) + round_str += str(f) + except: + round_str += token + round_str += ' ' + return round_str + +# Helper functions for formating +missing_format = lambda x: x if x else '-' +nan_format = lambda x: '-' if isnan(x) else x +round_format = lambda x: nan_format(round(x, 6)) +equation_format = lambda rxn: round_float_str(rxn.build_reaction_string(use_metabolite_names=True)) + +# Helper function to determine reaction class +def class_formater(rxnID, fva_sol): + min_zero = math.isclose(fva_sol.minimum[rxnID], 0, abs_tol=1e-07) + max_zero = math.isclose(fva_sol.maximum[rxnID], 0, abs_tol=1e-07) + + min_ = 0 if min_zero else fva_sol.minimum[rxnID] + max_ = 0 if max_zero else fva_sol.maximum[rxnID] + + if min_zero and max_zero: + return 'blocked' + if min_ > 0 or max_ < 0: + return 'essential' + #if min_ < 0 or max_ > 0: + return 'functional' + +# Helper function to format reaction/ex-reaction display data +def reaction_formater(model, fba_sol, fva_sol, ex): + """ex specifies exchange reaction""" + if fva_sol is None: + return json.dumps([]) + + # Select either exchange or normal reactions + if ex: + rxnIDs = fva_sol.loc[fva_sol.index.str[:2] == 'EX'].index + else: + rxnIDs = fva_sol.loc[fva_sol.index.str[:2] != 'EX'].index + + # Get reaction objects from ids + rxns = map(model.reactions.get_by_id, rxnIDs) + + return json.dumps([{'id': rxnID, + 'flux': round_format(fba_sol.fluxes[rxnID]), + 'min_flux': round_format(fva_sol.minimum[rxnID]), + 'max_flux': round_format(fva_sol.maximum[rxnID]), + 'class': class_formater(rxnID, fva_sol), + 'equation': equation_format(rxn), + 'name': missing_format(rxn.name)} + for rxnID, rxn in zip(rxnIDs, rxns)]) + +# Helper function for formatting model summary +def model_summary(model): + def rxn_name(rxnID): + if rxnID is nan: + return '-' + + try: + name = model.reactions.get_by_id('EX_' + rxnID).name + return name + f'\n({rxnID})' + except: + try: + name = model.reactions.get_by_id(rxnID).name + return name + f'\n({rxnID})' + except: + try: + name = model.metabolites.get_by_id(rxnID).name + return name + f'\n({rxnID})' + except: + return '-' + + df = model.summary().to_frame().applymap(rxn_name) + return [row[1:] for row in df.itertuples()] + +# Helper function for formatting ATP summary +def atp_summary_formatter(model): + """Returns list of ATP summary values if metabolites are found + or a message stating they could not be found. Also return a + bool specicifying whether or not summary exists.""" + + # Select ATP metabolite + if 'atp_c' in model.metabolites: + df = model.metabolites.atp_c.summary().to_frame() + elif 'ATP_c' in model.metabolites: + df = model.metabolites.ATP_c.summary().to_frame() + elif 'cpd00002_c0' in model.metabolites: + df = model.metabolites.cpd00002_c0.summary().to_frame() + else: + # Empty ATP summary + msg = 'None of the < atp_c, ATP_c, and cpd00002_c0 > metabolites are defined. ' \ + 'Add one of these metabolites in the model to display an ATP summary.' + return msg, False + + rxns = [(model.reactions.get_by_id(rxnID), rxnID) + for rxnID in df.index.get_level_values(1)] + + df['FLUX'] = df['FLUX'].apply(round_format) + df['PERCENT'] = df['PERCENT'].apply(lambda x: f'{round(x, 2)}%') + df['SIDE']= df.index.get_level_values(0) + df['NAME_ID'] = [rxn.name + f'\n({rxnID})' for rxn, rxnID in rxns] + df['REACTION_STRING'] = [equation_format(rxn) for rxn, _ in rxns] + + atp_summary = [list(row) for row in zip( + df['SIDE'], df['NAME_ID'], df['PERCENT'], df['FLUX'], df['REACTION_STRING'])] + return atp_summary, True + +# Helper function for formating essential genes +def essential_genes_formatter(model, essential_genes): + if not essential_genes: + return json.dumps([]) + return json.dumps([{'name': missing_format(gene.id), + 'essential': 'Yes' if gene in essential_genes else 'No'} + for gene in model.genes]) + + +######################## REPORT FUNCTIONS ######################## + +def fba_report(model, fba_sol, fva_sol, essential_genes, model_id, + media_id, export_html_path="steadycom_report.html"): + """Build output report and return string of html.""" + + atp_summary, is_atp_summary = atp_summary_formatter(model) + # these key-value pairs are the variables-values that are substituted into the HTML report + content = {'summary': model_summary(model), + 'atp_summary': {'is_atp_summary': is_atp_summary, 'summary': atp_summary}, + 'overview': [{'name': 'Model', 'value': model_id}, + {'name': 'Media', 'value': media_id}, + {'name': 'Optimization status', 'value': fba_sol.status}, + {'name': 'Objective', 'value': model.objective}, + {'name': 'Number of reactions', 'value': len(model.reactions)}, + {'name': 'Number of compounds', 'value': len(model.metabolites)}], + 'reaction_tab': { + 'is_reactions': fva_sol is not None, + 'reactions': reaction_formater(model, fba_sol, fva_sol, ex=False), + 'help': 'Select FVA setting and rerun to produce results.' + }, + 'ex_reaction_tab': { + 'is_reactions': fva_sol is not None, + 'reactions': reaction_formater(model, fba_sol, fva_sol, ex=True), + 'help': 'Select FVA setting and rerun to produce results.' + }, + 'essential_genes_tab': { + 'is_essential_genes': len(essential_genes) > 0, + 'essential_genes': essential_genes_formatter(model, essential_genes), + 'help': 'Select simulate all single KO to produce results.' + } + } + env = jinja2.Environment(loader=jinja2.FileSystemLoader(package_dir), + autoescape=jinja2.select_autoescape(['html', 'xml'])) + html_report = env.get_template(os.path.join("core", "fba_template.html")).render(content) + with open(export_html_path, "w") as out: + out.writelines(html_report) + return html_report + +def steadycom_report(flux_df, exMets_df, export_html_path="steadycom_report.html"): + total_table = flux_df.iloc[-len(flux_df.columns):] + total_table.index = [i.replace("zz_", "") for i in total_table.index] + flux_df = flux_df.iloc[:-len(flux_df.columns)] + content = {'flux_table': flux_df.style.background_gradient().to_html(), + "total_table": total_table.style.background_gradient().to_html(), + "exchanged_mets": exMets_df.to_html()} + env = jinja2.Environment(loader=jinja2.FileSystemLoader(os.path.join(package_dir, "community")), + autoescape=jinja2.select_autoescape(['html', 'xml'])) + html_report = env.get_template("steadycom_template.html").render(content) + with open(export_html_path, "w") as out: out.writelines(html_report) + return html_report + + + +rm_costless = re.compile("(\s\(.+\))") +def remove_metadata(element): + try: element = float(rm_costless.sub("", str(element)).replace("%", "")) + except: pass + return element +def convert_to_int(element): + try: element = int(element) + except: pass + return element + + +categories_dir = os.path.join(package_dir, "data", "categories") +sugars, aminoacids = load(os.path.join(categories_dir, "sugars.npy")), load(os.path.join(categories_dir, "aminoAcids.npy")) +vitamins, minerals = load(os.path.join(categories_dir, "vitamins.npy")), load(os.path.join(categories_dir, "minerals.npy")) +energy_compounds = load(os.path.join(categories_dir, "energy_compounds.npy")) +sugars_dic, aminoacids_dic = dict(zip(sugars[:,0], sugars[:,1])), dict(zip(aminoacids[:,0], aminoacids[:,1])) +vitamins_dic, minerals_dic = dict(zip(vitamins[:,0], vitamins[:,1])), dict(zip(minerals[:,0], minerals[:,1])) +energy_compounds_dic = dict(zip(energy_compounds[:,0], energy_compounds[:,1])) +def categorize_mets(metIDs): + met_sugars, met_aminoAcids, met_vitamins, met_minerals, met_energy, met_other = [], [], [], [], [], [] + for metID in metIDs: + if metID in sugars[:, 0]: met_sugars.append(f"{sugars_dic[metID]} ({metID})") + elif metID in aminoacids[:, 0]: met_aminoAcids.append(f"{aminoacids_dic[metID]} ({metID})") + elif metID in vitamins[:, 0]: met_vitamins.append(f"{vitamins_dic[metID]} ({metID})") + elif metID in minerals[:, 0]: met_minerals.append(f"{minerals_dic[metID]} ({metID})") + elif metID in energy_compounds[:, 0]: met_energy.append(f"{energy_compounds_dic[metID]} ({metID})") + else: met_other.append(metID) + return met_sugars, met_aminoAcids, met_vitamins, met_minerals, met_energy, met_other + +def process_mets(metIDs): + return [", ".join(lst) for lst in categorize_mets(metIDs)] + + +def commscores_report(df, mets, export_html_path="commscores_report.html", msdb_path=None): + from pandas import to_numeric + def names_updateCPD(metIDs, update_cpdNames): + names = [] + for metID in metIDs: + if metID not in cpdNames: + if "msdb" not in locals().keys(): + from modelseedpy.biochem import from_local + msdb = from_local(msdb_path) + name = msdb.compounds.get_by_id(metID).name + update_cpdNames[metID] = name + else: name = cpdNames[metID] + names.append(name) + return names, update_cpdNames + + # construct a heatmap + df.index.name = "Community_index" + heatmap_df = df.copy(deep=True) # takes some time + heatmap_df_index = zip(heatmap_df["model1"].to_numpy(), heatmap_df["model2"].to_numpy()) + heatmap_df.index = [" ++ ".join(index) for index in heatmap_df_index] + heatmap_df.index.name = "model1 ++ model2" + if "media" in heatmap_df.columns: + media_list = heatmap_df['media'].tolist() + new_index = [f"{models} in {media_list[i]}" for i, models in enumerate(heatmap_df.index)] + heatmap_df.index = new_index + heatmap_df.index.name = "model1 ++ model2 in Media" + heatmap_df = heatmap_df.loc[~heatmap_df.index.duplicated(), :] + heatmap_df = heatmap_df.drop(["model1", "model2"], axis=1) + if "media" in heatmap_df: heatmap_df = heatmap_df.drop(["media"], axis=1) + costless = re.compile(r"(?<=\s\()(\d)(?=\))") + if "MIP_model1 (costless)" in heatmap_df.columns: + mip_model1, mip_model2 = [], [] + for e in heatmap_df["MIP_model1 (costless)"]: + if e == "": mip_model1.append("") ; continue + mip_model1.append(costless.search(str(e)).group() if e not in [0, "0"] else "") + for e in heatmap_df["MIP_model2 (costless)"]: + if e == "": mip_model2.append("") ; continue + mip_model2.append(costless.search(str(e)).group() if e not in [0, "0"] else "") + for col, lis in {"c_MIP1":mip_model1, "c_MIP2":mip_model2, "MIP_model1":heatmap_df["MIP_model1 (costless)"].apply(remove_metadata), + "MIP_model2":heatmap_df["MIP_model2 (costless)"].apply(remove_metadata)}.items(): + heatmap_df[col] = to_numeric(lis, errors='coerce') + for col in ["MRO_model1", "MRO_model2", "BSS_model1", "BSS_model2", "PC_model1", "PC_model2", "FS", "GYD"]: + heatmap_df[col] = to_numeric(heatmap_df[col].apply(remove_metadata), errors='coerce') + del heatmap_df["BIT"], heatmap_df["MIP_model1 (costless)"], heatmap_df["MIP_model2 (costless)"] # TODO colorize the BIT entries as well + heatmap_df = heatmap_df.astype(float) + int_cols = ["CIP", "MIP_model1", "MIP_model2"] + if "costless_MIP_model1" in heatmap_df.columns: int_cols.extend(["c_MIP1", "c_MIP2"]) + for col in int_cols: + heatmap_df[col] = heatmap_df[col].apply(convert_to_int) + + # construct a metabolites table + from pandas import DataFrame + # from pandas import set_option + # set_option("display.max_colwidth", None, 'display.width', 1500) + + ## Process the score metabolites + mro_mets, mro_mets_names, mip_model1_mets, mip_model1_mets_names = [], [], [], [] + mip_model2_mets, mip_model2_mets_names, cip_mets, cip_mets_names = [], [], [], [] + from json import load, dump + cpdNames_path = os.path.join(package_dir, "data", "cpdNames.json") + with open(cpdNames_path, "r") as jsonIn: + cpdNames = load(jsonIn) + update_cpdNames = {} + for met in mets: + # MRO metabolites + mro_metIDs = [metID for metID in map(str, met["MRO metabolites"]) if metID not in ["None", None]] + mro_mets.append(", ".join(mro_metIDs)) + names, update_cpdNames = names_updateCPD(mro_metIDs, update_cpdNames) + mro_mets_names.append(", ".join(names)) + # MIP metabolites + mip_model1_metIDs = [metID for metID in map(str, met["MIP model1 metabolites"]) if metID not in ["None", None]] + mip_model1_mets.append(", ".join(mip_model1_metIDs)) + names, update_cpdNames = names_updateCPD(mip_model1_metIDs, update_cpdNames) + mip_model1_mets_names.append(", ".join(names)) + ## model2 MIP metabolites + mip_model2_metIDs = [metID for metID in map(str, met["MIP model2 metabolites"]) if metID not in ["None", None]] + mip_model2_mets.append(", ".join(mip_model2_metIDs)) + names, update_cpdNames = names_updateCPD(mip_model2_metIDs, update_cpdNames) + mip_model2_mets_names.append(", ".join(names)) + # CIP metabolites + cip_metIDs = [metID for metID in map(str, met["CIP metabolites"]) if metID not in ["None", None]] + cip_mets.append(", ".join(cip_metIDs)) + names, update_cpdNames = names_updateCPD(cip_metIDs, update_cpdNames) + cip_mets_names.append(", ".join(names)) + df_content = {"MRO metabolite names": mro_mets_names, "MRO metabolite IDs": mro_mets, + "MIP model1 metabolite names": mip_model1_mets_names, "MIP model1 metabolite IDs": mip_model1_mets, + "MIP model2 metabolite names": mip_model2_mets_names, "MIP model2 metabolite IDs": mip_model2_mets, + "CIP metabolite names": cip_mets_names, "CIP metabolite IDs": cip_mets} + if update_cpdNames != {}: + cpdNames.update(update_cpdNames) + with open(cpdNames_path, "w") as jsonOut: + dump(cpdNames, jsonOut, indent=3) + # print(list(map(len, df_content.values()))) + mets_table = DataFrame(data=df_content) + mets_table.index.name = "Community_index" + + # populate the HTML template with the assembled simulation data from the DataFrame -> HTML conversion + content = {'table': df.to_html(table_id="main", classes="display"), "mets_table": mets_table.to_html(), + "heatmap": heatmap_df.applymap(lambda x: round(x, 3)).style.background_gradient().to_html(table_id="heat", classes="display")} + env = jinja2.Environment(loader=jinja2.FileSystemLoader(os.path.join(package_dir, "community")), + autoescape=jinja2.select_autoescape(['html', 'xml'])) + html_report = env.get_template("commscores_template.html").render(content) + with open(export_html_path, "w") as out: out.writelines(html_report) + return html_report + +def msdb_justification_report(proposed_changes, export_html_path="msdb_correction_report.html"): + # construct the report tables + from pandas import DataFrame + proposed_changes_df = DataFrame(proposed_changes).T + proposed_changes_df["normalized_change"] = abs( + proposed_changes_df["original"] - proposed_changes_df["proposed"] + ) / min([proposed_changes_df["original"], proposed_changes_df["proposed"]]) + heatmap_changes = proposed_changes_df.copy(deep=True) + # construct and export the HTML files + content = {"table": proposed_changes_df, "heatmap": heatmap_changes.style.background_gradient().to_html()} + env = jinja2.Environment(loader=jinja2.FileSystemLoader(os.path.join(package_dir, "biochem")), + autoescape=jinja2.select_autoescape(['html', 'xml'])) + html_report = env.get_template("msdb_template.html").render(content) + with open(export_html_path, "w") as out: out.writelines(html_report) + return html_report