run_simulation_for_pam_mcpam function is now updated. can now used to…

… run full scale or core models

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")
 
 
 

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):

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):

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'