From bc4dc0c40695e7ef6319100681246454c8b80209 Mon Sep 17 00:00:00 2001 From: Titania Sugiarto Date: Fri, 6 Dec 2024 15:56:24 +0100 Subject: [PATCH] run_simulation_for_pam_mcpam function is now updated. can now used to run full scale or core models --- Examples/mcpam_generator.py | 23 +++++------- Scripts/mcpam_generation_uniprot_id.py | 2 +- Scripts/mcpam_simulations_analysis.py | 52 ++++++++------------------ src/PAModelpy/MembraneSector.py | 2 +- 4 files changed, 27 insertions(+), 52 deletions(-) diff --git a/Examples/mcpam_generator.py b/Examples/mcpam_generator.py index 0672ff6..926f317 100644 --- a/Examples/mcpam_generator.py +++ b/Examples/mcpam_generator.py @@ -1,27 +1,24 @@ from Scripts.mcpam_simulations_analysis import (run_pam_mcpam_core_with_optimized_kcats, - run_simulations_pam_mcpam_core, + run_simulations_pam_mcpam, set_up_ecolicore_pam, set_up_ecolicore_mcpam, compare_mu_for_different_sensitivities_ecolicore_pam, perform_single_gene_ko_for_all_genes, perform_and_plot_single_KO) -from cobra import Model, Reaction, Metabolite -import pandas as pd +from Scripts.mcpam_generation_uniprot_id import set_up_ecoli_mcpam, set_up_ecoli_pam import matplotlib.pyplot as plt; plt.rcdefaults() -import numpy as np -import matplotlib.pyplot as plt if __name__ == "__main__": - mcpam_core = set_up_ecolicore_mcpam(sensitivity=True) - mcpam_core.optimize() - print(mcpam_core.objective.value) - # pam_core = set_up_ecolicore_pam(sensitivity=True) - # - # models = [pam_core, mcpam_core] - # - # run_simulations_pam_mcpam_core(models) + mcpam = set_up_ecolicore_mcpam(sensitivity=False) + pam = set_up_ecolicore_pam(sensitivity=False) + pam.optimize() + mcpam.optimize() + print("PAM objective value: ", pam.objective.value) + print("mcPAM objective value: ", mcpam.objective.value) + models = [pam, mcpam] + run_simulations_pam_mcpam(models, type="core") diff --git a/Scripts/mcpam_generation_uniprot_id.py b/Scripts/mcpam_generation_uniprot_id.py index a484ada..ff926f3 100644 --- a/Scripts/mcpam_generation_uniprot_id.py +++ b/Scripts/mcpam_generation_uniprot_id.py @@ -359,7 +359,7 @@ def set_up_ecoli_mcpam(total_protein: Union[bool, float] = True, active_enzymes: unused_sector=unused_protein_sector, sensitivity=sensitivity, configuration = config, membrane_sector=membrane_sector ) - pamodel.solver = 'glpk' + return pamodel def _parse_enzyme_information_from_file(file_path:str): diff --git a/Scripts/mcpam_simulations_analysis.py b/Scripts/mcpam_simulations_analysis.py index 6c0504c..df2519b 100644 --- a/Scripts/mcpam_simulations_analysis.py +++ b/Scripts/mcpam_simulations_analysis.py @@ -202,35 +202,39 @@ def run_pam_mcpam_core_with_optimized_kcats(sensitivity:bool=True, return pam_core, mcpam_core -def run_simulations_pam_mcpam_core(models, print_area:bool=False): +def run_simulations_pam_mcpam(models, print_area:bool=False, type:str="full scale"): fontsize = 25 labelsize = 15 - #### Fluxes simulations for different glc uptake rates + # load phenotype data from excel file pt_data = pd.read_excel(os.path.join('Data', 'Ecoli_phenotypes','Ecoli_phenotypes_py_rev.xls'), sheet_name='Yields', index_col=None) + # Define the biomass name based on the used model + if type == "full scale": + biomass_name = 'BIOMASS_Ec_iML1515_core_75p37M' + max_area_list = np.linspace(0.1, 0.4, 4) + else: + biomass_name = 'BIOMASS_Ecoli_core_w_GAM' + max_area_list = np.linspace(0.01, 0.04, 4) + # extract reaction specific data rxn_to_pt = {} rxn_transform = { 'EX_ac_e': 'EX_ac_e', 'EX_co2_e': 'EX_co2_e', 'EX_o2_e': 'EX_o2_e', - 'BIOMASS_Ecoli_core_w_GAM':'BIOMASS_Ec_iML1515_core_75p37M', - # 'BIOMASS_Ec_iML1515_core_75p37M' :'BIOMASS_Ec_iML1515_core_75p37M' + biomass_name:'BIOMASS_Ec_iML1515_core_75p37M' } for rxn_id, pt_id in rxn_transform.items(): rxn_to_pt[rxn_id] = pt_data[['EX_glc__D_e', pt_id]].dropna().rename(columns={pt_id: rxn_id}) - glc_uptake_rates = np.linspace(0.5, 11.5, 20) + glc_uptake_rates = np.linspace(0.5, 14, 25) - # Initializing fluxes and concentrations for pam and mcpam core + # Initializing fluxes and concentrations for pam and mcpam fluxes_dict = {} concentrations_dict = {} - # Creating a list of different max area for mcpam - max_area_list = np.linspace(0.01, 0.04, 4) - for model, config in zip(models, ["PAM", "mcPAM"]): # disable pyruvate formate lyase (inhibited by oxygen) @@ -288,11 +292,11 @@ def run_simulations_pam_mcpam_core(models, print_area:bool=False): sns.set_palette(("colorblind")) # plot flux changes with glucose uptake - rxn_id = ['EX_ac_e', 'EX_co2_e', 'EX_o2_e', 'BIOMASS_Ecoli_core_w_GAM'] + rxn_id = ['EX_ac_e', 'EX_co2_e', 'EX_o2_e', biomass_name] ax_title = {'EX_ac_e': 'Acetate Secretion', 'EX_co2_e': 'CO2 Secretion', 'EX_o2_e': 'O2 uptake', - 'BIOMASS_Ecoli_core_w_GAM': 'Biomass Production'} + biomass_name: 'Biomass Production'} # rxn_id = ['EX_ac_e', 'EX_co2_e', 'EX_o2_e', 'BIOMASS_Ec_iML1515_core_75p37M'] fig, axs = plt.subplots(2, 2, dpi=90) for r, ax in zip(rxn_id, axs.flatten()): @@ -328,32 +332,6 @@ def run_simulations_pam_mcpam_core(models, print_area:bool=False): # show legend fig.subplots_adjust(top=0.85, wspace=0.5, hspace=0.5) - # Occupied, available area calculation for mcpam - mcpam_core = models[1] - mcpam_core.optimize() - occupied_area, available_area = mcpam_core.sectors.get_by_id("MembraneSector").calculate_occupied_membrane(mcpam_core) - - # Print the objective values for glc uptake rate 10 [mmol/gDW/h] - for model, config in zip(models, ["pam_core", "mcpam_core"]): - model.reactions.EX_glc__D_e.lower_bound = -10 - if config == "mcpam_core": - model.sectors.get_by_id("MembraneSector")._update_membrane_constraint(0.033, model) - print(f'{config} objective value for glc uptake rate 10 [mmol/gDW/h]:', model.objective.value) - print(f"mcPAM occupied area = {occupied_area} um2") - print(f"mcPAM available area = {available_area} um2") - - #Add parameter box - param_text = ( - f"mcPAM/PAM_core model \n" - f"Total protein: {mcpam_core.total_protein_fraction} g/g DW \n" - f"Max membrane area = {mcpam_core.membrane_sector.max_membrane_area * 100}% \n" - f"mcPAM occupied area = {occupied_area/available_area*100}%" - ) - - #Position the text box in the upper left corner of each subplot - fig.text(0.6, 2.75, 'param_text', transform=ax.transAxes, fontsize=fontsize, - verticalalignment='top', bbox=dict(facecolor='white', alpha=0.5, boxstyle="round,pad=0.3")) - plt.show() def perform_single_gene_ko_for_all_genes(model): diff --git a/src/PAModelpy/MembraneSector.py b/src/PAModelpy/MembraneSector.py index 36d4bc2..7259a93 100644 --- a/src/PAModelpy/MembraneSector.py +++ b/src/PAModelpy/MembraneSector.py @@ -15,7 +15,7 @@ def __init__( alpha_numbers_dict: {}, enzyme_location: {}, cog_class: {} = None, - max_area: float = 0.033, + max_area: float = 0.27, configuration=Config): self.id = 'MembraneSector'