8. Low Level API¶
Try it with Binder¶
Recall you can also try these examples with Binder.
1. core¶
Import core modules¶
In order to use contrabass library you have to import the following core module.
[1]:
from contrabass.core import CobraMetabolicModel
Download example model¶
[8]:
import urllib.request
FILE = "aureus.xml"
MODEL_URL = 'https://raw.githubusercontent.com/openCONTRABASS/CONTRABASS/main/docs/source/aureus.xml'
urllib.request.urlretrieve(MODEL_URL, FILE)
[8]:
('aureus.xml', <http.client.HTTPMessage at 0x7f21ef540510>)
Read metabolic model¶
Read input metabolic model in Systems Biology Markup Language (SBML) format. SBML is an XML-based standard for systems biology model’s exchange.
Allowed file formats are: * xml * json * yml
[7]:
model = CobraMetabolicModel("aureus.xml")
Get model info¶
The following methods allow us to print on the command line data from the model.
Model info¶
[3]:
model.print_model_info()
MODEL INFO
-------------------------------------------------------
MODEL: MODEL1507180070
REACTIONS: 743
METABOLITES: 655
GENES: 619
COMPARTMENTS: c
e
Metabolites¶
[4]:
model.print_metabolites()
MODEL: MODEL1507180070 - NUMBER OF METABOLITES: 655
METABOLITE | COMPARTMENT | REACTION ID
-------------------------------------------------------
10fthf_c | c | MTHFC
| | biomass_SA_7b
| | GARFT
| | biomass_SA_8a
| | FTHFL
| | biomass_SA_7a
| | AICART
12dgr_EC_c | c | biomass_SA_3a
| | biomass_SA_lipids_only
| | biomass_SA_2a
| | biomass_SA_3b
|
Reactions¶
[5]:
model.print_reactions()
MODEL: MODEL1507180070 - NUMBER OF REACTIONS: 743
REACTION ID | UPPER BOUND | LOWER BOUND | REACTION
-------------------------------------------------------
3M2OBLOXRD | 999999.0 | -999999.0 | 3mob_c + h_c + lpam_c <=> 2mpdhl_c + co2_c
3M2OPLOXRD | 999999.0 | -999999.0 | 3mop_c + h_c + lpam_c <=> 2mbdhl_c + co2_c
4M2OPLOXRD | 999999.0 | -999999.0 | 4mop_c + h_c + lpam_c <=> 3mbdhl_c + co2_c
6PGALSZ | 999999.0 | 0.0 | h2o_c + lac6p_c --> dgal6p_c + glc_DASH_D_c
6PHBG | 999999.0 | 0.0 | h2o_c + salc6p_c --> 2hymeph_c + g6p_c
ABTAr | 999999.0 | 0.0 | 4abut_c + akg_c --> glu_DASH_L_c + sucsal_c
ACACT1r | 999999.0
Genes¶
[6]:
model.print_genes()
MODEL: MODEL1507180070 - NUMBER OF GENES: 619
GENE ID | GENE NAME | REACTION ID | GPR RELATION
-------------------------------------------------------
SA0008 | SA0008 | HISDr | SA0008
SA0009 | SA0009 | SERTRS | SA0009
SA0011 | SA0011 | HSERTA | SA0011
SA0016 | SA0016 | ADSS | SA0016
SA0036 | SA0036 | GPDDA2 | SA0036 or SA0820 or SA1542 or SA0969 or SA0220
| | GPDDA5 | SA0036 or SA0820 or SA1542 or SA0969 or SA0220
| | GPDDA4 | SA0036 or SA0820 or SA1542 or SA0969 or SA0220
|
Dead-end metabolites¶
Dead-end metabolites are those metabolites which are not consumed or not produced by any reaction of a given compartment of the model, including exchange reactions.
Finding Dead-end metabolites¶
find_dem(). Method dem() returns a python dict with the following content:[7]:
model.find_dem()
model.dem()
[7]:
{'c': [<Metabolite uppg1_c at 0x7f288f6b0110>,
<Metabolite mi1p_DASH_D_c at 0x7f288f67a110>,
<Metabolite malm_c at 0x7f288f6741d0>,
<Metabolite nal2a6o_c at 0x7f288f67e250>,
<Metabolite ssaltpp_c at 0x7f288f69e390>,
<Metabolite trnaala_c at 0x7f288f6a43d0>,
<Metabolite trnaarg_c at 0x7f288f6a44d0>,
<Metabolite adcobdam_c at 0x7f2858f96510>,
<Metabolite trnaasn_c at 0x7f288f6a4590>,
<Metabolite 2ombz_c at 0x7f2858fe2610>,
<Metabolite trnaasp_c at 0x7f288f6a4650>,
<Metabolite trnacys_c at 0x7f288f6a4750>,
<Metabolite pala_SA_c at 0x7f28b43a8810>,
<Metabolite 2pglyc_c at 0x7f2858fe28d0>,
<Metabolite trnagly_c at 0x7f288f6a48d0>,
<Metabolite sucr_c at 0x7f288f69e8d0
Dead-end metabolites can be printed on the command line with print_dem() method.
[8]:
model.print_dem()
MODEL: MODEL1507180070 - NUMBER OF DEM: 2 - COMPARTMENT: ALL
METABOLITE | COMPARTMENT | REACTION ID
-------------------------------------------------------
uppg1_c | c | UPPDC2
mi1p_DASH_D_c | c | MI1PP
malm_c | c | PYRZAM
nal2a6o_c | c | NALN6
ssaltpp_c | c | SHCHCS2
trnaala_c | c | ALATRS
trnaarg_c | c | ARGTRS
adcobdam_c | c | ADCYRS
trnaasn_c | c | ASNTRS
2ombz_c | c | URFGTT
trnaasp_c | c | ASPTRS
trnacys_c | c | CYSTRS
pala_SA_c | c
Printed dead-end metabolites can also be limited to a specific compartment. Available compartments are given by model.compartments() method.
[9]:
model.print_dem(compartment="e")
MODEL: MODEL1507180070 - NUMBER OF DEM: 2 - COMPARTMENT: e
METABOLITE | COMPARTMENT | REACTION ID
-------------------------------------------------------
zn2_e | e | EX_zn2_e
| | ZNabc
spmd_e | e | EX_spmd_e
| | SPMDabc
salcn_e | e | EX_salcn_e
| | SALCpts
tre_e | e | EX_tre_e
| | TREpts
mn2_e | e | EX_mn2_e
| | MNabc
| | MNt2
ura_e | e | URAt2
|
Removing dead-end metabolites¶
Method remove_dem() removes dead-end metabolites from the model. Once dead-end metabolites are deleted, some reactions might not produce a metabolite anymore. This method deletes also these reactions and loops again until no dead-end metabolite is found:
function remove_dead_end_metabolites(model)
do while previous_metabolites != current_metabolites:
previous_metabolites = model metabolites
delete all metabolites in find_dead_end_metabolites(model)
for reaction that produced or consumed dead-end metabolites:
if reaction produces or consumes 0 metabolites [and is not exchange nor demand]:
delete reaction on model
current_metabolites = model metabolites
return model
Reactions that are deleted with this method can be asjusted with the following params:
keep_all_exchange_demand_reactions :
True: (default): A reaction is considered exchange or demand if it initially does not produce or consume any metabolite. This means the reaction represents an input or output of the metabolic network. If the parameter is true, none of these reactions will be removed during the execution.
False: If a reaction is a cobra boundary reaction (calculated with heuristics) that reaction will not be deleted during the execution. Note that exchange or demand reactions that are not identified by the cobrapy framework as boundary reactions will be deleted.
delete_exchange : Caution: this parameter is only intended for closed models without exchange or demand reactions, that is, models without reactions that do not produce or do not consume one metabolites on one of its sides. Note that most models do not match this property. Note also that this type of model will always yield 0 with Flux Balance Analysis.
True: all the reactions that produce or consume 0 metabolites are deleted whether they are exchange/demand or not.
False: (default) deleted according to ‘keep_all_incomplete_reactions’ param.
[10]:
print("Metabolites: ", len(model.metabolites()), "\nReactions: ", len(model.reactions()))
Metabolites: 655
Reactions: 743
[11]:
model.remove_dem()
[12]:
print("Metabolites: ", len(model.metabolites()), "\nReactions: ", len(model.reactions()))
Metabolites: 486
Reactions: 739
Dead reactions¶
Dead reactions are those reactions with upper and lower flux equal to zero.
Finding Dead reactions¶
Dead reactions of the model are calculated with method dead_reactions(), which returns a list of cobra.core.reactions representing dead reactions.
[13]:
model.dead_reactions()
[13]:
[<Reaction SA_biomass_1a at 0x7f288f110ad0>,
<Reaction biomass_SA_2a at 0x7f288f06c750>,
<Reaction biomass_SA_2b at 0x7f288f06c2d0>,
<Reaction biomass_SA_3a at 0x7f288f06c810>,
<Reaction biomass_SA_3b at 0x7f288f06cf90>,
<Reaction biomass_SA_4a at 0x7f288eff60d0>,
<Reaction biomass_SA_5a at 0x7f288eff6050>,
<Reaction biomass_SA_6a at 0x7f288eff6fd0>,
<Reaction biomass_SA_6b at 0x7f288f06cd50>,
<Reaction biomass_SA_7a at 0x7f288f01de90>,
<Reaction biomass_SA_7b at 0x7f288f01df10>,
<Reaction biomass_SA_lipids_only at 0x7f288f06cc10>,
<Reaction biomass_SA_nuc_only at 0x7f288f06c7d0>,
<Reaction biomass_SA_only_AA at 0x7f288f06c690>]
Chokepoint reactions¶
Chokepoint reactions are those reactions that are the only consumer or producer of a given metabolite that is not a dead-end metabolite.
Finding chokepoint reactions¶
find_chokepoints().chokepoints() returns a list of tuples of types (cobra.core.reactions, cobra.core.metabolites) representing chokepoint reactions and the metabolite of which they are the only consumer/producer.[14]:
# Read initial model again
model = CobraMetabolicModel("aureus.xml")
model.find_chokepoints()
model.chokepoints()
[14]:
[(<Reaction PAPA_SA at 0x7f288ec42750>,
<Metabolite 12dgr_SA_c at 0x7f288f683790>),
(<Reaction PROD2 at 0x7f288f4637d0>, <Metabolite 1pyr5c_c at 0x7f2858f9da50>),
(<Reaction DHDPRy at 0x7f288ee48a50>,
<Metabolite 23dhdp_c at 0x7f288f67a150>),
(<Reaction DHDPS at 0x7f288ee48a90>, <Metabolite 23dhdp_c at 0x7f288f67a150>),
(<Reaction DHAD1 at 0x7f288ee48cd0>, <Metabolite 23dhmb_c at 0x7f288f67a350>),
(<Reaction KARA1i at 0x7f288edae690>,
<Metabolite 23dhmb_c at 0x7f288f67a350>),
(<Reaction DHAD2 at 0x7f288ee48910>, <Metabolite 23dhmp_c at 0x7f288f67a750>),
(<Reaction KARA2i at 0x7f288edae890>,
<Metabolite 23dhmp_c at 0x7f288f67a750>),
(<Reaction DHPPDA at 0x7f288ee4bcd0>,
<Met
Dead reactions can be excluded from the computation of chokepoints with parameter exclude_dead_reactions=True. This is strongly suggested as dead reactions are considered forward reactions by default an this can lead to misinterpretations.
[15]:
model.find_chokepoints(exclude_dead_reactions=True)
model.chokepoints()
[15]:
[(<Reaction PAPA_SA at 0x7f288ec42750>,
<Metabolite 12dgr_SA_c at 0x7f288f683790>),
(<Reaction PROD2 at 0x7f288f4637d0>, <Metabolite 1pyr5c_c at 0x7f2858f9da50>),
(<Reaction DHDPRy at 0x7f288ee48a50>,
<Metabolite 23dhdp_c at 0x7f288f67a150>),
(<Reaction DHDPS at 0x7f288ee48a90>, <Metabolite 23dhdp_c at 0x7f288f67a150>),
(<Reaction DHAD1 at 0x7f288ee48cd0>, <Metabolite 23dhmb_c at 0x7f288f67a350>),
(<Reaction KARA1i at 0x7f288edae690>,
<Metabolite 23dhmb_c at 0x7f288f67a350>),
(<Reaction DHAD2 at 0x7f288ee48910>, <Metabolite 23dhmp_c at 0x7f288f67a750>),
(<Reaction KARA2i at 0x7f288edae890>,
<Metabolite 23dhmp_c at 0x7f288f67a750>),
(<Reaction DHPPDA at 0x7f288ee4bcd0>,
<Met
Chokepoint reactions can also be printed on the command line with print_chokepoints.
[16]:
model.print_chokepoints()
MODEL: MODEL1507180070 - NUMBER OF CHOKEPOINTS: 448
METABOLITE ID | METABOLITE NAME | REACTION ID | REACTION NAME
------------------------------------------------------------
12dgr_SA_c | 1,2-Daicylglycerol (Saureus) | PAPA_SA | Phosphatidate phosphatase
1pyr5c_c | 1-Pyrroline-5-carboxylate | PROD2 | Proline dehydrogenase
23dhdp_c | 2,3-Dihydrodipicolinate | DHDPRy | dihydrodipicolinate reductase (NADPH)
23dhdp_c | 2,3-Dihydrodipicolinate | DHDPS | dihydrodipicolinate synthase
23dhmb_c | (R)-2,3-Dihydroxy-3-methylbutanoate | DHAD1 |
Flux Balance Analysis¶
Flux Balance Analysis (FBA) is a mathematical procedure to estimate the fluxes of reactions in a metabolic model.
Method get_growth() calculates the objective value (growth rate) that maximizes the objective function. This method uses cobra.core.model.slim_optimize and was obtained from cobra.flux_analysis.deletion.
The objective value obtained with FBA can be accessed with objective_value().
[17]:
model.get_growth()
model.objective_value()
[17]:
0.1580502916027849
The objective function that is maximized during FBA can be accessed with objective()
[18]:
model.objective()
[18]:
'1.0*biomass_SA_8a - 1.0*biomass_SA_8a_reverse_4ce77'
This objective function can be changed by another reaction of the model with set_objective(). This method receives the id of the reaction that will be set as the new objective value.
[19]:
model.set_objective("DHAD1")
model.objective()
[19]:
'1.0*DHAD1 - 1.0*DHAD1_reverse_39dca'
Flux Variability Analysis¶
Flux Vatiability Analysis (FVA) is a mathematical procedure used to calculate the ‘’minimum and maximum flux for reactions in the network while maintaining some state of the network, e.g., supporting 90% of maximal possible biomass production rate’’.
fva() runs Flux Variability Analysis on the model. This method runs cobra.flux_analysis.variability and as so, it allows the same parameters (see the previous link):Method fva() returns an error list if there was an error while running FVA or an empty list [] otherwise.
[20]:
model.fva(threshold=0.95)
[20]:
[]
[21]:
model.fva(threshold=0.95, verbose=True)
FLUX VARIABILITY ANALYSIS: MODEL1507180070
REACTION: 3-Methyl-2-oxobutanoate:lipoamide oxidoreductase(decarboxylating and acceptor-2-methylpropanoylating)
fva ranges: [ 0.0 , 0.0 ]
REACTION: 3-Methyl-2-oxopentanoate:lipoamide oxidoreductase(decarboxylating and acceptor-2-methylpropanoylating)
fva ranges: [ 0.0 , 0.0 ]
REACTION: 4-Methyl-2-oxopentanoate:lipoamide oxidoreductase(decarboxylating and acceptor-2-methylpropanoylating)
fva ranges: [ 0.0 , 0.0 ]
REACTION: 6-phospho-beta-galactosidase
fva ranges: [ 0.0 , 0.0 ]
REACTION: 6-phospho-beta-glucosidase
fva ranges: [ 0.0 , 0.0 ]
RE
The result obtained with FVA can be accessed with get_fva(). This method returns a list of tuples: (cobra.core.reaction, [float] maximum flux, [float] minimum flux)
[22]:
model.get_fva()
[22]:
[(<Reaction 3M2OBLOXRD at 0x7f288eeab950>, 0.0, 0.0),
(<Reaction 3M2OPLOXRD at 0x7f288eeab990>, 0.0, 0.0),
(<Reaction 4M2OPLOXRD at 0x7f288eeabc10>, 0.0, 0.0),
(<Reaction 6PGALSZ at 0x7f288eeabad0>, 0.0, 0.0),
(<Reaction 6PHBG at 0x7f288eeabb10>, 0.0, 0.0),
(<Reaction ABTAr at 0x7f288eead890>, 0.0, 0.0),
(<Reaction ACACT1r at 0x7f288eeadd50>,
0.00032511320071648854,
-2065.6249999999704),
(<Reaction ACACT2r at 0x7f288eeadcd0>, 0.0, -2065.624999999973),
(<Reaction ACACT3r at 0x7f288eead710>, 0.0, -2065.624999999973),
(<Reaction ACACT4r at 0x7f288eead150>, 0.0, -2065.624999999973),
(<Reaction ACACT5r at 0x7f288eeadd10>, 0.0, -2065.624999999973),
(<Reaction ACACT6r at 0x7f288eead8
Updating flux of reactions with FVA¶
contrabass adds an extra parameter to this method which is update_flux. If True the method fva() updates the fluxes of reactions with the maximum and minimum flux values obtained with FVA (see the example below).
[23]:
# Print initial reactions of the model with initial flux values.
model.print_reactions()
MODEL: MODEL1507180070 - NUMBER OF REACTIONS: 743
REACTION ID | UPPER BOUND | LOWER BOUND | REACTION
-------------------------------------------------------
3M2OBLOXRD | 999999.0 | -999999.0 | 3mob_c + h_c + lpam_c <=> 2mpdhl_c + co2_c
3M2OPLOXRD | 999999.0 | -999999.0 | 3mop_c + h_c + lpam_c <=> 2mbdhl_c + co2_c
4M2OPLOXRD | 999999.0 | -999999.0 | 4mop_c + h_c + lpam_c <=> 3mbdhl_c + co2_c
6PGALSZ | 999999.0 | 0.0 | h2o_c + lac6p_c --> dgal6p_c + glc_DASH_D_c
6PHBG | 999999.0 | 0.0 | h2o_c + salc6p_c --> 2hymeph_c + g6p_c
ABTAr | 999999.0 | 0.0 | 4abut_c + akg_c --> glu_DASH_L_c + sucsal_c
ACACT1r | 999999.0
[24]:
# Update reactions flux values with FVA
model.fva(update_flux=True)
[24]:
[]
[25]:
# Print reactions of the model with constrained reactions flux values.
model.print_reactions()
MODEL: MODEL1507180070 - NUMBER OF REACTIONS: 743
REACTION ID | UPPER BOUND | LOWER BOUND | REACTION
-------------------------------------------------------
3M2OBLOXRD | 0.0 | 0.0 | 3mob_c + h_c + lpam_c --> 2mpdhl_c + co2_c
3M2OPLOXRD | 0.0 | 0.0 | 3mop_c + h_c + lpam_c --> 2mbdhl_c + co2_c
4M2OPLOXRD | 0.0 | 0.0 | 4mop_c + h_c + lpam_c --> 3mbdhl_c + co2_c
6PGALSZ | 0.0 | 0.0 | h2o_c + lac6p_c --> dgal6p_c + glc_DASH_D_c
6PHBG | 0.0 | 0.0 | h2o_c + salc6p_c --> 2hymeph_c + g6p_c
ABTAr | 0.0 | 0.0 | 4abut_c + akg_c --> glu_DASH_L_c + sucsal_c
ACACT1r | -3.26642
Essential genes¶
Essential genes are those genes that cause a zero growth rate when knocked out.
Finding essential genes¶
Essential genes are calculated with the method find_essential_genes_1(). This method returns an error list if there was an error during the computing or an empty list [] otherwise.
[26]:
model.find_essential_genes_1()
[26]:
[]
If essential genes were computed with the method above, they can be accessed with essential_genes(). This method return a list of cobra.core.gene with essential genes.
[27]:
model.essential_genes()
[27]:
{<Gene SA0183 at 0x7f2858fe2b10>,
<Gene SA0512 at 0x7f288eff6310>,
<Gene SA0728 at 0x7f288eea6890>,
<Gene SA0729 at 0x7f288eea6950>,
<Gene SA0823 at 0x7f288ee87d10>,
<Gene SA0910 at 0x7f288ee89890>,
<Gene SA0911 at 0x7f288ee89950>,
<Gene SA0912 at 0x7f288ee89a10>,
<Gene SA0913 at 0x7f288ee89ad0>,
<Gene SA0937 at 0x7f288eea5650>,
<Gene SA0938 at 0x7f288eea5710>,
<Gene SA0994 at 0x7f288ee8d350>,
<Gene SA0995 at 0x7f288ee8d410>,
<Gene SA0996 at 0x7f288ee8d4d0>,
<Gene SA1088 at 0x7f288ee8b4d0>,
<Gene SA1089 at 0x7f288ee8b590>,
<Gene SA1244 at 0x7f288ee91ad0>,
<Gene SA1245 at 0x7f288ee91b90>,
<Gene SA1255 at 0x7f288ee91dd0>,
<Gene SA1521 at 0x7f288ee93410>,
<Gene SA1585 at 0x7
Essential reactions¶
Essential reactions are those genes that cause a zero growth rate when knocked out.
Finding essential reactions¶
Essential reactions are calculated with the method compute_essential_reactions(). This method returns an error list if there was an error during the computing, or an empty list [] otherwise.
[31]:
model.compute_essential_reactions()
If essential reactions were computed with the method above, they can be accessed with essential_reactions(). This method return a dict with keys cobra.core.reaction with all the reactions of the model, and values float with the result of computing FBA with the reaction knocked-out.
[32]:
model.essential_reactions()
[32]:
[<Reaction CO2t at 0x7f288ee40750>,
<Reaction KARA1i at 0x7f288edae690>,
<Reaction DHAD1 at 0x7f288ee48cd0>,
<Reaction PDHcr at 0x7f288ec4cc10>,
<Reaction EX_pro_DASH_L_e at 0x7f288edf5350>,
<Reaction SUCOAS at 0x7f288eead690>,
<Reaction PFK at 0x7f288ec4cf50>,
<Reaction FBA at 0x7f288edfb3d0>,
<Reaction AKGDa at 0x7f288ee34c90>,
<Reaction GLCpts at 0x7f288ee25590>,
<Reaction PGI at 0x7f288f537490>,
<Reaction FUM at 0x7f288ee0a650>,
<Reaction EX_co2_e at 0x7f288ee5f850>,
<Reaction PROD2 at 0x7f288f4637d0>,
<Reaction EX_glc_DASH_D_e at 0x7f288ee67550>,
<Reaction ACLS at 0x7f288ee84d50>,
<Reaction AKGDb at 0x7f288ee34d50>,
<Reaction EX_o2_e at 0x7f288ee6dcd0>,
<Reaction TPI a
Essential genes reactions¶
Essential genes reactions are those reactions that are knocked-out when an essential gene is knocked-out.
Finding essential genes reactions¶
Essential genes reactions are calculated with the method find_essential_genes_reactions(). This method returns an error list if there was an error during the computing or an empty list [] otherwise.
[33]:
model.find_essential_genes_reactions()
[33]:
[]
If essential genes reactions were computed with the method above, they can be accessed with essential_genes_reactions(). This method returns a dict with keys cobra.core.reaction with all the reactions of the model, and values a list of cobra.core.genes with the essential genes that cause the knock-out of the reaction.
[34]:
model.essential_genes_reactions()
[34]:
{<Reaction PGK at 0x7f288f53c950>: [<Gene SA0728 at 0x7f288eea6890>],
<Reaction CYTBD at 0x7f288ee42e90>: [<Gene SA0910 at 0x7f288ee89890>,
<Gene SA0911 at 0x7f288ee89950>,
<Gene SA0912 at 0x7f288ee89a10>,
<Gene SA0913 at 0x7f288ee89ad0>,
<Gene SA0937 at 0x7f288eea5650>,
<Gene SA0938 at 0x7f288eea5710>],
<Reaction PROD2 at 0x7f288f4637d0>: [<Gene SA1585 at 0x7f288ee97110>],
<Reaction TPI at 0x7f288eb48190>: [<Gene SA0729 at 0x7f288eea6950>],
<Reaction DHAD2 at 0x7f288ee48910>: [<Gene SA1858 at 0x7f288ee83a10>],
<Reaction DHAD1 at 0x7f288ee48cd0>: [<Gene SA1858 at 0x7f288ee83a10>],
<Reaction AKGDb at 0x7f288ee34d50>: [<Gene SA1244 at 0x7f288ee91ad0>],
<Reaction ILETA at 0x7f2
2. core.Facade core.FacadeUtils¶
Modules FacadeUtils and Facade from contrabass provide specific methods for the generation of reports from the computation of critical reactions taking into account flux constraints.
Import FacadeUtils modules¶
Import the module
[35]:
from contrabass.core import Facade
from contrabass.core import FacadeUtils
Generate vulnerabilities summary spreadsheet¶
run_summary_model generates a spreadsheet file from a model file with the following data:Method declaration:
run_summary_model(self, model_path, print_f, arg1, arg2, objective=None, fraction=1.0)
model_path: Path of the SBML metabolic model file.print_f: Callback function to inform of the computation progression. The function must have a declaration like the following:custom_function_name(message, arg1, arg2)arg1: First parameter to pass to the callback function if any.arg2: Second parameter to pass to the callback function if any.objective: Reactions id to be used as objective function.fraction: Fraction of optmimum to be used with FVA.[36]:
def callback_print_ignore(message, arg1, arg2):
pass
def callback_print_logger(message, arg1, arg2):
logger = arg1
logger.info(message)
def callback_print(message, arg1, arg2):
print(arg1 + message)
facadeUtils = FacadeUtils()
spreadsheet = facadeUtils.run_summary_model("aureus.xml", callback_print, "LOG:", None, fraction=0.95)
facadeUtils.save_spreadsheet("output.xls", spreadsheet)
LOG:Reading model...
LOG:Generating models...
LOG:Searching Dead End Metabolites (D.E.M.)...
LOG:Searching chokepoint reactions...
LOG:Searching essential reactions...
LOG:Searching essential genes...
LOG:Searching essential genes reactions...
LOG:Removing Dead End Metabolites (D.E.M.)...
LOG:Searching essential reactions...
LOG:Searching new chokepoint reactions...
LOG:Searching essential genes...
LOG:Searching essential genes reactions...
LOG:Running Flux Variability Analysis...
LOG:Searching Dead End Metabolites (D.E.M.)...
LOG:Searching new chokepoint reactions...
LOG:Searching essential genes...
LOG:Searching essential genes reactions...
LOG:Searching essential reactions...
LOG:Removing
Generate growth-dependent reactions¶
The method generate_growth_dependent_report of Facade computes and returns growth dependent reactions.
The output can be saved in a spreadsheet file or html report. The resulting files include a report of the effect that varying values of optimal growth with FVA has on the folloing sets: - Reversible Reactions (RR) - Non-reversible Reactions (NR) - Dead Reactions (DR) - Chokepoint reactions (CP)
generate_growth_dependent_report(self, config)config: contrabass.core.utils.GrowthDependentCPConfigcontrabass.core.utils.GrowthDependentCPConfigmodel_path: Path of the SBML metabolic model file.print_f: Callback function to inform of the computation progression. The function must have a declaration like the following: custom_function_name(message, arg1, arg2)arg1: First parameter to pass to the function if any.arg2: Second parameter to pass to the function if any.objective: Reactions id to be used as objective function.fraction: Fraction of optmimum to be used with FVA.output_path_spreadsheet: Output path to save spreadsheet report. If None no file is generated.output_path_html: Output path to save html report. If None no file is generated.[37]:
# import config class
from contrabass.core.utils.GrowthDependentCPConfig import GrowthDependentCPConfig
def callback_print(message, arg1, arg2):
print(message)
config = GrowthDependentCPConfig()
config.print_f = callback_print
config.model_path = "aureus.xml"
config.output_path_spreadsheet = "output.xls"
config.output_path_html = "report.html"
facade = Facade()
facade.generate_growth_dependent_report(config)
Reading model...
Computing essential reactions...
Computing optimal growth essential reactions
Computing growth dependent essential reactions
Running Flux Variability Analysis with fraction: 0.0
Running Flux Variability Analysis with fraction: 0.1
Running Flux Variability Analysis with fraction: 0.2
Running Flux Variability Analysis with fraction: 0.3
Running Flux Variability Analysis with fraction: 0.4
Running Flux Variability Analysis with fraction: 0.5
Running Flux Variability Analysis with fraction: 0.6
Running Flux Variability Analysis with fraction: 0.7
Running Flux Variability Analysis with fraction: 0.8
Running Flux Variability Analysis with fraction: 0.9
Running Flux Variability Ana
[ ]: