amici.sbml_import.SbmlImporter¶
-
class
amici.sbml_import.
SbmlImporter
(sbml_source, show_sbml_warnings=False, from_file=True)[source]¶ Class to generate AMICI C++ files for a model provided in the Systems Biology Markup Language (SBML).
- Variables
show_sbml_warnings – indicates whether libSBML warnings should be displayed
symbols – dict carrying symbolic definitions
sbml_reader –
The libSBML sbml reader
Warning
Not storing this may result in a segfault.
sbml_doc –
document carrying the sbml definition
Warning
Not storing this may result in a segfault.
sbml – SBML model to import
compartments – dict of compartment ids and compartment volumes
stoichiometric_matrix – stoichiometric matrix of the model
flux_vector – reaction kinetic laws
_local_symbols – model symbols for sympy to consider during sympification see locals`argument in `sympy.sympify
species_assignment_rules – Assignment rules for species. Key is symbolic identifier and value is assignment value
compartment_assignment_rules – Assignment rules for compartments. Key is symbolic identifier and value is assignment value
parameter_assignment_rules – assignment rules for parameters, these parameters are not permissible for sensitivity analysis
initial_assignments – initial assignments for parameters, these parameters are not permissible for sensitivity analysis
sbml_parser_settings – sets behaviour of SBML Formula parsing
-
__init__
(sbml_source, show_sbml_warnings=False, from_file=True)[source] Create a new Model instance.
- Parameters
sbml_source (
typing.Union
[str
,libsbml.Model
]) – Either a path to SBML file where the model is specified, or a model string as created by sbml.sbmlWriter( ).writeSBMLToString() or an instance of libsbml.Model.show_sbml_warnings (
bool
) – Indicates whether libSBML warnings should be displayed.from_file (
bool
) – Whether sbml_source is a file name (True, default), or an SBML string
Methods Summary
__init__
(sbml_source[, show_sbml_warnings, …])Create a new Model instance.
add_d_dt
(d_dt, variable, variable0, name)Creates or modifies species, to implement rate rules for compartments and species, respectively.
add_local_symbol
(key, value)Add local symbols with some sanity checking for duplication which would indicate redefinition of internals, which SBML permits, but we don’t.
Check possible events in the model, as AMICI does currently not support
Check whether all required SBML features are supported.
is_assignment_rule_target
(element)Checks if an element has a valid assignment rule in the specified model.
is_rate_rule_target
(element)Checks if an element has a valid assignment rule in the specified model.
process_conservation_laws
(ode_model, …)Find conservation laws in reactions and species.
sbml2amici
([model_name, output_dir, …])Generate and compile AMICI C++ files for the model provided to the constructor.
Methods
-
__init__
(sbml_source, show_sbml_warnings=False, from_file=True)[source]¶ Create a new Model instance.
- Parameters
sbml_source (
typing.Union
[str
,libsbml.Model
]) – Either a path to SBML file where the model is specified, or a model string as created by sbml.sbmlWriter( ).writeSBMLToString() or an instance of libsbml.Model.show_sbml_warnings (
bool
) – Indicates whether libSBML warnings should be displayed.from_file (
bool
) – Whether sbml_source is a file name (True, default), or an SBML string
-
add_d_dt
(d_dt, variable, variable0, name)[source]¶ Creates or modifies species, to implement rate rules for compartments and species, respectively.
- Parameters
d_dt (
sympy.core.expr.Expr
) – The rate rule (or, right-hand side of an ODE).variable (
sympy.core.symbol.Symbol
) – The subject of the rate rule.variable0 (
typing.Union
[float
,sympy.core.expr.Expr
]) – The initial value of the variable.name (
str
) – Species name, only applicable if this function generates a new species
- Return type
-
add_local_symbol
(key, value)[source]¶ Add local symbols with some sanity checking for duplication which would indicate redefinition of internals, which SBML permits, but we don’t.
- Parameters
key (
str
) – local symbol keyvalue (
sympy.core.expr.Expr
) – local symbol value
-
check_event_support
()[source]¶ Check possible events in the model, as AMICI does currently not support
delays in events
priorities of events
events fired at initial time
Furthermore, event triggers are optional (e.g., if an event is fired at initial time, no trigger function is necessary). In this case, warn that this event will have no effect.
- Return type
-
check_support
()[source]¶ Check whether all required SBML features are supported. Also ensures that the SBML contains at least one reaction, or rate rule, or assignment rule, to produce change in the system over time.
- Return type
-
is_assignment_rule_target
(element)[source]¶ Checks if an element has a valid assignment rule in the specified model.
- Parameters
element (
libsbml.SBase
) – SBML variable- Return type
- Returns
boolean indicating truth of function name
-
is_rate_rule_target
(element)[source]¶ Checks if an element has a valid assignment rule in the specified model.
- Parameters
element (
libsbml.SBase
) – SBML variable- Return type
- Returns
boolean indicating truth of function name
-
process_conservation_laws
(ode_model, volume_updates)[source]¶ Find conservation laws in reactions and species.
- Parameters
ode_model – ODEModel object with basic definitions
volume_updates – List with updates for the stoichiometric matrix accounting for compartment volumes
- Returns volume_updates_solver
List (according to reduced stoichiometry) with updates for the stoichiometric matrix accounting for compartment volumes
- Return type
-
sbml2amici
(model_name=None, output_dir=None, observables=None, constant_parameters=None, sigmas=None, noise_distributions=None, verbose=40, assume_pow_positivity=False, compiler=None, allow_reinit_fixpar_initcond=True, compile=True, compute_conservation_laws=True, simplify=<function SbmlImporter.<lambda>>, log_as_log10=True, generate_sensitivity_code=True, **kwargs)[source]¶ Generate and compile AMICI C++ files for the model provided to the constructor.
The resulting model can be imported as a regular Python module (if compile=True), or used from Matlab or C++ as described in the documentation of the respective AMICI interface.
Note that this generates model ODEs for changes in concentrations, not amounts unless the hasOnlySubstanceUnits attribute has been defined for a particular species.
Sensitivity analysis for local parameters is enabled by creating global parameters _{reactionId}_{localParameterName}.
- Parameters
model_name (
typing.Optional
[str
]) – name of the model/model directoryoutput_dir (
typing.Optional
[str
]) – seeamici.ode_export.ODEExporter.set_paths()
observables (
typing.Optional
[typing.Dict
[str
,typing.Dict
[str
,str
]]]) – dictionary( observableId:{‘name’:observableName (optional), ‘formula’:formulaString)}) to be added to the modelconstant_parameters (
typing.Optional
[typing.Iterable
[str
]]) – list of SBML Ids identifying constant parameterssigmas (
typing.Optional
[typing.Dict
[str
,typing.Union
[str
,float
]]]) – dictionary(observableId: sigma value or (existing) parameter name)noise_distributions (
typing.Optional
[typing.Dict
[str
,typing.Union
[str
,typing.Callable
]]]) – dictionary(observableId: noise type). If nothing is passed for some observable id, a normal model is assumed as default. Either pass a noise type identifier, or a callable generating a custom noise string.verbose (
typing.Union
[int
,bool
]) – verbosity level for logging,True
/False
default tologging.Error
/logging.DEBUG
assume_pow_positivity (
bool
) – if set toTrue
, a special pow function is used to avoid problems with state variables that may become negative due to numerical errorscompiler (
typing.Optional
[str
]) – distutils/setuptools compiler selection to build the python extensionallow_reinit_fixpar_initcond (
bool
) – seeamici.ode_export.ODEExporter
compile (
bool
) – IfTrue
, compile the generated Python package, ifFalse
, just generate code.compute_conservation_laws (
bool
) – if set toTrue
, conservation laws are automatically computed and applied such that the state-jacobian of the ODE right-hand-side has full rank. This option should be set toTrue
when using the newton algorithm to compute steadystate sensitivities.simplify (
typing.Callable
) – seeODEModel._simplify
log_as_log10 (
bool
) – IfTrue
, log in the SBML model will be parsed aslog10
(default), ifFalse
, log will be parsed as natural logarithmln
generate_sensitivity_code (
bool
) – IfFalse
, the code required for sensitivity computation will not be generated
- Return type