SBML.jl — load systems biology models from SBML files

This package provides a straightforward way to load model- and simulation-relevant information from SBML files.

The library provides a single function readSBML to load a SBML.Model:

julia> using SBML
julia> mdl = readSBML("Ec_core_flux1.xml")
SBML.Model(…)

julia> mdl.compartments
2-element Array{String,1}:
 "Extra_organism"
 "Cytosol"

There are several functions to help you with using the data in the usual COBRA-style workflows, such as getS:

julia> metabolites, reactions, S = getS(mdl)
julia> metabolites
77-element Array{String,1}:
 "M_succoa_c"
 "M_ac_c"
 "M_etoh_c"
  ⋮

julia> S
77×77 SparseArrays.SparseMatrixCSC{Float64,Int64} with 308 stored entries:
  [60,  1]  =  -1.0
  [68,  1]  =  1.0
  [1 ,  2]  =  1.0
  [6 ,  2]  =  -1.0
  ⋮
  [23, 76]  =  1.0
  [56, 76]  =  -1.0
  [30, 77]  =  -1.0
  [48, 77]  =  1.0

julia> Matrix(S)
77×77 Array{Float64,2}:
 0.0   1.0  0.0  0.0  0.0  0.0  0.0  …  0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  1.0  0.0  0.0  0.0  0.0     0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0  -1.0  0.0  0.0  0.0  0.0  0.0  …  0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0   0.0  1.0  -1.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0   0.0  0.0   0.0  0.0  0.0  0.0
 0.0   0.0  0.0  0.0  0.0  0.0  0.0     0.0  -1.0  0.0   0.0  0.0  0.0  0.0
 ⋮                         ⋮         ⋱  ⋮                          ⋮    

Function reference

Data structures

SBML.MaybeType
Maybe{X}

Type shortcut for "X or nothing" or "nullable X" in javaspeak. Name got inspired by our functional friends.

source
SBML.MathType

A simplified representation of MathML-specified math AST

source
SBML.MathApplyType

Function application ("call by name", no tricks allowed) in mathematical expression

source
SBML.MathConstType

A constant identified by name (usually something like pi, e or true) in mathematical expression

source
SBML.MathTimeType

A special value representing the current time of the simulation, with a special name.

source
SBML.MathValType

A literal value (usually a numeric constant) in mathematical expression

source
SBML.ModelType

Structure that collects the model-related data. Contains parameters, units, compartments, species and reactions and gene_products, and additional notes and annotation (also present internally in some of the data fields). The contained dictionaries are indexed by identifiers of the corresponding objects.

source
SBML.ReactionType

Reaction with stoichiometry that assigns reactants and products their relative consumption/production rates (accessible in field stoichiometry), lower/upper bounds (in tuples lb and ub, with unit names), and objective coefficient (oc). Also may contains notes and annotation.

source
SBML.SpeciesType

Species metadata – contains a human-readable name, a compartment identifier, formula, charge, and additional notes and annotation.

source
SBML.UnitPartType

Part of a measurement unit definition that corresponds to the SBML definition of Unit. For example, the unit "per square megahour", Mh^(-2), is written as:

SBML.UnitPart("second",  # base SI unit, this says we are measuring time
         -2,        # exponent, says "per square"
         6,         # log-10 scale of the unit, says "mega"
         1/3600)    # second-to-hour multiplier

Compound units (such as "volt-amperes" and "dozens of yards per ounce") are built from multiple UnitParts; see the definition of field units in SBML.Model.

source

Base functions

SBML.SBMLVersionMethod
function SBMLVersion()

Get the version of the used SBML library in Julia version format.

source
SBML.extractModelMethod

" function extractModel(mdl::VPtr)::SBML.Model

Take the SBMLModel_t pointer and extract all information required to make a valid SBML.Model structure.

source
SBML.getAssociationMethod
function getAssociation(x::VPtr)::GeneProductAssociation

Convert a pointer to SBML FbcAssociation_t to the GeneProductAssociation tree structure.

source
SBML.get_optional_stringMethod
get_optional_string(x::VPtr, fn_sym)::Maybe{String}

Like get_string, but returns nothing instead of throwing an exception.

This is used to get notes and annotations and several other things (see get_notes, get_annotations)

source
SBML.get_stringMethod
get_string(x::VPtr, fn_sym)::Maybe{String}

C-call the SBML function fn_sym with a single parameter x, interpret the result as a string and return it, or throw exception in case the pointer is NULL.

source
SBML.readSBMLFunction
function readSBML(fn::String, sbml_conversion = model->nothing)::SBML.Model

Read the SBML from a XML file in fn and return the contained SBML.Model.

The sbml_conversion is a function that does an in-place modification of the single parameter, which is the C pointer to the loaded SBML document (C type SBMLDocument*). Several functions for doing that are prepared, including set_level_and_version, libsbml_convert, and convert_simplify_math.

Example

m = readSBML("my_model.xml", doc -> begin
    set_level_and_version(3, 1)(doc)
    convert_simplify_math(doc)
end)
source

libsbml representation converters

The converters are intended to be used as parameters of readSBML.

SBML.libsbml_convertMethod
libsbml_convert(conversion_options::Vector{Pair{String, Dict{String, String}}})

A converter that runs the SBML conversion routine, with specified conversion options. The argument is a vector of pairs to allow specifying the order of conversions.

source
SBML.libsbml_convertMethod
libsbml_convert(converter::String; kwargs...)

Quickly construct a single run of a libsbml converter from keyword arguments.

Example

readSBML("example.xml", libsbml_convert("stripPackage", package="layout"))
source

Data helpers

SBML.check_errorsMethod
check_errors(success::Integer, doc::Ptr{Cvoid}, error::Exception)

If success is a 0-valued Integer (a logical false), then call get_error_messages to show the error messages reported by SBML in the doc document and throw the error if they are more than 1. success is typically the value returned by an SBML C function operating on doc which returns a boolean flag to signal a successful operation.

source
SBML.extensive_kinetic_mathMethod
extensive_kinetic_math(m::SBML.Model, formula::SBML.Math)

Convert a SBML math formula to "extensive" kinetic laws, where the references to species that are marked as not having only substance units are converted from amounts to concentrations.

Handling of units in the conversion process is ignored in this version.

source
SBML.getLBsMethod
getLBs(m::SBML.Model)::Vector{Tuple{Float64,String}}

Extract a vector of lower bounds of reaction rates from the model. All bounds are accompanied with the unit of the corresponding value (this behavior is based on SBML specification).

source
SBML.getOCsMethod
getOCs(m::SBML.Model)::Vector{Float64}

Extract the vector of objective coefficients of each reaction.

source
SBML.getSMethod
function getS(m::SBML.Model; zeros=spzeros)::Tuple{Vector{String},Vector{String},AbstractMatrix{Float64}}

Extract the vector of species (aka metabolite) identifiers, vector of reaction identifiers, and the (dense) stoichiometry matrix from an existing SBML.Model. Returns a tuple with these values.

The matrix is sparse by default (initially constructed by SparseArrays.spzeros). You can fill in a custom empty matrix constructed to argument zeros; e.g. running with zeros=zeros will produce a dense matrix.

source
SBML.getUBsMethod
getUBs(m::SBML.Model)::Vector{Tuple{Float64,String}}

Likewise to getLBs, extract the upper bounds.

source
SBML.get_error_messagesMethod
get_error_messages(doc::Ptr{Cvoid}, error::Exception)

Show the error messages reported by SBML in the doc document and throw the error if they are more than 1.

source
SBML.initial_amountsMethod
initial_amounts(m::SBML.Model; convert_concentrations = false)

Return initial amounts for each species as a generator of pairs species_name => initial_amount; the amount is set to nothing if not available. If convert_concentrations is true and there is information about initial concentration available together with compartment size, the result is computed from the species' initial concentration.

In the current version, units of the measurements are completely ignored.

Example

# get the initial amounts as dictionary
Dict(initial_amounts(model, convert_concentrations = true))

# remove the empty entries
Dict(k => v for (k,v) in initial_amounts(model) if !isnothing(v))
source

Math and Symbolics.jl compatibility

SBML.default_symbolics_constantsConstant
const default_symbolics_constants::Dict{String, Any}

A dictionary of default constants filled in place of SBML Math constants in the symbolics conversion.

source
SBML.default_symbolics_mappingConstant
default_symbolics_mapping :: Dict{String,Any}

Default mapping of SBML function names to Julia functions, represented as a dictionary from Strings (SBML names) to anything evaluable as Julia&Symbolics functions, such as symbols and expressions.

The default mapping only contains the basic SBML functions that are unambiguously represented in Julia; it is supposed to be extended by the user if more functions need to be supported.

source
Base.convertMethod
Base.convert(
    ::Type{Num},
    x::SBML.Math;
    mapping = default_symbolics_mapping,
    convert_time = (x::SBML.MathTime) -> Num(Variable(Symbol(x.id))).val,
    convert_const = (x::SBML.MathConst) -> Num(default_symbolics_constants[x.id]),
)

Convert SBML.Math to Num type from Symbolics package. The conversion of functions can be customized by supplying a custom mapping; if nothing is supplied, default_symbolics_mapping that translates basic functions to their Julia equivalents is assumed.

Translation of MathLambda is not supported by Symbolics.

MathTime is handled specially, the function from the argument convert_time is called to possibly specify any desired behavior. By default, it just creates a variable with the same name as the time variable name stored in SBML.

source

Internal math helpers

SBML.ast_isMethod
ast_is(ast::VPtr, what::Symbol)::Bool

Helper for quickly recognizing kinds of ASTs

source
SBML.parse_mathMethod
parse_math(ast::VPtr)::Math

This attempts to parse out a decent Julia-esque (Math AST from a pointer to ASTNode_t.

source