amici-dev / amici

Advanced Multilanguage Interface to CVODES and IDAS

Advanced Multilanguage Interface for CVODES and IDAS

About

AMICI provides a multi-language (Python, C++, Matlab) interface for the SUNDIALS solvers CVODES (for ordinary differential equations) and IDAS (for algebraic differential equations). AMICI allows the user to read differential equation models specified as SBML or PySB and automatically compiles such models into Python modules, C++ libraries or Matlab .mex simulation files.

In contrast to the (no longer maintained) sundialsTB Matlab interface, all necessary functions are transformed into native C++ code, which allows for a significantly faster simulation.

Beyond forward integration, the compiled simulation file also allows for forward sensitivity analysis, steady state sensitivity analysis and adjoint sensitivity analysis for likelihood-based output functions.

The interface was designed to provide routines for efficient gradient computation in parameter estimation of biochemical reaction models, but it is also applicable to a wider range of differential equation constrained optimization problems.

Current build status

Features

• SBML import
• PySB import
• Generation of C++ code for model simulation and sensitivity computation
• Access to and high customizability of CVODES and IDAS solver
• Python, C++, Matlab interface
• Sensitivity analysis
  • forward
  • steady state
  • adjoint
  • first- and second-order
• Pre-equilibration and pre-simulation conditions
• Support for discrete events and logical operations

Interfaces & workflow

The AMICI workflow starts with importing a model from either SBML (Matlab, Python), PySB (Python), or a Matlab definition of the model (Matlab-only). From this input, all equations for model simulation are derived symbolically and C++ code is generated. This code is then compiled into a C++ library, a Python module, or a Matlab .mex file and is then used for model simulation.

Getting started

The AMICI source code is available at https://github.com/AMICI-dev/AMICI/. To install AMICI, first read the installation instructions for Python, C++ or Matlab. There are also instructions for using AMICI inside containers.

To get you started with Python-AMICI, the best way might be checking out this Jupyter notebook .

To get started with Matlab-AMICI, various examples are available in matlab/examples/.

Comprehensive documentation is available at https://amici.readthedocs.io/en/latest/.

Any contributions to AMICI are welcome (code, bug reports, suggestions for improvements, ...).

Getting help

In case of questions or problems with using AMICI, feel free to post an issue on GitHub. We are trying to get back to you quickly.

Projects using AMICI

There are several tools for parameter estimation offering good integration with AMICI:

• pyPESTO: Python library for optimization, sampling and uncertainty analysis
• pyABC: Python library for parallel and scalable ABC-SMC (Approximate Bayesian Computation - Sequential Monte Carlo)
• parPE: C++ library for parameter estimation of ODE models offering distributed memory parallelism with focus on problems with many simulation conditions.

Publications

Citeable DOI for the latest AMICI release:

There is a list of publications using AMICI. If you used AMICI in your work, we are happy to include your project, please let us know via a GitHub issue.

When using AMICI in your project, please cite

• Fröhlich, F., Weindl, D., Schälte, Y., Pathirana, D., Paszkowski, Ł., Lines, G.T., Stapor, P. and Hasenauer, J., 2021. AMICI: High-Performance Sensitivity Analysis for Large Ordinary Differential Equation Models. Bioinformatics, btab227, DOI:10.1093/bioinformatics/btab227.

@article{frohlich2020amici,
  title={AMICI: High-Performance Sensitivity Analysis for Large Ordinary Differential Equation Models},
  author={Fr{\"o}hlich, Fabian and Weindl, Daniel and Sch{\"a}lte, Yannik and Pathirana, Dilan and Paszkowski, {\L}ukasz and Lines, Glenn Terje and Stapor, Paul and Hasenauer, Jan},
  journal = {Bioinformatics},
  year = {2021},
  month = {04},
  issn = {1367-4803},
  doi = {10.1093/bioinformatics/btab227},
  note = {btab227},
  eprint = {https://academic.oup.com/bioinformatics/advance-article-pdf/doi/10.1093/bioinformatics/btab227/36866220/btab227.pdf},
}

When presenting work that employs AMICI, feel free to use one of the icons in documentation/gfx/, which are available under a CC0 license:

amici's Issues

warning for creation of symbolic variables

Currently matlab produces the warning:

Warning: Support of character vectors that are not valid variable names or define a number will be removed in a future release. To create
symbolic expressions, first create symbolic variables and then use operations on them.

This will require changing lots of code and is likely going to break many things in one of the next MATLAB releases.

undefined externals compiling model to mexw64

In amiwrap and compileC, the following externals are undefined:
__imp___iob_func
__imp_fprintf
__imp_vsprintf
__imp_printf
In addition, sprint is causing a warning, even though stdio.h is imported.
The details are in the attachment.
error-2015-11-14.docx

memory leak on integration errors

ami_if(false,...,...) and ami_if(true,...,...) not behaving as expected

This is due to heaviside(0) = 1/2. We need to explicitely check for those cases in ami_if.

automatically add and remove auxiliary folder from path when calling amiwrap

interface to FATODE

FATODE (http://people.cs.vt.edu/~asandu/Software/FATODE/index.html) provides explicit and implicit Runge-Kutta methods with support for forward and adjoint sensitivities and rootfinding. I expect that RK methods perform much better for models with many events and for adjoint sensitivities with many timepoints as repeated restarting hurts multi-step methods a lot.

return status not properly initialised

example_7 and example_8 are broken

example_7 and example_8 are DAE models which are not extensively tested yet.

const model dimensions

make fields in udata struct const and make init_modeldims the constructor.

Rewrite in C++

Moving all MATLAB code to C++ would provide many benefits
a) speedup of symbolic processing
b) better accessibility from other toolboxes
c) easier continuous integration

possible alternatives for the matlab symbolic toolbox are ginac (http://www.ginac.de) and reduce (http://www.reduce-algebra.com). reduce is much more powerful compared to ginac, but written in lisp 😞. For code generation in ginac we could use the routines from http://www.higgs.de/~stefanw/sector_decomposition/

installAMICI references non-existing directory

installAMICI
Warning: Name is nonexistent or not a directory:$HOME/src/AMICI.github/auxiliary/datahash
In path (line 109)
In addpath (line 88)
In installAMICI (line 4)

memory leak on integration failure

For the ribosomal transfection model I encounter memory leaks caused by
N_VCloneEmpty_Serial and setupAMI after integrations failures at the discontinuity.

rewrite generateM

https://github.com/ICB-DCM/AMICI/blob/master/%40amimodel/generateM.m is shitty copypasta zombiecodewhich repeatedly caused trouble now. Rewrite or move to C part ASAP.

faster model compilation for large models

for large models recompilation can be slow due to necessary preprocessing of symbolic equations. we can avoid that by saving the model struct as a mat which carries the hash of the syms file in its name. Then we always load the (processed) model unless the syms file changes.

Separate output parameter and dynamic parameters

Currently, sensitivities are computed for all model parameters and not only for those that have nonzero state sensitivities. This does not matter too much for small systems but can become problematic for larger systems.

amidata and amioption are slow

especially for small models, the passing/initialisation of those object can take considerable amount of time. implement checks for the class to avoid the call to the constructor.

• make amioption fast
• make amidata fast

evaluation of jacobian and xdot with adjoint sensitivities

those should be evaluated at the end of the time interval and not at the beginning for adjoint sensitivities

Amiwrap.m: File-existance test necessary?

Is the following check in amiwrap.m necessary?

    if(exist(symfun,'file') ~= 2)
        error(['"' symfun '" must be the name of a matlab function in the matlab path. Please check whether the folder containing "' symfun '" is in the matlab path. AMICI currently does not support absolute or relative paths in its input arguments.'])
    end

Amiwrap expects a function name, but checks if the file exists. I think this is an unnecessary limitation.

For example, it limits the use of package functions such as:
'''
amiwrap('bla', 'package.bla', pwd);
'''

I suggest replacing it by try/catch.

missing tests

We want tests for the following things

• correct return flags on errors
• correct simulations results
• correct sensitivities
• correct generation of model specific functions
• memory leaks via valgrind

forward and adjoint options

these model options are currently not functional. Either fix them or remove the functionality completely (the differences is relatively small anyways).

fix llh output for event objective

• change getEventObjective to use fJz
• change symbolic expression of Jz to be event-output-specific (like time-resolved objective Jy)
• implement routines for .sllh for forward sensitivities, check adjoint sensitivities (they have never been tested
• implement additional example (neuron model from bioinformatics paper is probably a good fit)
• add second order testcases

update CVODES and SuiteSparse

For both toolboxes new versions are available

sllh field for forward sensitivities

Actually all the code is in place to fill this field, yet is currently not done for forward sensitivities.

Support for SBRML and csv import export

It would be really great to have an outport/export function from amidata to csv/xls data aswell as SBRML. The simulation result could then also be saved in the same format.

amioptions is a handle

this can lead to unexpected behaviour. the superclass should be changed from hgetset (handle class) to some value class with similar properties (SolverOptions maybe?)

Accessor macros

The current accessor macros in [eur]data.h are causing trouble when including header files with function prototypes having arguments with the same names as the accessor macros (not unlikely with names as "p" or "k"). This becomes a pain in the ass when including external libraries and using AMICI-generated C-code in bigger projects.
Although tedious, what do you think about prefixing the accessor functions with "AMI(CI)_" and/or write them in uppercase to be clearly identifiable as preprocessor macros?

Commit number in simulation routine

It would be great to have a comment in the simulation routine that indicates the git commit number that was used to generate the simulation routine. This would increase reproducibility and help debugging.

condition assignments are not handled correctly

SBML events that make assignments to species with boundaryCondition="true" or constant="true" attributes are not handled correctly as they should only be triggered when the trigger transitions from false to true. Currently, they are also triggered when the trigger transitions from true to false.

This requires substantial rewriting of code as the current handling via heaviside functions cannot be fixed but needs to be replaced by events which can make changes constants. Alternatively we could also keep them as states, which might come at some additional cost but would generally simplify SBML import.

This is not a fringe SBML feature, but can be expected to frequently occur in models that are relevant to us.

improve code for functions that involve sparse matrix vector multiplications

The efficency/size of several functions could be dramatically reduced if we move from actual expressions to Sparse Matrix Vector Multiplication. This would improve vectorisation (this gives us the speedup) and provide more compact formulas (this gives us the size reduction). Also code generation would be faster as it would no longer be necessary to compute respective expressions.
Balancing formula size, readability and execution time is tricky, so think carefully where to apply this.

sol.srz, sol.s2rz and sol.s2z are wrong for second order forward sensitivities

😭

Sensitivities not returned by model simulation

The simulation of a model is not returning SX or SY. This might be a little bug, or just a user error. Details are in the attachment.
error-2015-11-17.docx

Memory deallocation in case of errors

The code for memory cleanup is quite messy at the moment. Especially in the case of error I am not too sure that we do not run into memory leaks. Here the structure should definitely be improved (set flags whether certain parts of the memory have been allocated?)

Beginners problems with amiwrap ad installToolbox

Here are my notes, sorry for the length. My workaround for the installToolbox is at the end.

What does “@” mean in “@amimodel”? According to the warning below, this might be a “method directory”, but I can’t find out what that is.

Does installToolbox.m need to run in AMICI directory?
addpath(genpath('auxiliary'))
Warning: Method directories not allowed in MATLAB path: @amimodel

In path (line 109)
In addpath (line 86)
In installToolbox (line 2)
In MakeModels (line 4)

Error using amiwrap (line 23)
symfun must be the name of a matlab function in the matlab path
Error in MakeModels (line 5)
amiwrap('M1Trivial','M1Trivial');
“symfun” has no meaning for user. Documentation calls it “example_model_syms” Simply “parameter 2” would mean more.

Now trying
amiwrap('M1Trivial','models/M1Trivial.m');

Error using eval
Undefined function or variable 'models'.
Error in amimodel (line 93)
model = eval(symfun);
Error in amiwrap (line 55)
model = amimodel(symfun,modelname);
Error in MakeModels (line 5)
amiwrap('M1Trivial','models/M1Trivial.m');

Now I tried renaming “M1Trivial.m” to “M1Trivial_syms.m”.
Same problems.

Now trying:
addpath('AMICI','AMICI/models');
cd('AMICI');
installToolbox;
cd('models');
model=M1Trivial_syms;
cd('..');
amiwrap('M1Trivial',model);
Error:
Error using amiwrap (line 20)
symfun must be a string
Error in MakeModels (line 9)
amiwrap('M1Trivial',model);
That disagrees with the document

Now I have put my make script in the models folder and I am getting much further. There is a problem with an undeclared identifier M_PI.
Error using mex
amiwrap.c
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\include/cvodewrap.h(107): warning C4098: 'AMIFree': 'void' function
returning a value
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(36): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(713): warning C4554: '|': check operator precedence for
possible error; use parentheses to clarify precedence
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(723): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(724): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(731): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(732): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(738): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(739): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(746): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(747): warning C4267: '=': conversion from 'size_t' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(916): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\src/amici.c(920): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(84): warning C4244: '=': conversion from 'double' to 'int',
possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(110): warning C4244: '=': conversion from 'double' to 'int',
possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(115): warning C4244: '=': conversion from 'double' to 'int',
possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(125): warning C4244: '=': conversion from 'double' to 'int',
possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(130): warning C4244: '=': conversion from 'double' to 'int',
possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(136): warning C4244: '=': conversion from 'double' to 'int',
possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(198): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(204): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(267): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(309): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(314): error C2065: 'M_PI': undeclared identifier
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(416): warning C4244: 'function': conversion from 'double' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(417): warning C4244: 'function': conversion from 'double' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(436): warning C4244: 'function': conversion from 'double' to
'int', possible loss of data
C:\Data\Tom\Research\Raedler\mRNA-Helmholtz\Tom\AMICI\amiwrap.c(437): warning C4244: 'function': conversion from 'double' to
'int', possible loss of data
Error in amimodel/compileC (line 381)
eval(['mex ' COPT ' ' CLIBS ' -output ' fullfile(wrap_path,'models',this.modelname,[prefix '_' this.modelname])
' ' fullfile(wrap_path,cstr) includesstr objectsstr ])
Error in amiwrap (line 76)
model.compileC();
Error in MakeModels (line 13)
amiwrap('M1Trivial','M1Trivial');

Workaround for installToolbox:
modelPath=fileparts(mfilename('fullpath'));
cd('..');
amiciPath=pwd;
addpath(amiciPath,modelPath);
% installToolbox;
addpath(fullfile(amiciPath,'auxiliary'));
addpath(fullfile(amiciPath,'auxiliary','datahash'));
% addpath(fullfile(amiciPath,'@amimodel'));
addpath(fullfile(amiciPath,'symbolic'));
cd(modelPath);

Move application of chain rule factor from .m to .c

Include field to parameter type in UserData, apply chain rule factor in amici.cpp, remove respective code from simulate___.m

compver

the compver property in amimodel is cumbersome, replace it by checksum over generateC and @amimodel methods.

amimodel.generateCstandalone

is this missing or deprecated?

Xdata structs

We should turn these structs into proper C++ classes with Functions which parse hdf5 or MATLAB rhs. This will make maintenance easier and improve readability/interpretability of code. We should probably also get rid of these accessor macros.

Note that udata and tempdata are used in C code (model specific functions).

Error when building deltaqB with 'Xcode Clang++'.

Hello,

I am using AMICI with Matlab R_2016b and Xcode version 8.2.1 for a DAE model. When I call: amiwrap(modelname,’example_model_syms’,dir,o2flag) to compile my model,
Mex completes successfully until I get the following errors:

"deltaqB | Building with 'Xcode Clang++'.
Error using mex
/Users/harry/GitHub/AMICI/models/MEC/MEC_deltaqB.cpp:17:5: error:
use of undeclared identifier 'ip'
for(ip = 0; ip<np; ip++) {
^
/Users/harry/GitHub/AMICI/models/MEC/MEC_deltaqB.cpp:17:13: error:
use of undeclared identifier 'ip'
for(ip = 0; ip<np; ip++) {
^
/Users/harry/GitHub/AMICI/models/MEC/MEC_deltaqB.cpp:17:20: error:
use of undeclared identifier 'ip'
for(ip = 0; ip<np; ip++) {
^
/Users/harry/GitHub/AMICI/models/MEC/MEC_deltaqB.cpp:18:15: error:
use of undeclared identifier 'ip'
switch (plist[ip]) {
^
4 errors generated.

Error in amimodel/compileC (line 379)
eval(['mex ' DEBUG COPT ...

Error in amiwrap (line 90)
model.compileC();"

Every time the MEC_deltaqB.cpp file is generated, it is missing the line:
int ip;
I am not sure why this line is not being generated. It may be a bug or a user error.

Best,
Harry Dudley

Event firing at initial timepoint

Event triggers that are 0 at timepoint t=t0 should fire their respective events before simulation starts. 6595513 only throws errors in trivial examples, but should be removed once fixed.

continuous integration and testing

we now have support for continuous integration (via travis-ci). currently this only builds a pre-generated version of model_dirac, but this should be extended to other models. Moreover we should also implement tests similar to amiexamples.

Conditional jump or move depends on unitialised value(s)

occurs for transfection model, first in initHeaviside. Probably what is causing more problems later on.

Error in example_dirac on staging

Error using symengine
The dimensions do not match.

Error in sym/privBinaryOp (line 937)
Csym = mupadmex(op,args{1}.s, args{2}.s, varargin{:});

Error in * (line 270)
X = privBinaryOp(A, B, 'symobj::mtimes');

Error in amifun/getSyms (line 499)
this.sym = diff(model.fun.root.sym,'t') +
model.fun.drootdx.sym*model.fun.xdot.sym_noopt;

Error in amimodel/getFun (line 48)
[fun,this] = fun.getSyms(this);

Error in amimodel/checkDeps (line 38)
this.getFun([],deps{id});

Error in amimodel/getFun (line 22)
cflag = this.checkDeps([],fun.deps);

Error in amimodel/checkDeps (line 38)
this.getFun([],deps{id});

Error in amimodel/getFun (line 25)
cflag = this.checkDeps([],fun.deps);

Error in amimodel/checkDeps (line 38)
this.getFun([],deps{id});

Error in amimodel/getFun (line 25)
cflag = this.checkDeps([],fun.deps);

Error in amimodel/parseModel (line 154)
this.getFun(HTable,funs{ifun});

Error in amiwrap (line 73)
model.parseModel();

Error in example_dirac (line 7)
amiwrap('model_dirac','model_dirac_syms',exdir)

====================================
For some reason, model.fun.xdot.sym_noopt is empty...

adjoint sensitivies for models with events

this has not been tested thoroughly yet and might fail

report error test failures and convergence test failures

error test failures -> bad linear solver/jacobian
convergence test failures -> ill-posed rhs

should be included in sol struct for more detailed system analysis

segmentation fault for second order sensitivities with event observables

title says it all 😞

Adjoint sensitivities for second order sensitivities

At the moment it is only possibly to compute (directional) second order sensitivities and in a forward-forward manner. Yet, forward-adjoint methods would be much more efficient.

very small doc issue

AMICI_guide.pdf
2.1.2 Options
.maxsteps should be 1e4, not 1e-8.
Should noforward and noadjoint be in the table?

Sensitivites are not working correctly for logic functions

maximal parameter number

we should define a maximal parameter number for every model somewhere and check against that.

Model recompilation

Sometimes the model compilation does not work properly, I suppose that this happens as some dependencies of symbolic variables are not correctly set and do not get updated. This will lead to crashing/incorrect simulations. The same goes for compilations that were stopped during the generation of .c files, AMICI should somehow detect and fix this. An easy workaround is deleting the respective folder in the models directory to force the recompilations.

Ideally the dependencies should be checked via a unit test that compiles only parts of the model.