Types

Base types

COBREXA._defaultMethod
_default(d, x::Union{Nothing, T} where T) -> Any

Fold the Maybe{T} down to T by defaulting.

source
COBREXA.AnnotationsType
Annotations = Dict{String,Vector{String}}

Dictionary used to store (possible multiple) standardized annotations of something, such as a Metabolite and a Reaction.

Example

Annotations("PubChem" => ["CID12345", "CID54321"])
source
COBREXA.GeneAssociationType
GeneAssociation = Vector{Vector{String}}

An association to genes, represented as a logical formula in a positive disjunctive normal form (DNF). (The 2nd-level vectors of strings are connected by "and" to form conjunctions, and the 1st-level vectors of these conjunctions are connected by "or" to form the DNF.)

source
COBREXA.MetabolicModelType
abstract type MetabolicModel end

A helper supertype of everything usable as a linear-like model for COBREXA functions.

If you want your model type to work with COBREXA, add the MetabolicModel as its supertype, and implement the accessor functions. Accessors reactions, metabolites, stoichiometry, bounds and objective must be implemented; others are not mandatory and default to safe "empty" values.

source
COBREXA.NotesType
Notes = Dict{String,Vector{String}}

Free-form notes about something (e.g. a Gene), categorized by "topic".

source

Model types and contents

COBREXA.CoreModelType
mutable struct CoreModel <: MetabolicModel

A "bare bones" core linear optimization problem of the form, with reaction and metabolite names.

min c^T x
s.t. S x = b
      xₗ ≤ x ≤ xᵤ

Fields

  • S::SparseArrays.SparseMatrixCSC{Float64, Int64}

  • b::SparseArrays.SparseVector{Float64, Int64}

  • c::SparseArrays.SparseVector{Float64, Int64}

  • xl::Vector{Float64}

  • xu::Vector{Float64}

  • rxns::Vector{String}

  • mets::Vector{String}

  • grrs::Vector{Union{Nothing, Vector{Vector{String}}}}

source
Base.convertMethod
convert(_::Type{CoreModel}, m::MetabolicModel) -> MetabolicModel

Make a CoreModel out of any compatible model type.

source
COBREXA.balanceMethod
balance(a::CoreModel) -> SparseArrays.SparseVector{Float64, Int64}

CoreModel target flux balance.

source
COBREXA.boundsMethod
bounds(a::CoreModel) -> Tuple{Vector{Float64}, Vector{Float64}}

CoreModel flux bounds.

source
COBREXA.genesMethod
genes(a::CoreModel) -> Vector{String}

Collect all genes contained in the CoreModel. The call is expensive for large models, because the vector is not stored and instead gets rebuilt each time this function is called.

source
COBREXA.objectiveMethod
objective(a::CoreModel) -> SparseArrays.SparseVector{Float64, Int64}

CoreModel objective vector.

source
COBREXA.stoichiometryMethod
stoichiometry(a::CoreModel) -> SparseArrays.SparseMatrixCSC{Float64, Int64}

CoreModel stoichiometry matrix.

source
COBREXA.CoreCouplingType
mutable struct CoreCoupling{M} <: ModelWrapper

A matrix-based wrap that adds reaction coupling matrix to the inner model. A flux x feasible in this model must satisfy:

    cₗ ≤ C x ≤ cᵤ

Fields

  • lm::Any

  • C::SparseArrays.SparseMatrixCSC{Float64, Int64}

  • cl::Vector{Float64}

  • cu::Vector{Float64}

source
COBREXA.CoreModelCoupledType
const CoreModelCoupled = CoreCoupling{CoreModel}

A matrix-based linear model with additional coupling constraints in the form:

    cₗ ≤ C x ≤ cᵤ

Internally, the model is implemented using CoreCoupling that contains a single CoreModel.

source
Base.convertMethod
convert(::Type{CoreCoupling{M}}, mm::MetabolicModel; clone_coupling) -> MetabolicModel

Make a CoreCoupling out of any compatible model type.

source
COBREXA.couplingMethod
coupling(a::CoreCoupling) -> SparseArrays.SparseMatrixCSC{Float64, Int64}

Coupling constraint matrix for a CoreCoupling.

source
COBREXA.coupling_boundsMethod
coupling_bounds(a::CoreCoupling) -> Tuple{Vector{Float64}, Vector{Float64}}

Coupling bounds for a CoreCoupling.

source
COBREXA.FluxSummaryType
struct FluxSummary

A struct used to store summary information about the solution of a constraint based analysis result.

Fields

  • biomass_fluxes::OrderedCollections.OrderedDict{String, Float64}

  • import_fluxes::OrderedCollections.OrderedDict{String, Float64}

  • export_fluxes::OrderedCollections.OrderedDict{String, Float64}

  • unbounded_fluxes::OrderedCollections.OrderedDict{String, Float64}

source
COBREXA.flux_summaryMethod
flux_summary(flux_result::Union{Nothing, Dict{String, Float64}}; exclude_exchanges, exchange_prefixes, biomass_strings, exclude_biomass, small_flux_bound, large_flux_bound, keep_unbounded) -> FluxSummary

Summarize a dictionary of fluxes into small, useful representation of the most important information contained. Useful for pretty-printing and quickly exploring the results. Internally this function uses looks_like_biomass_reaction and looks_like_exchange_reaction. The corresponding keyword arguments passed to these functions. Use this if your model has non-standard ids for reactions. Fluxes smaller than small_flux_bound are not stored, while fluxes larger than large_flux_bound are only stored if keep_unbounded is true.

Example

julia> sol = flux_dict(flux_balance_analysis(model, Tulip.Optimizer))
julia> fr = flux_summary(sol)
Biomass:
  BIOMASS_Ecoli_core_w_GAM: 0.8739
Import:
  EX_o2_e:     -21.7995
  EX_glc__D_e: -10.0
  EX_nh4_e:    -4.7653
  EX_pi_e:     -3.2149
Export:
  EX_h_e:      17.5309
  EX_co2_e:    22.8098
  EX_h2o_e:    29.1758
source
COBREXA.FluxVariabilitySummaryType
struct FluxVariabilitySummary

Stores summary information about the result of a flux variability analysis.

Fields

  • biomass_fluxes::Dict{String, Vector{Union{Nothing, Float64}}}

  • exchange_fluxes::Dict{String, Vector{Union{Nothing, Float64}}}

source
COBREXA.flux_variability_summaryMethod
flux_variability_summary(flux_result::Tuple{Dict{String, Dict{String, Float64}}, Dict{String, Dict{String, Float64}}}; exclude_exchanges, exchange_prefixes, biomass_strings, exclude_biomass) -> FluxVariabilitySummary

Summarize a dictionary of flux dictionaries obtained eg. from flux_variability_analysis_dict. The simplified summary representation is useful for pretty-printing and easily showing the most important results.

Internally this function uses looks_like_biomass_reaction and looks_like_exchange_reaction. The corresponding keyword arguments are passed to these functions. Use this if your model has an uncommon naming of reactions.

Example

julia> sol = flux_variability_analysis_dict(model, Gurobi.Optimizer; bounds = objective_bounds(0.99))
julia> flux_res = flux_variability_summary(sol)
Biomass                     Lower bound   Upper bound
  BIOMASS_Ecoli_core_w_GAM: 0.8652        0.8652
Exchange
  EX_h2o_e:                 28.34         28.34
  EX_co2_e:                 22.0377       22.0377
  EX_o2_e:                  -22.1815      -22.1815
  EX_h_e:                   17.3556       17.3556
  EX_glc__D_e:              -10.0         -10.0
  EX_nh4_e:                 -4.8448       -4.8448
  EX_pi_e:                  -3.2149       -3.2149
  EX_for_e:                 0.0           0.0
  ...                       ...           ...
source
COBREXA.GeneType
Gene() -> Gene
Gene(id; name, notes, annotations) -> Gene

A convenient constructor for a Gene.

source
COBREXA.GeneType
mutable struct Gene

Fields

  • id::String

  • name::Union{Nothing, String}

  • notes::Dict{String, Vector{String}}

  • annotations::Dict{String, Vector{String}}

source
COBREXA.HDF5ModelType
mutable struct HDF5Model <: MetabolicModel

A model that is stored in HDF5 format. The model data is never really pulled into memory, but instead mmap'ed as directly as possible into the Julia structures. This makes reading the HDF5Models extremely fast, at the same time the (uncached) HDF5Models can be sent around efficiently among distributed nodes just like Serialized models, provided the nodes share a common storage.

All HDF5Models must have the backing disk storage. To create one, use save_h5_model or save_model with .h5 file extension. To create a temporary model that behaves like a model "in memory", save it to a temporary file. For related reasons, you can not use convert models to HDF5Model format, because the conversion would impliy having the model saved somewhere.

Fields

  • h5::Union{Nothing, HDF5.File}

  • filename::String

source
COBREXA.IsozymeType
mutable struct Isozyme

Information about isozyme composition and activity.

Fields

  • gene_product_count::Dict{String, Int64}

  • kcat_forward::Float64

  • kcat_reverse::Float64

source
COBREXA.JSONModelType
struct JSONModel <: MetabolicModel

A struct used to store the contents of a JSON model, i.e. a model read from a file ending with .json. These model files typically store all the model data in arrays of JSON objects (represented in Julia as vectors of dictionaries).

Usually, not all of the fields of the input JSON can be easily represented when converting to other models, care should be taken to avoid losing information.

Direct work with the json structure is not very efficient; the model structure therefore caches some of the internal structure in the extra fields. The single-parameter JSONModel constructor creates these caches correctly from the json. The model structure is designed as read-only, and changes in json invalidate the cache.

Example

model = load_json_model("some_model.json")
model.json # see the actual underlying JSON
reactions(model) # see the list of reactions

Fields

  • json::Dict{String, Any}

  • rxn_index::Dict{String, Int64}

  • rxns::Vector{Any}

  • met_index::Dict{String, Int64}

  • mets::Vector{Any}

  • gene_index::Dict{String, Int64}

  • genes::Vector{Any}

source
COBREXA.boundsMethod
bounds(model::JSONModel) -> Tuple{Vector, Vector}

Get the bounds for reactions, assuming the information is stored in .lower_bound and .upper_bound.

source
COBREXA.gene_nameMethod
gene_name(model::JSONModel, gid::String) -> Any

Return the name of gene with ID gid.

source
COBREXA.genesMethod
genes(model::JSONModel) -> Vector

Extract gene names from a JSON model.

source
COBREXA.metabolite_formulaMethod
metabolite_formula(model::JSONModel, mid::String) -> Union{Nothing, Dict{String, Int64}}

Parse and return the metabolite .formula

source
COBREXA.metabolitesMethod
metabolites(model::JSONModel) -> Vector

Extract metabolite names (stored as .id) from JSON model.

source
COBREXA.objectiveMethod
objective(model::JSONModel) -> SparseArrays.SparseVector{_A, Int64} where _A

Collect .objective_coefficient keys from model reactions.

source
COBREXA.reaction_gene_associationMethod
reaction_gene_association(model::JSONModel, rid::String) -> Union{Nothing, Vector{Vector{String}}}

Parses the .gene_reaction_rule from reactions.

source
COBREXA.reactionsMethod
reactions(model::JSONModel) -> Vector

Extract reaction names (stored as .id) from JSON model.

source
COBREXA.stoichiometryMethod
stoichiometry(model::JSONModel) -> SparseArrays.SparseMatrixCSC{Float64, Int64}

Get the stoichiometry. Assuming the information is stored in reaction object under key .metabolites.

source
COBREXA.MATModelType
struct MATModel <: MetabolicModel

Wrapper around the models loaded in dictionaries from the MATLAB representation.

Fields

  • mat::Dict{String, Any}
source
Base.convertMethod
convert(_::Type{MATModel}, m::MetabolicModel) -> MetabolicModel

Convert any metabolic model to MATModel.

source
COBREXA.balanceMethod
balance(m::MATModel) -> Any

Extracts balance from the MAT model, defaulting to zeroes if not present.

source
COBREXA.boundsMethod
bounds(m::MATModel) -> Tuple{Any, Any}

Extracts bounds from the MAT file, saved under lb and ub.

source
COBREXA.coupling_boundsMethod
coupling_bounds(m::MATModel) -> Tuple{Any, Any}

Extracts the coupling constraints. Currently, there are several accepted ways to store these in MATLAB models; this takes the constraints from vectors cl and cu.

source
COBREXA.genesMethod
genes(m::MATModel) -> Any

Extracts the possible gene list from genes key.

source
COBREXA.metabolite_chargeMethod
metabolite_charge(m::MATModel, mid::String) -> Union{Nothing, Int64}

Extract metabolite charge from metCharge or metCharges.

source
COBREXA.metabolite_formulaMethod
metabolite_formula(m::MATModel, mid::String) -> Union{Nothing, Dict{String, Int64}}

Extract metabolite formula from key metFormula or metFormulas.

source
COBREXA.metabolitesMethod
metabolites(m::MATModel) -> Vector{String}

Extracts metabolite names from mets key in the MAT file.

source
COBREXA.objectiveMethod
objective(m::MATModel) -> Any

Extracts the objective from the MAT model (defaults to zeroes).

source
COBREXA.reaction_gene_associationMethod
reaction_gene_association(m::MATModel, rid::String) -> Union{Nothing, Vector{Vector{String}}}

Extracts the associations from grRules key, if present.

source
COBREXA.reactionsMethod
reactions(m::MATModel) -> Vector{String}

Extracts reaction names from rxns key in the MAT file.

source
COBREXA.balanceMethod
balance(a::MetabolicModel) -> SparseArrays.SparseVector{Float64, Int64}

Get the sparse balance vector of a model.

source
COBREXA.boundsMethod
bounds(a::MetabolicModel) -> Tuple{Vector{Float64}, Vector{Float64}}

Get the lower and upper solution bounds of a model.

source
COBREXA.couplingMethod
coupling(a::MetabolicModel) -> SparseArrays.SparseMatrixCSC{Float64, Int64}

Get a matrix of coupling constraint definitions of a model. By default, there is no coupling in the models.

source
COBREXA.coupling_boundsMethod
coupling_bounds(a::MetabolicModel) -> Tuple{Vector{Float64}, Vector{Float64}}

Get the lower and upper bounds for each coupling bound in a model, as specified by coupling. By default, the model does not have any coupling bounds.

source
COBREXA.fluxesMethod
fluxes(a::MetabolicModel) -> Vector{String}

In some models, the reactions that correspond to the columns of stoichiometry matrix do not fully represent the semantic contents of the model; for example, fluxes may be split into forward and reverse reactions, reactions catalyzed by distinct enzymes, etc. Together with reaction_flux (and n_fluxes) this specifies how the flux is decomposed into individual reactions.

By default (and in most models), fluxes and reactions perfectly correspond.

source
COBREXA.gene_annotationsMethod
gene_annotations(a::MetabolicModel, gene_id::String) -> Any

Return standardized names that identify the corresponding gene or product. The dictionary assigns vectors of possible identifiers to identifier system names, e.g. "PDB" => ["PROT01"].

source
COBREXA.gene_nameMethod
gene_name(model::MetabolicModel, gid::String) -> Any

Return the name of gene with ID gid.

source
COBREXA.gene_notesMethod
gene_notes(model::MetabolicModel, gene_id::String) -> Any

Return the notes associated with the gene gene_id in model.

source
COBREXA.genesMethod
genes(a::MetabolicModel) -> Vector{String}

Return identifiers of all genes contained in the model. By default, there are no genes.

In SBML, these are usually called "gene products" but we write genes for simplicity.

source
COBREXA.metabolite_annotationsMethod
metabolite_annotations(a::MetabolicModel, metabolite_id::String) -> Any

Return standardized names that may help to reliably identify the metabolite. The dictionary assigns vectors of possible identifiers to identifier system names, e.g. "ChEMBL" => ["123"] or "PubChem" => ["CID123", "CID654645645"].

source
COBREXA.metabolite_chargeMethod
metabolite_charge(model::MetabolicModel, metabolite_id::String) -> Any

Return the charge associated with metabolite metabolite_id in model. Returns nothing if charge not present.

source
COBREXA.metabolite_compartmentMethod
metabolite_compartment(model::MetabolicModel, metabolite_id::String) -> Any

Return the compartment of metabolite metabolite_id in model if it is assigned. If not, return nothing.

source
COBREXA.metabolite_formulaMethod
metabolite_formula(model::MetabolicModel, metabolite_id::String) -> Any

Return the formula of metabolite metabolite_id in model. Return nothing in case the formula is not known or irrelevant.

source
COBREXA.metabolite_notesMethod
metabolite_notes(model::MetabolicModel, metabolite_id::String) -> Any

Return the notes associated with metabolite reaction_id in model.

source
COBREXA.metabolitesMethod
metabolites(a::MetabolicModel) -> Vector{String}

Return a vector of metabolite identifiers in a model. The vector precisely corresponds to the rows in stoichiometry matrix.

As with reactionss, some metabolites in models may be virtual, representing purely technical equality constraints.

source
COBREXA.n_genesMethod
n_genes(a::MetabolicModel) -> Int64

Return the number of genes in the model (as returned by genes). If you just need the number of the genes, this may be much more efficient than calling genes and measuring the array.

source
COBREXA.objectiveMethod
objective(a::MetabolicModel) -> SparseArrays.SparseVector{Float64, Int64}

Get the objective vector of the model. Analysis functions, such as flux_balance_analysis, are supposed to maximize dot(objective, x) where x is a feasible solution of the model.

source
COBREXA.precache!Method
precache!(a::MetabolicModel)

Do whatever is feasible to get the model into a state that can be read from as-quickly-as-possible. This may include e.g. generating helper index structures and loading delayed parts of the model from disk. The model should be modified "transparently" in-place. Analysis functions call this right before applying modifications or converting the model to the optimization model using make_optimization_model; usually on the same machine where the optimizers (and, generally, the core analysis algorithms) will run. The calls are done in a good hope that the performance will be improved.

By default, it should be safe to do nothing.

source
COBREXA.reaction_annotationsMethod
reaction_annotations(a::MetabolicModel, reaction_id::String) -> Any

Return standardized names that may help identifying the reaction. The dictionary assigns vectors of possible identifiers to identifier system names, e.g. "Reactome" => ["reactomeID123"].

source
COBREXA.reaction_fluxMethod
reaction_flux(a::MetabolicModel) -> Any

Retrieve a sparse matrix that describes the correspondence of a solution of the linear system to the fluxes (see fluxes for rationale). Returns a sparse matrix of size (n_reactions(a), n_fluxes(a)). For most models, this is an identity matrix.

source
COBREXA.reaction_gene_associationMethod
reaction_gene_association(a::MetabolicModel, reaction_id::String) -> Any

Returns the sets of genes that need to be present so that the reaction can work (technically, a DNF on gene availability, with positive atoms only).

For simplicity, nothing may be returned, meaning that the reaction always takes place. (in DNF, that would be equivalent to returning [[]].)

source
COBREXA.reaction_nameMethod
reaction_name(model::MetabolicModel, rid::String) -> Any

Return the name of reaction with ID rid.

source
COBREXA.reaction_notesMethod
reaction_notes(model::MetabolicModel, reaction_id::String) -> Any

Return the notes associated with reaction reaction_id in model.

source
COBREXA.reaction_stoichiometryMethod
reaction_stoichiometry(m::MetabolicModel, rid::String) -> Any

Return the stoichiometry of reaction with ID rid in the model. The dictionary maps the metabolite IDs to their stoichiometric coefficients.

source
COBREXA.reaction_subsystemMethod
reaction_subsystem(model::MetabolicModel, reaction_id::String) -> Any

Return the subsystem of reaction reaction_id in model if it is assigned. If not, return nothing.

source
COBREXA.reactionsMethod
reactions(a::MetabolicModel) -> Vector{String}

Return a vector of reaction identifiers in a model. The vector precisely corresponds to the columns in stoichiometry matrix.

For technical reasons, the "reactions" may sometimes not be true reactions but various virtual and helper pseudo-reactions that are used in the metabolic modeling, such as metabolite exchanges, separate forward and reverse reactions, supplies of enzymatic and genetic material and virtual cell volume, etc. To simplify the view of the model contents use reaction_flux.

source
COBREXA.stoichiometryMethod
stoichiometry(a::MetabolicModel) -> SparseArrays.SparseMatrixCSC{Float64, Int64}

Get the sparse stoichiometry matrix of a model. A feasible solution x of a model m is defined as satisfying the equations:

  • stoichiometry(m) * x .== balance(m)
  • x .>= lbs
  • y .<= ubs
  • `(lbs, ubs) == bounds(m)
source
COBREXA.MetaboliteType
Metabolite() -> Metabolite
Metabolite(id; name, formula, charge, compartment, notes, annotations) -> Metabolite

A constructor for Metabolites.

source
COBREXA.MetaboliteType
mutable struct Metabolite

Fields

  • id::String

  • name::Union{Nothing, String}

  • formula::Union{Nothing, String}

  • charge::Union{Nothing, Int64}

  • compartment::Union{Nothing, String}

  • notes::Dict{String, Vector{String}}

  • annotations::Dict{String, Vector{String}}

source
COBREXA.ReactionType
Reaction(id::String, metabolites) -> Reaction
Reaction(id::String, metabolites, dir; default_bound) -> Reaction

Convenience constructor for Reaction. The reaction equation is specified using metabolites, which is a dictionary mapping metabolite ids to stoichiometric coefficients. The direcion of the reaction is set through dir which can take :bidirectional, :forward, and :reverse as values. Finally, the default_bound is the value taken to mean infinity in the context of constraint based models, often this is set to a very high flux value like 1000.

source
COBREXA.ReactionType
Reaction() -> Reaction
Reaction(id; name, metabolites, lb, ub, grr, subsystem, notes, annotations, objective_coefficient) -> Reaction

A constructor for Reaction that only takes a reaction id and assigns default/uninformative values to all the fields that are not explicitely assigned.

source
COBREXA.ReactionType
mutable struct Reaction

A structure for representing a single reaction in a StandardModel.

Fields

  • id::String

  • name::Union{Nothing, String}

  • metabolites::Dict{String, Float64}

  • lb::Float64

  • ub::Float64

  • grr::Union{Nothing, Vector{Vector{String}}}

  • subsystem::Union{Nothing, String}

  • notes::Dict{String, Vector{String}}

  • annotations::Dict{String, Vector{String}}

  • objective_coefficient::Float64

source
COBREXA.ReactionStatusType
mutable struct ReactionStatus

Used for concise reporting of modeling results.

Fields

  • already_present::Bool

  • index::Int64

  • info::String

source
COBREXA.SBMLModelType
struct SBMLModel <: MetabolicModel

Construct the SBML model and add the necessary cached indexes, possibly choosing an active objective.

source
COBREXA.SBMLModelType
struct SBMLModel <: MetabolicModel

Thin wrapper around the model from SBML.jl library. Allows easy conversion from SBML to any other model format.

Fields

  • sbml::SBML.Model

  • reaction_ids::Vector{String}

  • reaction_idx::Dict{String, Int64}

  • metabolite_ids::Vector{String}

  • metabolite_idx::Dict{String, Int64}

  • gene_ids::Vector{String}

  • active_objective::String

source
Base.convertMethod
convert(_::Type{SBMLModel}, mm::MetabolicModel) -> MetabolicModel

Convert any metabolic model to SBMLModel.

source
COBREXA.balanceMethod
balance(model::SBMLModel) -> SparseArrays.SparseVector{Float64, Int64}

Balance vector of a SBMLModel. This is always zero.

source
COBREXA.boundsMethod
bounds(model::SBMLModel) -> Tuple{Vector{Float64}, Vector{Float64}}

Get the lower and upper flux bounds of model SBMLModel. Throws DomainError in case if the SBML contains mismatching units.

source
COBREXA.gene_annotationsMethod
gene_annotations(model::SBMLModel, gid::String) -> Dict{String, Vector{String}}

Return the annotations of gene with ID gid.

source
COBREXA.gene_nameMethod
gene_name(model::SBMLModel, gid::String) -> Union{Nothing, String}

Return the name of gene with ID gid.

source
COBREXA.gene_notesMethod
gene_notes(model::SBMLModel, gid::String) -> Dict{String, Vector{String}}

Return the notes about gene with ID gid.

source
COBREXA.metabolite_annotationsMethod
metabolite_annotations(model::SBMLModel, mid::String) -> Dict{String, Vector{String}}

Return the annotations of metabolite with ID mid.

source
COBREXA.metabolite_nameMethod
metabolite_name(model::SBMLModel, mid::String) -> Union{Nothing, String}

Return the name of metabolite with ID mid.

source
COBREXA.metabolite_notesMethod
metabolite_notes(model::SBMLModel, mid::String) -> Dict{String, Vector{String}}

Return the notes about metabolite with ID mid.

source
COBREXA.reaction_annotationsMethod
reaction_annotations(model::SBMLModel, rid::String) -> Dict{String, Vector{String}}

Return the annotations of reaction with ID rid.

source
COBREXA.reaction_nameMethod
reaction_name(model::SBMLModel, rid::String) -> Union{Nothing, String}

Return the name of reaction with ID rid.

source
COBREXA.reaction_notesMethod
reaction_notes(model::SBMLModel, rid::String) -> Dict{String, Vector{String}}

Return the notes about reaction with ID rid.

source
COBREXA.SerializedType
mutable struct Serialized{M} <: ModelWrapper

A meta-model that represents a model that is serialized on the disk. The internal model will be loaded on-demand by using any accessor, or by calling precache! directly.

Fields

  • m::Union{Nothing, M} where M

  • filename::String

source
COBREXA.precache!Method
precache!(model::Serialized)

Load the Serialized model from disk in case it's not alreadly loaded.

source
COBREXA.unwrap_modelMethod
unwrap_model(m::Serialized) -> Union{Nothing, M} where M

Unwrap the serialized model (precaching it transparently).

source
COBREXA.StandardModelType
mutable struct StandardModel <: MetabolicModel

StandardModel is used to store a constraint based metabolic model with meta-information. Meta-information is defined as annotation details, which include gene-reaction-rules, formulas, etc.

This model type seeks to keep as much meta-information as possible, as opposed to CoreModel and CoreModelCoupled, which keep the bare neccessities only. When merging models and keeping meta-information is important, use this as the model type. If meta-information is not important, use the more efficient core model types. See CoreModel and CoreModelCoupled for comparison.

In this model, reactions, metabolites, and genes are stored in ordered dictionaries indexed by each struct's id field. For example, model.reactions["rxn1_id"] returns a Reaction where the field id equals "rxn1_id". This makes adding and removing reactions efficient.

Note that the stoichiometric matrix (or any other core data, e.g. flux bounds) is not stored directly as in CoreModel. When this model type is used in analysis functions, these core data structures are built from scratch each time an analysis function is called. This can cause performance issues if you run many small analysis functions sequentially. Consider using the core model types if performance is critical.

See also: Reaction, Metabolite, Gene

Example

model = load_model(StandardModel, "my_model.json")
keys(model.reactions)

Fields

  • id::String

  • reactions::OrderedCollections.OrderedDict{String, Reaction}

  • metabolites::OrderedCollections.OrderedDict{String, Metabolite}

  • genes::OrderedCollections.OrderedDict{String, Gene}

source
Base.convertMethod
convert(_::Type{StandardModel}, model::MetabolicModel) -> MetabolicModel

Convert any MetabolicModel into a StandardModel. Note, some data loss may occur since only the generic interface is used during the conversion process.

source
COBREXA.balanceMethod
balance(model::StandardModel) -> SparseArrays.SparseVector{Float64, Int64}

Return the balance of the linear problem, i.e. b in Sv = 0 where S is the stoichiometric matrix and v is the flux vector.

source
COBREXA.boundsMethod
bounds(model::StandardModel) -> Tuple{Vector{Float64}, Vector{Float64}}

Return the lower and upper bounds, respectively, for reactions in model. Order matches that of the reaction ids returned in reactions().

source
COBREXA.gene_annotationsMethod
gene_annotations(model::StandardModel, id::String) -> Dict{String, Vector{String}}

Return the annotation associated with gene id in model. Return an empty Dict if not present.

source
COBREXA.gene_nameMethod
gene_name(m::StandardModel, gid::String) -> Union{Nothing, String}

Return the name of gene with ID id.

source
COBREXA.gene_notesMethod
gene_notes(model::StandardModel, id::String) -> Dict{String, Vector{String}}

Return the notes associated with gene id in model. Return an empty Dict if not present.

source
COBREXA.genesMethod
genes(model::StandardModel) -> Vector{String}

Return a vector of gene id strings in model.

source
COBREXA.lower_boundsMethod
lower_bounds(model::StandardModel) -> Vector{Float64}

Return the lower bounds for all reactions in model in sparse format.

source
COBREXA.metabolite_annotationsMethod
metabolite_annotations(model::StandardModel, id::String) -> Dict{String, Vector{String}}

Return the annotation associated with metabolite id in model. Return an empty Dict if not present.

source
COBREXA.metabolite_chargeMethod
metabolite_charge(model::StandardModel, id::String) -> Union{Nothing, Int64}

Return the charge associated with metabolite id in model. Return nothing if not present.

source
COBREXA.metabolite_compartmentMethod
metabolite_compartment(model::StandardModel, id::String) -> Union{Nothing, String}

Return compartment associated with metabolite id in model. Return nothing if not present.

source
COBREXA.metabolite_formulaMethod
metabolite_formula(model::StandardModel, id::String) -> Union{Nothing, Dict{String, Int64}}

Return the formula of reaction id in model. Return nothing if not present.

source
COBREXA.metabolite_nameMethod
metabolite_name(m::StandardModel, mid::String) -> Union{Nothing, String}

Return the name of metabolite with ID id.

source
COBREXA.metabolite_notesMethod
metabolite_notes(model::StandardModel, id::String) -> Dict{String, Vector{String}}

Return the notes associated with metabolite id in model. Return an empty Dict if not present.

source
COBREXA.metabolitesMethod
metabolites(model::StandardModel) -> Vector{String}

Return a vector of metabolite id strings contained in model. The order of metabolite strings returned here matches the order used to construct the stoichiometric matrix.

source
COBREXA.n_genesMethod
n_genes(model::StandardModel) -> Int64

Return the number of genes in model.

source
COBREXA.n_reactionsMethod
n_reactions(model::StandardModel) -> Int64

Return the number of reactions contained in model.

source
COBREXA.objectiveMethod
objective(model::StandardModel) -> SparseArrays.SparseVector{Float64, Int64}

Return sparse objective vector for model.

source
COBREXA.reaction_annotationsMethod
reaction_annotations(model::StandardModel, id::String) -> Dict{String, Vector{String}}

Return the annotation associated with reaction id in model. Return an empty Dict if not present.

source
COBREXA.reaction_gene_associationMethod
reaction_gene_association(model::StandardModel, id::String) -> Union{Nothing, Vector{Vector{String}}}

Return the gene reaction rule in string format for reaction with id in model. Return nothing if not available.

source
COBREXA.reaction_nameMethod
reaction_name(m::StandardModel, rid::String) -> Union{Nothing, String}

Return the name of reaction with ID id.

source
COBREXA.reaction_notesMethod
reaction_notes(model::StandardModel, id::String) -> Dict{String, Vector{String}}

Return the notes associated with reaction id in model. Return an empty Dict if not present.

source
COBREXA.reaction_subsystemMethod
reaction_subsystem(model::StandardModel, id::String) -> Union{Nothing, String}

Return the subsystem associated with reaction id in model. Return nothing if not present.

source
COBREXA.reactionsMethod
reactions(model::StandardModel) -> Vector{String}

Return a vector of reaction id strings contained in model. The order of reaction ids returned here matches the order used to construct the stoichiometric matrix.

source
COBREXA.stoichiometryMethod
stoichiometry(model::StandardModel) -> SparseArrays.SparseMatrixCSC{Float64, Int64}

Return the stoichiometric matrix associated with model in sparse format.

source
COBREXA.upper_boundsMethod
upper_bounds(model::StandardModel) -> Vector{Float64}

Return the upper bounds for all reactions in model in sparse format. Order matches that of the reaction ids returned in reactions().

source

Model type wrappers

COBREXA.ExpressionLimitedModelType
struct ExpressionLimitedModel <: ModelWrapper

ExpressionLimitedModel follows the methodology of the E-Flux algorithm to constraint the flux through the reactions in order to simulate the limited expression of genes (and thus the limited availability of the gene products).

Use make_expression_limited_model or with_expression_limits to construct the models.

E-Flux algorithm is closer described by: Colijn, Caroline, Aaron Brandes, Jeremy Zucker, Desmond S. Lun, Brian Weiner, Maha R. Farhat, Tan-Yun Cheng, D. Branch Moody, Megan Murray, and James E. Galagan. "Interpreting expression data with metabolic flux models: predicting Mycobacterium tuberculosis mycolic acid production." PLoS computational biology 5, no. 8 (2009): e1000489.

Fields

  • relative_expression::Dict{String, Float64}

    Relative gene expression w.r.t. to some chosen reference; the normalization and scale of the values should match the expectations of the bounding_function.

source
COBREXA.GeckoModelType
struct GeckoModel <: ModelWrapper

A model with complex enzyme concentration and capacity bounds, as described in Sánchez, Benjamín J., et al. "Improving the phenotype predictions of a yeast genome-scale metabolic model by incorporating enzymatic constraints." Molecular systems biology 13.8 (2017): 935.

Use make_gecko_model or with_gecko to construct this kind of model.

The model wraps another "internal" model, and adds following modifications:

  • enzymatic reactions with known enzyme information are split into multiple forward and reverse variants for each isozyme,
  • reaction coupling is added to ensure the groups of isozyme reactions obey the global reaction flux bounds from the original model,
  • gene concentrations specified by each reaction and its gene product stoichiometry, can constrained by the user to reflect measurements, such as from mass spectrometry,
  • additional coupling is added to simulate total masses of different proteins grouped by type (e.g., membrane-bound and free-floating proteins), which can be again constrained by the user (this is slightly generalized from original GECKO algorithm, which only considers a single group of indiscernible proteins).

The structure contains fields columns that describe the contents of the stoichiometry matrix columns, coupling_row_reaction, coupling_row_gene_product and coupling_row_mass_group that describe correspondence of the coupling rows to original model and determine the coupling bounds (note: the coupling for gene product is actually added to stoichiometry, not in coupling), and inner, which is the original wrapped model. The objective of the model includes also the extra columns for individual genes, as held by coupling_row_gene_product.

Implementation exposes the split reactions (available as reactions(model)), but retains the original "simple" reactions accessible by fluxes. The related constraints are implemented using coupling and coupling_bounds.

Fields

  • objective::SparseArrays.SparseVector{Float64, Int64}

  • columns::Vector{COBREXA._gecko_reaction_column}

  • coupling_row_reaction::Vector{Int64}

  • coupling_row_gene_product::Vector{Tuple{Int64, Tuple{Float64, Float64}}}

  • coupling_row_mass_group::Vector{COBREXA._gecko_capacity}

  • inner::MetabolicModel

source
COBREXA._gecko_capacityType
struct _gecko_capacity

A helper struct that contains the gene product capacity terms organized by the grouping type, e.g. metabolic or membrane groups etc.

Fields

  • group_id::String

  • gene_product_idxs::Vector{Int64}

  • gene_product_molar_masses::Vector{Float64}

  • group_upper_bound::Float64

source
COBREXA._gecko_reaction_columnType
struct _gecko_reaction_column

A helper type for describing the contents of GeckoModels.

Fields

  • reaction_idx::Int64

  • isozyme_idx::Int64

  • direction::Int64

  • reaction_coupling_row::Int64

  • lb::Float64

  • ub::Float64

  • gene_product_coupling::Vector{Tuple{Int64, Float64}}

source
COBREXA.balanceMethod
balance(model::GeckoModel) -> Any

Return the balance of the reactions in the inner model, concatenated with a vector of zeros representing the enzyme balance of a GeckoModel.

source
COBREXA.couplingMethod
coupling(model::GeckoModel) -> Any

Return the coupling of GeckoModel. That combines the coupling of the wrapped model, coupling for split (arm) reactions, and the coupling for the total enzyme capacity.

source
COBREXA.genesMethod
genes(model::GeckoModel) -> Any

Return the gene ids of genes that have enzymatic constraints associated with them.

source
COBREXA.n_genesMethod
n_genes(model::GeckoModel) -> Int64

Return the number of genes that have enzymatic constraints associated with them.

source
COBREXA.n_reactionsMethod
n_reactions(model::GeckoModel) -> Int64

Returns the number of all irreversible reactions in model as well as the number of gene products that take part in enzymatic reactions.

source
COBREXA.objectiveMethod
objective(model::GeckoModel) -> SparseArrays.SparseVector{Float64, Int64}

Return the objective of the GeckoModel. Note, the objective is with respect to the internal variables, i.e. reactions(model), which are the unidirectional reactions and the genes involved in enzymatic reactions that have kinetic data.

source
COBREXA.reactionsMethod
reactions(model::GeckoModel) -> Any

Returns the internal reactions in a GeckoModel (these may be split to forward- and reverse-only parts with different isozyme indexes; reactions IDs are mangled accordingly with suffixes).

source
COBREXA.stoichiometryMethod
stoichiometry(model::GeckoModel) -> Any

Return a stoichiometry of the GeckoModel. The enzymatic reactions are split into unidirectional forward and reverse ones, each of which may have multiple variants per isozyme.

source
COBREXA.SMomentModelType
struct SMomentModel <: ModelWrapper

An enzyme-capacity-constrained model using sMOMENT algorithm, as described by Bekiaris, Pavlos Stephanos, and Steffen Klamt, "Automatic construction of metabolic models with enzyme constraints" BMC bioinformatics, 2020.

Use make_smoment_model or with_smoment to construct the models.

The model is constructed as follows:

  • stoichiometry of the original model is retained as much as possible, but enzymatic reations are split into forward and reverse parts (marked by a suffix like ...#forward and ...#reverse),
  • coupling is added to simulate a virtual metabolite "enzyme capacity", which is consumed by all enzymatic reactions at a rate given by enzyme mass divided by the corresponding kcat,
  • the total consumption of the enzyme capacity is constrained to a fixed maximum.

The SMomentModel structure contains a worked-out representation of the optimization problem atop a wrapped MetabolicModel, in particular the separation of certain reactions into unidirectional forward and reverse parts, an "enzyme capacity" required for each reaction, and the value of the maximum capacity constraint. Original coupling in the inner model is retained.

In the structure, the field columns describes the correspondence of stoichiometry columns to the stoichiometry and data of the internal wrapped model, and total_enzyme_capacity is the total bound on the enzyme capacity consumption as specified in sMOMENT algorithm.

This implementation allows easy access to fluxes from the split reactions (available in reactions(model)), while the original "simple" reactions from the wrapped model are retained as fluxes. All additional constraints are implemented using coupling and coupling_bounds.

Fields

  • columns::Vector{COBREXA._smoment_column}

  • total_enzyme_capacity::Float64

  • inner::MetabolicModel

source
COBREXA._smoment_columnType
struct _smoment_column

A helper type that describes the contents of SMomentModels.

Fields

  • reaction_idx::Int64

  • direction::Int64

  • lb::Float64

  • ub::Float64

  • capacity_required::Float64

source
COBREXA.couplingMethod
coupling(model::SMomentModel) -> Any

Return the coupling of SMomentModel. That combines the coupling of the wrapped model, coupling for split reactions, and the coupling for the total enzyme capacity.

source
COBREXA.reactionsMethod
reactions(model::SMomentModel) -> Vector

Returns the internal reactions in a SMomentModel (these may be split to forward- and reverse-only parts; reactions IDs are mangled accordingly with suffixes).

source
COBREXA.stoichiometryMethod
stoichiometry(model::SMomentModel) -> Any

Return a stoichiometry of the SMomentModel. The enzymatic reactions are split into unidirectional forward and reverse ones.

source