reconstruction.ecoli.dataclasses.process.metabolism
SimulationData for metabolism process
TODO:
improved estimate of ILE/LEU abundance or some external data point
implement L1-norm minimization for AA concentrations
find concentration for PI[c]
add (d)NTP byproduct concentrations
- class reconstruction.ecoli.dataclasses.process.metabolism.ConcentrationUpdates(concDict, relative_changes, equilibriumReactions, exchange_data_dict, all_metabolite_ids, linked_metabolites)[source]
Bases:
object
- _exchange_flux_present(exchange_data)[source]
Caches the presence of exchanges in each media condition based on exchange_data to set concentrations in concentrations_based_on_nutrients().
- Parameters:
exchange_data (dict[str, Any]) –
dictionary of exchange data for all media conditions with keys:
{importUnconstrainedExchangeMolecules (dict[str, set[str]]): exchange molecules (with location tag) for each media key that do not have an upper bound on their flux, importConstrainedExchangeMolecules (dict[str, dict[str, float with mol/mass/time units]]): constrained molecules (with location tag) for each media key with upper bound flux constraints}
- Returns:
Sets of molecules IDs (with location tags) that can be imported for each media ID
- Return type:
- reconstruction.ecoli.dataclasses.process.metabolism.DRY_MASS_UNITS = 1 [fg]
Units for dry mass
- exception reconstruction.ecoli.dataclasses.process.metabolism.InvalidReactionDirectionError[source]
Bases:
Exception
- reconstruction.ecoli.dataclasses.process.metabolism.K_CAT_UNITS = 1.0 [1/s]
Units for k cat values
- reconstruction.ecoli.dataclasses.process.metabolism.METABOLITE_CONCENTRATION_UNITS = 1.0 [mol/L]
Units for metabolite concentrations
- class reconstruction.ecoli.dataclasses.process.metabolism.Metabolism(raw_data, sim_data)[source]
Bases:
object
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
sim_data (SimulationDataEcoli) – Simulation data object
- solver
solver ID, should match a value in modular_fba.py, set by
_set_solver_values()
- Type:
- kinetic_objective_weight
weighting for the kinetic objective, 1-weighting for the homeostatic objective, set by
_set_solver_values()
- Type:
- kinetic_objective_weight_in_range
weighting for deviations from the kinetic target within min and max ranges, set by
_set_solver_values()
- Type:
- secretion_penalty_coeff
penalty on secretion fluxes, set by
_set_solver_values()
- Type:
- metabolite_charge
mapping of metabolite IDs to charge, set by
_add_metabolite_charge()
- concentration_updates
- kinetic_constraint_reactions
- kinetic_constraint_enzymes
- kinetic_constraint_substrates
- _kcats
- _saturations
- _enzymes
- constraint_is_kcat_only
- _compiled_enzymes
- _compiled_saturation
- reaction_stoich
- maintenance_reaction
- reaction_catalysts
- catalyst_ids
- reactions_with_catalyst
- catalysis_matrix_I
- catalysis_matrix_J
- catalysis_matrix_V
- use_all_constraints
- constraints_to_disable
- base_reaction_ids
- reaction_id_to_base_reaction_id
- amino_acid_export_kms
- transport_reactions
transport reaction IDs in the metabolic network (includes reverse reactions and reactions with kinetic constraints), set by
_build_transport_reactions()
- aa_synthesis_pathways
data for allosteric inhibition of amino acid pathways indexed by amino acid ID with location tag and nested dictionary with the following keys:
{'enzymes' (str): limiting/regulated enzyme ID in synthesis pathway with location tag, 'kcat_data' (units.Unum): kcat associated with enzyme reaction with units of 1/time, 'ki' (Tuple[units.Unum, units.Unum]]): lower and upper limits of KI associated with enzyme reaction with units of mol/volume}
Set by
_build_amino_acid_pathways()
- KI_aa_synthesis
KI for each AA for synthesis portion of supply (in units of
METABOLITE_CONCENTRATION_UNITS
), set byset_phenomological_supply_constants()
- Type:
- KM_aa_export
KM for each AA for export portion of supply (in units of
METABOLITE_CONCENTRATION_UNITS
), set byset_phenomological_supply_constants()
- Type:
- fraction_supply_rate
fraction of AA supply that comes from a base synthesis rate, set by
set_phenomological_supply_constants()
- Type:
- fraction_import_rate
fraction of AA supply that comes from AA import if nutrients are present, set by
set_phenomological_supply_constants()
- Type:
- ppgpp_synthesis_reaction
reaction ID for ppGpp synthesis (catalyzed by RelA and SpoT), set by
_build_ppgpp_reactions()
- Type:
- ppgpp_degradation_reaction
reaction ID for ppGpp degradation (catalyzed by SpoT), set by
_build_ppgpp_reactions()
- Type:
- ppgpp_reaction_names
names of reaction involved in ppGpp, set by
_build_ppgpp_reactions()
- ppgpp_reaction_metabolites
names of metabolites in ppGpp reactions, set by
_build_ppgpp_reactions()
- ppgpp_reaction_stoich
2D array with metabolites on rows and reactions on columns containing the stoichiometric coefficient, set by
_build_ppgpp_reactions()
- Type:
- aa_to_exporters
dictonary that maps aa to transporters involved in export reactions. Set by
set_mechanistic_export_constants()
.
- aa_to_exporters_matrix
correlation matrix. Columns correspond to exporting enzymes and rows to amino acids. Set by
set_mechanistic_export_constants()
.- Type:
- aa_exporter_names
names of all exporters. Set by
set_mechanistic_export_constants()
.- Type:
- aa_export_kms
kms corresponding to generic transport/enzyme reactions for each AA in concentration units. Set by
set_mechanistic_export_constants()
.- Type:
- export_kcats_per_aa
kcats corresponding to generic export reactions for each AA. Units in counts/second. Set by
set_mechanistic_export_constants()
andset_mechanistic_export_constants()
.- Type:
- aa_to_importers
dictonary that maps aa to transporters involved in import reactions. Set by
set_mechanistic_uptake_constants()
.
- aa_to_importers_matrix
correlation matrix. Columns correspond to importing enzymes and rows to amino acids. Set by
set_mechanistic_export_constants()
.- Type:
- aa_importer_names
names of all importers. Set by
set_mechanistic_export_constants()
.- Type:
- import_kcats_per_aa
kcats corresponding to generic import reactions for each AA. Units in counts/second. Set by
set_mechanistic_export_constants()
.- Type:
- aa_enzymes
enzyme ID with location tag for each enzyme that can catalyze an amino acid pathway with
enzyme_to_amino_acid
mapping these to each amino acid. Set byset_mechanistic_supply_constants()
.- Type:
- aa_kcats_fwd
forward kcat value for each synthesis pathway in units of
K_CAT_UNITS
, ordered by amino acid molecule group. Set byset_mechanistic_supply_constants()
.- Type:
- aa_kcats_rev
reverse kcat value for each synthesis pathway in units of
K_CAT_UNITS
, ordered by amino acid molecule group. Set byset_mechanistic_supply_constants()
.- Type:
- aa_kis
KI value for each synthesis pathway in units of
METABOLITE_CONCENTRATION_UNITS
, ordered by amino acid molecule group. Will be inf if there is no inhibitory control. Set byset_mechanistic_supply_constants()
.- Type:
- aa_upstream_kms
KM value associated with the amino acid that feeds into each synthesis pathway in units of
METABOLITE_CONCENTRATION_UNITS
, ordered by amino acid molecule group. Will be 0 if there is no upstream amino acid considered. Set byset_mechanistic_supply_constants()
.
- aa_reverse_kms
KM value associated with the amino acid in each synthesis pathway in units of
METABOLITE_CONCENTRATION_UNITS
, ordered by amino acid molecule group. Will be inf if the synthesis pathway is not reversible. Set byset_mechanistic_supply_constants()
.- Type:
- aa_upstream_mapping
index of the upstream amino acid that feeds into each synthesis pathway, ordered by amino acid molecule group. Set by
set_mechanistic_supply_constants()
.- Type:
- enzyme_to_amino_acid
relationship mapping from aa_enzymes to amino acids (n enzymes, m amino acids). Will contain a 1 if the enzyme associated with the row can catalyze the pathway for the amino acid associated with the column. Set by
set_mechanistic_supply_constants()
.- Type:
- aa_forward_stoich
relationship mapping from upstream amino acids to downstream amino acids (n upstream, m downstream). Will contain a -1 if the amino acid associated with the row is required for synthesis of the amino acid associated with the column. Set by
set_mechanistic_supply_constants()
.- Type:
- aa_reverse_stoich
relationship mapping from upstream amino acids to downstream amino acids (n downstream, m upstream). Will contain a -1 if the amino acid associated with the row is produced through a reverse reaction from the amino acid associated with the column. Set by
set_mechanistic_supply_constants()
.- Type:
- aa_import_kis
inhibition constants for amino acid import based on the internal amino acid concentration
- Type:
- specific_import_rates
import rates expected in rich media conditions for each amino acid normalized by dry cell mass in units of
K_CAT_UNITS
/DRY_MASS_UNITS
, ordered by amino acid molecule group. Set byset_mechanistic_supply_constants()
.- Type:
- max_specific_import_rates
max import rates for each amino acid without including internal concentration inhibition normalized by dry cell mass in units of
K_CAT_UNITS
/DRY_MASS_UNITS
, ordered by amino acid molecule group. Set byset_mechanistic_supply_constants()
.- Type:
- _add_metabolite_charge(raw_data)[source]
Compiles all metabolite charges.
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
- Attributes set:
- _build_amino_acid_pathways(raw_data, sim_data)[source]
Creates mapping between enzymes and amino acid pathways with allosteric inhibition feedback from the amino acid.
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
sim_data (SimulationDataEcoli) – Simulation data object
- Attributes set:
- _build_biomass(raw_data, sim_data)[source]
Calculates metabolite concentration targets.
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
sim_data (SimulationDataEcoli) – Simulation data object
- Attributes set:
- _build_linked_metabolites(raw_data, conc_dict)[source]
Calculates ratio between linked metabolites to keep it constant throughout a simulation.
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
conc_dict (dict[str, Unum]) – Mapping of metabolite IDs to homeostatic concentrations calculated by
_build_biomass()
- Returns:
Mapping from a linked metabolite to its lead metabolite and concentration ratio to be maintained:
{'lead' (str): metabolite to link the concentration to, 'ratio' (float): ratio to multiply the lead concentration by}
- Return type:
- _build_metabolism(raw_data, sim_data)[source]
Build the matrices/vectors for metabolism (FBA) Reads in and stores reaction and kinetic constraint information
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
sim_data (SimulationDataEcoli) – Simulation data object
- Attributes set:
- _build_ppgpp_reactions(raw_data, sim_data)[source]
Creates structures for ppGpp reactions for use in polypeptide_elongation.
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
sim_data (SimulationDataEcoli) – Simulation data object
- _build_transport_reactions(raw_data, sim_data)[source]
Creates list of transport reactions that are included in the reaction network.
- Parameters:
raw_data (KnowledgeBaseEcoli) – Raw data object
sim_data (SimulationDataEcoli) – Simulation data object
- Attributes set:
- static _construct_default_saturation_equation(mets, kms, kis, known_mets)[source]
- Parameters:
- Returns:
Saturation equation with metabolites to replace delimited by double quote (e.g. “metabolite”)
- Return type:
- static _extract_custom_constraint(constraint, reactant_tags, product_tags, known_mets)[source]
- Parameters:
values defining a kinetic constraint:
{'customRateEquation' (str): mathematical representation of rate (must contain 'kcat*E'), 'customParameterVariables' (dict[str, str]): mapping of variable names in the rate equation to metabolite IDs without location tags (must contain 'E'), 'customParameterConstants' (list[str]): constant strings in the rate equation that correspond to values (must contain 'kcat'), 'customParameterConstantValues' (list[float]): values for each of the constant strings, 'Temp' (float or ''): temperature of measurement}
reactant_tags (dict[str, str]) – mapping of molecule IDs without a location tag to molecule IDs with a location tag for all reactants
product_tags (dict[str, str]) – mapping of molecule IDs without a location tag to molecule IDs with a location tag for all products
known_mets (set[str]) – molecule IDs with a location tag for molecules with known concentrations
- Returns:
2-element tuple containing
kcats: temperature adjusted kcat value, in units of 1/s
saturation: saturation equation with metabolites to replace delimited by double quote (eg. “metabolite”)
- Return type:
- _is_transport_rxn(stoich)[source]
Determines if the metabolic reaction with a given stoichiometry is a transport reactions that transports metabolites between different compartments. A metabolic reaction is considered to be a transport reaction if the substrate set and the product share the same metabolite tagged into different compartments.
- static _lambdify_constraints(constraints)[source]
Creates str representations of kinetic terms to be used to create kinetic constraints that are returned with getKineticConstraints().
- Parameters:
constraints (dict[str, Any]) –
valid kinetic constraints for each reaction:
{reaction ID: { 'enzyme': enzyme catalyst (str), 'kcat': kcat values (list[float]), 'saturation': saturation equations (list[str]) }}
- Returns:
7-element tuple containing
rxns: sorted reaction IDs for reactions with a kinetic constraint
enzymes: sorted enzyme IDs for enzymes that catalyze a kinetic reaction
substrates: sorted substrate IDs for substrates that are needed for kinetic saturation terms
all_kcats: (n rxns, 3) min, mean and max kcat value for each reaction
all_saturations: sympy str representation of a list of saturation terms (eg. ‘[s[0] / (1 + s[0]), 2 / (2 + s[1])]’)
all_enzymes: sympy str representation of enzymes for each reaction (e.g. ‘[e[0], e[2], e[1]]’)
constraint_is_kcat_only: True if reaction only has kcat values and no saturation terms
- Return type:
tuple[list[str], list[str], list[str], ndarray[Any, dtype[float64]], str, str, ndarray[Any, dtype[bool_]]]
- static _replace_enzyme_reactions(constraints, stoich, rxn_catalysts, reversible_rxns, rxn_id_to_compiled_id)[source]
Modifies reaction IDs in data structures to duplicate reactions with kinetic constraints and multiple enzymes.
- Parameters:
constraints (dict[tuple[str, str], dict[str, list[Any]]]) –
valid kinetic constraints for each reaction/enzyme pair:
{(reaction ID, enzyme with location tag): { 'kcat': kcat values (list[float]), 'saturation': saturation equations (list[str]) }}
stoich (dict[str, dict[str, int]]) –
stoichiometry of metabolites for each reaction (if None, data is loaded from raw_data and sim_data):
{reaction ID: {metabolite ID with location tag: stoichiometry}}
rxn_catalysts (dict[str, list[str]]) –
enzyme catalysts for each reaction with known catalysts, likely a subset of reactions in stoich (if None, data is loaded from raw_data and sim_data):
{reaction ID: enzyme IDs with location tag}
reversible_rxns (list[str]) – reaction IDs for reactions that have a reverse complement, does not have reverse tag
rxn_id_to_compiled_id (dict[str, str]) – mapping from reaction IDs to the IDs of the original reactions they were derived from
- Returns:
5-element tuple containing
new_constraints: valid kinetic constraints for each reaction:
{reaction ID: { 'enzyme': enzyme catalyst (str), 'kcat': kcat values (list[float]), 'saturation': saturation equations (list[str]) }}
stoich: stoichiometry of metabolites for each reaction with updated reactions for enzyme catalyzed kinetic reactions:
{reaction ID: {metabolite ID with location tag: stoichiometry}}
rxn_catalysts: enzyme catalysts for each reaction with known catalysts, likely a subset of reactions in stoich with updated reactions for enzyme catalyzed kinetic reactions:
{reaction ID: enzyme IDs with location tag}
reversible_rxns: reaction IDs for reactions that have a reverse complement with updated reactions for enzyme catalyzed kinetic reactions, does not have reverse tag
rxn_id_to_compiled_id: mapping from reaction IDs to the IDs of the original reactions they were derived from, with updated reactions for enzyme catalyzed kinetic reactions
- Return type:
tuple[dict[str, Any], dict[str, dict[str, int]], dict[str, list[str]], list[str], dict[str, str]]
- _set_solver_values(constants)[source]
Sets values to be used in the FBA solver.
- Attributes set:
- Parameters:
constants (Constants)
- aa_supply_scaling(aa_conc, aa_present)[source]
Called during polypeptide_elongation process Determine amino acid supply rate scaling based on current amino acid concentrations.
- Parameters:
aa_conc (Unum) – internal concentration for each amino acid (ndarray[float])
aa_present (Unum) – whether each amino acid is in the external environment or not (ndarray[bool])
- Returns:
Scaling for the supply of each amino acid with higher supply rate if >1, lower supply rate if <1
- Return type:
- amino_acid_export(aa_transporters_counts, aa_conc, mechanistic_uptake)[source]
Calculate the rate of amino acid export.
- Parameters:
- Returns:
Rate of export for each amino acid (unitless but represents counts of amino acid per second)
- Return type:
- amino_acid_import(aa_in_media, dry_mass, internal_aa_conc, aa_transporters_counts, mechanisitic_uptake)[source]
Calculate the rate of amino acid uptake.
- Parameters:
aa_in_media (ndarray[Any, dtype[bool_]]) – bool for each amino acid being present in current media
dry_mass (Unum) – current dry mass of the cell, with mass units
internal_aa_conc (Unum | ndarray[Any, dtype[float64]]) – internal concentrations of amino acids
aa_transporters_counts (ndarray[Any, dtype[int64]]) – counts of each transporter
mechanisitic_uptake (bool) – if true, the uptake is calculated based on transporters
- Returns:
Rate of uptake for each amino acid (unitless but represents counts of amino acid per second)
- Return type:
- amino_acid_synthesis(counts_per_aa_fwd, counts_per_aa_rev, aa_conc)[source]
Calculate the net rate of synthesis for amino acid pathways (can be negative with reverse reactions).
- Parameters:
counts_per_aa_fwd (ndarray[Any, dtype[int64]]) – counts for enzymes in forward reactions for each amino acid
counts_per_aa_rev (ndarray[Any, dtype[int64]]) – counts for enzymes in loss reactions for each amino acid
aa_conc (Unum | ndarray[Any, dtype[float64]]) – concentrations of each amino acid with mol/volume units
- Returns:
3-element tuple containing
synthesis: net rate of synthesis for each amino acid pathway. array is unitless but represents counts of amino acid per second
forward_fraction: saturated fraction for forward reactions
loss_fraction: saturated fraction for loss reactions
- Return type:
tuple[ndarray[Any, dtype[float64]], ndarray[Any, dtype[float64]], ndarray[Any, dtype[float64]]]
Note
Currently does not match saturation terms used in calc_kcats since it assumes only a reverse or degradation KM exists for simpler calculations
- exchange_constraints(exchangeIDs, coefficient, targetUnits, media_id, unconstrained, constrained, concModificationsBasedOnCondition=None)[source]
Called during Metabolism process Returns the homeostatic objective concentrations based on the current nutrients Returns levels for external molecules available to exchange based on the current nutrients
- static extract_kinetic_constraints(raw_data, sim_data, stoich=None, catalysts=None, known_metabolites=None)[source]
Load and parse kinetic constraint information from raw_data
- Parameters:
raw_data (KnowledgeBaseEcoli) – knowledge base data
sim_data (SimulationDataEcoli) – simulation data
stoich (dict[str, dict[str, int]] | None) –
stoichiometry of metabolites for each reaction (if
None
, data is loaded fromraw_data
andsim_data
):{reaction ID: {metabolite ID with location tag: stoichiometry}}
catalysts (dict[str, list[str]] | None) –
enzyme catalysts for each reaction with known catalysts, likely a subset of reactions in
stoich
(ifNone
, data is loaded fromraw_data
andsim_data
:{reaction ID: enzyme IDs with location tag}
known_metabolites (Set[str] | None) – metabolites with known concentrations
- Returns:
Valid kinetic constraints for each reaction/enzyme pair:
{(reaction ID, enzyme with location tag): { 'kcat': kcat values (list[float]), 'saturation': saturation equations (list[str]) }}
- Return type:
- static extract_reactions(raw_data, sim_data)[source]
Extracts reaction data from raw_data to build metabolism reaction network with stoichiometry, reversibility and enzyme catalysts.
- Parameters:
raw_data (KnowledgeBaseEcoli) – knowledge base data
sim_data (SimulationDataEcoli) – simulation data
- Returns:
5-element tuple containing
base_rxn_ids: list of base reaction IDs from which reaction IDs were derived from
reaction_stoich: stoichiometry of metabolites for each reaction:
{reaction ID: {metabolite ID with location tag: stoichiometry}}
reversible_reactions: reaction IDs for reactions that have a reverse complement, does not have reverse tag
reaction_catalysts: enzyme catalysts for each reaction with known catalysts, likely a subset of reactions in stoich:
{reaction ID: enzyme IDs with location tag}
rxn_id_to_base_rxn_id: mapping from reaction IDs to the IDs of the base reactions they were derived from:
{reaction ID: base ID}
- Return type:
tuple[list[str], dict[str, dict[str, int]], list[str], dict[str, list[str]], dict[str, str]]
- get_aa_to_transporters_mapping_data(sim_data, export=False)[source]
Creates a dictionary that maps amino acids with their transporters. Based on this dictionary, it creates a correlation matrix with rows as AA and columns as transporters.
- Parameters:
sim_data (SimulationDataEcoli) – simulation data
export (bool) – if True, the parameters calculated are for mechanistic export instead of uptake
- Returns:
3-element tuple containing
aa_to_transporters: dictonary that maps aa to transporters involved in transport reactions
aa_to_transporters_matrix: correlation matrix. Columns correspond to transporter enzymes and rows to amino acids
aa_transporters_names: names of all transporters
- Return type:
tuple[dict[str, list], ndarray[Any, dtype[float64]], ndarray[Any, dtype[str_]]]
- get_kinetic_constraints(enzymes, substrates)[source]
Allows for dynamic code generation for kinetic constraint calculation for use in Metabolism process. Inputs should be unitless but the order of magnitude should match the kinetics parameters (umol/L/s).
If trying to pickle sim_data object after function has been called, _compiled_enzymes and _compiled_saturation might not be able to be pickled. See __getstate__(), __setstate__() comments on PR 111 to address.
Returns np.array of floats of the kinetic constraint target for each reaction with kinetic parameters
- Parameters:
enzymes (Unum) – concentrations of enzymes associated with kinetic constraints (mol / volume units)
substrates (Unum) – concentrations of substrates associated with kinetic constraints (mol / volume units)
- Returns:
Array of dimensions (n reactions, 3) where each row contains the min, mean and max kinetic constraints for each reaction with kinetic constraints (mol / volume / time units)
- Return type:
Unum
- get_pathway_enzyme_counts_per_aa(enzyme_counts)[source]
Get the counts of enzymes for forward and reverse reactions in the amino acid synthesis network based on all of the enzymes used in the network. Useful to get the counts to pass to amino_acid_synthesis() from counts based on self.aa_enzymes.
- Parameters:
enzyme_counts (ndarray) – counts of all enzymes in the amino acid network
- Returns:
- 2-element tuple containing
counts_per_aa_fwd: counts of enzymes for the forward reaction for each amino acid
counts_per_aa_rev: counts of enzymes for the reverse reaction for each amino acid
- Return type:
- static match_reaction(stoich, catalysts, rxn_to_match, enz, mets, direction=None)[source]
Matches a given reaction (rxn_to_match) to reactions that exist in stoich given that enz is known to catalyze the reaction and mets are reactants in the reaction. Can perform a fuzzy reaction match since rxn_to_match just needs to be part of the actual reaction name to match specific instances of a reaction. (eg. rxn_to_match=”ALCOHOL-DEHYDROG-GENERIC-RXN” can match “ALCOHOL-DEHYDROG-GENERIC-RXN-ETOH/NAD//ACETALD/NADH/PROTON.30.”).
- Parameters:
stoich (dict[str, dict[str, int]]) –
stoichiometry of metabolites for each reaction:
{reaction ID: {metabolite ID with location tag: stoichiometry}}
catalysts (dict[str, list[str]]) –
enzyme catalysts for each reaction with known catalysts, likely a subset of reactions in stoich:
{reaction ID: enzyme IDs with location tag}
rxn_to_match (str) – reaction ID from kinetics to match to existing reactions
enz (str) – enzyme ID with location tag
mets (list[str]) – metabolite IDs with no location tag from kinetics
direction (str | None) – reaction directionality,
'forward'
or'reverse'
orNone
- Returns:
Matched reaction IDs in stoich
- Return type:
- set_mechanistic_export_constants(sim_data, cell_specs, basal_container)[source]
Calls get_aa_to_transporters_mapping_data() for AA export, which calculates the total amount of export transporter counts per AA. Kcats are calculated using the same exchange rates as for uptake and transporter counts. Missing KMs are calculated based on present KMs. This is done by calculating the average factor for KMs compared to estimated concentration (av_factor = sum(KM / concentration) / n_aa_with_kms). ** KM = av_factor * concentration
- Parameters:
sim_data (SimulationDataEcoli) – simulation data
cell_specs (dict[str, dict]) – mapping from condition to calculated cell properties
basal_container (ndarray) – average initial bulk molecule counts in the basal condition (structured Numpy array, see Bulk Molecules)
- set_mechanistic_supply_constants(sim_data, cell_specs, basal_container, with_aa_container)[source]
Sets constants to determine amino acid supply during translation. Used with amino_acid_synthesis() and amino_acid_import() during simulations but supply can alternatively be determined phenomologically. This approach is more detailed and should better respond to environmental changes and perturbations but has more variability related to gene expression and regulation.
- Parameters:
sim_data (SimulationDataEcoli) – simulation data
cell_specs (dict[str, dict]) – mapping from condition to calculated cell properties
basal_container (ndarray) – average initial bulk molecule counts in the basal condition (structured Numpy array, see Bulk Molecules)
with_aa_container (ndarray) – average initial bulk molecule counts in the
with_aa
condition
- Sets class attributes:
Assumptions:
Only one reaction is limiting in an amino acid pathway (typically the first and one with KI) and the kcat for forward or reverse directions will apply to all enzymes that can catalyze that step
kcat for reverse and degradation reactions is the same (each amino acid only has reverse or degradation at this point but that could change with modifications to the amino_acid_pathways flat file)
- set_mechanistic_uptake_constants(sim_data, cell_specs, with_aa_container)[source]
Based on the matrix calculated in get_aa_to_transporters_mapping_data(), we calculate the total amount of transporter counts per AA.
- Parameters:
sim_data (SimulationDataEcoli) – simulation data
cell_specs (dict[str, dict]) – mapping from condition to calculated cell properties
with_aa_container (ndarray) – average initial bulk molecule counts in the
with_aa
condition (structured Numpy array, see Bulk Molecules)
- set_phenomological_supply_constants(sim_data)[source]
Sets constants to determine amino acid supply during translation. Used with aa_supply_scaling() during simulations but supply can alternatively be determined mechanistically. This approach may require manually adjusting constants (fraction_supply_inhibited and fraction_supply_exported) but has less variability related to gene expression and regulation.
- Parameters:
sim_data (SimulationDataEcoli) – simulation data
- Attributes set:
- Assumptions:
Each internal amino acid concentration in ‘minimal_plus_amino_acids’ media is not lower than in ‘minimal’ media
- TODO (Travis):
Better handling of concentration assumption
- reconstruction.ecoli.dataclasses.process.metabolism.amino_acid_export_jit(aa_transporters_counts, aa_conc, mechanistic_uptake, aa_to_exporters_matrix, export_kcats_per_aa, aa_export_kms)
- reconstruction.ecoli.dataclasses.process.metabolism.amino_acid_import_jit(aa_in_media, dry_mass, internal_aa_conc, aa_transporters_counts, mechanisitic_uptake, aa_import_kis, aa_to_importers_matrix, import_kcats_per_aa, max_specific_import_rates)
- reconstruction.ecoli.dataclasses.process.metabolism.amino_acid_synthesis_jit(counts_per_aa_fwd, counts_per_aa_rev, aa_conc, aa_upstream_kms, aa_kis, aa_reverse_kms, aa_degradation_kms, aa_forward_stoich, aa_kcats_fwd, aa_reverse_stoich, aa_kcats_rev)
- reconstruction.ecoli.dataclasses.process.metabolism.np_apply_along_axis(func1d, axis, arr)
- reconstruction.ecoli.dataclasses.process.metabolism.np_prod(array, axis)