AMICI Python example "steadystate"

This is an example using the [model_steadystate_scaled.sbml] model to demonstrate and test SBML import and AMICI Python interface.


In [1]:
# SBML model we want to import
sbml_file = 'model_steadystate_scaled.xml'
# Name of the model that will also be the name of the python module
model_name = 'model_steadystate_scaled'
# Directory to which the generated model code is written
model_output_dir = model_name

import libsbml
import importlib
import amici
import os
import sys
import numpy as np
import matplotlib.pyplot as plt

The example model

Here we use libsbml to show the reactions and species described by the model (this is independent of AMICI).


In [2]:
sbml_reader = libsbml.SBMLReader()
sbml_doc = sbml_reader.readSBML(sbml_file)
sbml_model = sbml_doc.getModel()
dir(sbml_doc)

print('Species: ', [s.getId() for s in sbml_model.getListOfSpecies()])

print('\nReactions:')
for reaction in sbml_model.getListOfReactions():
    reactants = ' + '.join(['%s %s'%(int(r.getStoichiometry()) if r.getStoichiometry() > 1 else '', r.getSpecies()) for r in reaction.getListOfReactants()])
    products  = ' + '.join(['%s %s'%(int(r.getStoichiometry()) if r.getStoichiometry() > 1 else '', r.getSpecies()) for r in reaction.getListOfProducts()])
    reversible = '<' if reaction.getReversible() else ''
    print('%3s: %10s %1s->%10s\t\t[%s]' % (reaction.getId(), 
                        reactants,
                        reversible,
                         products,
                        libsbml.formulaToL3String(reaction.getKineticLaw().getMath())))


Species:  ['x1', 'x2', 'x3']

Reactions:
 r1:       2 x1  ->        x2		[p1 * x1^2]
 r2:   x1 +  x2  ->        x3		[p2 * x1 * x2]
 r3:         x2  ->      2 x1		[p3 * x2]
 r4:         x3  ->  x1 +  x2		[p4 * x3]
 r5:         x3  ->          		[k0 * x3]
 r6:             ->        x1		[p5]

Importing an SBML model, compiling and generating an AMICI module

Before we can use AMICI to simulate our model, the SBML model needs to be translated to C++ code. This is done by amici.SbmlImporter.


In [3]:
# Create an SbmlImporter instance for our SBML model
sbml_importer = amici.SbmlImporter(sbml_file)

In this example, we want to specify fixed parameters, observables and a $\sigma$ parameter. Unfortunately, the latter two are not part of the SBML standard. However, they can be provided to amici.SbmlImporter.sbml2amici as demonstrated in the following.

Constant parameters

Constant parameters, i.e. parameters with respect to which no sensitivities are to be computed (these are often parameters specifying a certain experimental condition) are provided as a list of parameter names.


In [4]:
constantParameters = ['k0']

Observables

Specifying observables is beyond the scope of SBML. Here we define them manually.

If you are looking for a more scalable way for defining observables, ideally inside the SBML file, then checkout https://github.com/ICB-DCM/PEtab. There use a convention based on SBML's AssignmentRule to specify Model outputs within the SBML file.


In [5]:
# Define observables
observables = {
    'observable_x1': {'name': '', 'formula': 'x1'}, 
    'observable_x2': {'name': '', 'formula': 'x2'}, 
    'observable_x3': {'name': '', 'formula': 'x3'}, 
    'observable_x1_scaled': {'name': '', 'formula': 'scaling_x1 * x1'}, 
    'observable_x2_offsetted': {'name': '', 'formula': 'offset_x2 + x2'}, 
    'observable_x1withsigma': {'name': '', 'formula': 'x1'}
}

$\sigma$ parameters

To specify measurement noise as a parameter, we simply provide a dictionary with (preexisting) parameter names as keys and a list of observable names as values to indicate which sigma parameter is to be used for which observable.


In [6]:
sigmas = {'observable_x1withsigma': 'observable_x1withsigma_sigma'}

Generating the module

Now we can generate the python module for our model. amici.SbmlImporter.sbml2amici will symbolically derive the sensitivity equations, generate C++ code for model simulation, and assemble the python module.


In [7]:
sbml_importer.sbml2amici(model_name, 
                         model_output_dir, 
                         verbose=False,
                         observables=observables,
                         constantParameters=constantParameters,
                         sigmas=sigmas)

Importing the module and loading the model

If everything went well, we need to add the previously selected model output directory to our PYTHON_PATH and are then ready to load newly generated model:


In [8]:
sys.path.insert(0, os.path.abspath(model_output_dir))
model_module = importlib.import_module(model_name)

And get an instance of our model from which we can retrieve information such as parameter names:


In [9]:
model = model_module.getModel()

print("Model parameters:", model.getParameterIds())
print("Model outputs:   ", model.getObservableIds())
print("Model states:    ", model.getStateIds())


Model parameters: ('p1', 'p2', 'p3', 'p4', 'p5', 'scaling_x1', 'offset_x2', 'observable_x1withsigma_sigma')
Model outputs:    ('observable_x1', 'observable_x2', 'observable_x3', 'observable_x1_scaled', 'observable_x2_offsetted', 'observable_x1withsigma')
Model states:     ('x1', 'x2', 'x3')

Using the python-generated AMICI model from Matlab

It is also possible to use a Python-AMICI imported model from Matlab. You might want to do this because:

  • you don't have the Symbolic Math Toolbox available
  • model generation using the AMICI Python interface can be significantly faster
  • python SBML import has some extra features that are not available through the matlab interface
  • model equations will differ between Python- and Matlab-import

In [10]:
print('''To use python-generated model in Matlab, ensure you added AMICI to your matlab path and run:

modelName = '{model_name}';
modelDir = '{model_output_dir}';
amimodel.compileAndLinkModel(modelName, modelDir, [], [], [], []);
amimodel.generateMatlabWrapper({nx}, {ny}, {np}, {nk}, {nz}, {o2flag}, [], [ modelDir '/simulate_' modelName '.m'], modelName, 'lin', 1, 1);
'''.format(model_name=model_name, 
           model_output_dir=os.path.abspath(model_output_dir),
           nx = model.nxtrue_rdata,
           ny = model.nytrue,
           np = model.np(),
           nk = model.nk(),
           nz = model.nz,
           o2flag = model.o2mode
          ))


To use python-generated model in Matlab, ensure you added AMICI to your matlab path and run:

modelName = 'model_steadystate_scaled';
modelDir = '/home/dweindl/src/AMICI-devel/python/examples/example_steadystate/model_steadystate_scaled';
amimodel.compileAndLinkModel(modelName, modelDir, [], [], [], []);
amimodel.generateMatlabWrapper(3, 6, 8, 1, 0, 0, [], [ modelDir '/simulate_' modelName '.m'], modelName, 'lin', 1, 1);

This will use the matlab compiler to generate a mex file and will create a Matlab wrapper script. You have to run this only once after generating the Python model. Afterwards you can use the model from Matlab just as if it was generted using amiwrap.m directly.

Running simulations and analyzing results

After importing the model, we can run simulations using amici.runAmiciSimulation. This requires a Model instance and a Solver instance. Optionally you can provide measurements inside an ExpData instance, as shown later in this notebook.


In [11]:
# Create Model instance
model = model_module.getModel()

# set timepoints for which we want to simulate the model
model.setTimepoints(np.linspace(0, 60, 60)) 

# Create solver instance
solver = model.getSolver()

# Run simulation using default model parameters and solver options
rdata = amici.runAmiciSimulation(model, solver)

In [13]:
print('Simulation was run using model default parameters as specified in the SBML model:')
print(model.getParameters())


Simulation was run using model default parameters as specified in the SBML model:
(1.0, 0.5, 0.4, 2.0, 0.1, 2.0, 3.0, 0.2)

Simulation results are provided as numpy.ndarrays in the returned dictionary:


In [14]:
#np.set_printoptions(threshold=8, edgeitems=2)
for key, value in rdata.items():
    print('%12s: ' % key, value)


         ptr:  <amici.amici.ReturnDataPtr; proxy of <Swig Object of type 'std::unique_ptr< amici::ReturnData > *' at 0x7fdc43ec3cc0> >
           t:  [ 0.          1.01694915  2.03389831  3.05084746  4.06779661  5.08474576
  6.10169492  7.11864407  8.13559322  9.15254237 10.16949153 11.18644068
 12.20338983 13.22033898 14.23728814 15.25423729 16.27118644 17.28813559
 18.30508475 19.3220339  20.33898305 21.3559322  22.37288136 23.38983051
 24.40677966 25.42372881 26.44067797 27.45762712 28.47457627 29.49152542
 30.50847458 31.52542373 32.54237288 33.55932203 34.57627119 35.59322034
 36.61016949 37.62711864 38.6440678  39.66101695 40.6779661  41.69491525
 42.71186441 43.72881356 44.74576271 45.76271186 46.77966102 47.79661017
 48.81355932 49.83050847 50.84745763 51.86440678 52.88135593 53.89830508
 54.91525424 55.93220339 56.94915254 57.96610169 58.98305085 60.        ]
           x:  [[0.1        0.4        0.7       ]
 [0.57995051 0.73365809 0.0951589 ]
 [0.55996496 0.71470091 0.0694127 ]
 [0.5462855  0.68030366 0.06349394]
 [0.53561883 0.64937432 0.05923555]
 [0.52636487 0.62259567 0.05568686]
 [0.51822014 0.59943346 0.05268079]
 [0.51103767 0.57935661 0.05012037]
 [0.5047003  0.56191593 0.04793052]
 [0.49910666 0.54673518 0.0460508 ]
 [0.49416809 0.53349812 0.04443206]
 [0.48980688 0.52193768 0.043034  ]
 [0.48595477 0.51182733 0.04182339]
 [0.48255177 0.50297415 0.04077267]
 [0.47954512 0.49521321 0.03985882]
 [0.47688834 0.48840307 0.03906254]
 [0.4745405  0.48242201 0.03836756]
 [0.47246548 0.47716503 0.0377601 ]
 [0.47063147 0.4725413  0.03722844]
 [0.46901037 0.46847203 0.03676259]
 [0.4675774  0.46488882 0.03635397]
 [0.46631066 0.46173208 0.03599523]
 [0.46519082 0.45894988 0.03568002]
 [0.46420083 0.45649685 0.03540285]
 [0.4633256  0.45433333 0.03515899]
 [0.46255181 0.45242458 0.0349443 ]
 [0.46186769 0.45074017 0.0347552 ]
 [0.46126283 0.44925338 0.03458856]
 [0.46072804 0.44794076 0.03444166]
 [0.46025521 0.44678169 0.03431212]
 [0.45983714 0.44575805 0.03419785]
 [0.4594675  0.44485389 0.03409701]
 [0.45914066 0.44405515 0.03400802]
 [0.45885167 0.44334948 0.03392946]
 [0.45859615 0.44272595 0.03386009]
 [0.45837021 0.44217497 0.03379883]
 [0.45817044 0.44168806 0.03374473]
 [0.45799379 0.44125773 0.03369693]
 [0.4578376  0.44087739 0.03365471]
 [0.45769949 0.44054121 0.0336174 ]
 [0.45757737 0.44024405 0.03358444]
 [0.45746939 0.43998138 0.03355531]
 [0.45737391 0.43974917 0.03352956]
 [0.45728949 0.4395439  0.03350681]
 [0.45721483 0.43936243 0.0334867 ]
 [0.45714882 0.43920199 0.03346892]
 [0.45709046 0.43906015 0.03345321]
 [0.45703884 0.43893475 0.03343932]
 [0.45699321 0.43882387 0.03342704]
 [0.45695285 0.43872585 0.03341618]
 [0.45691717 0.43863918 0.03340658]
 [0.45688562 0.43856255 0.0333981 ]
 [0.45685772 0.43849479 0.0333906 ]
 [0.45683305 0.43843488 0.03338397]
 [0.45681123 0.43838191 0.0333781 ]
 [0.45679194 0.43833508 0.03337292]
 [0.45677489 0.43829366 0.03336833]
 [0.4567598  0.43825705 0.03336428]
 [0.45674647 0.43822467 0.0333607 ]
 [0.45673468 0.43819604 0.03335753]]
          x0:  [0.1 0.4 0.7]
        x_ss:  [nan nan nan]
          sx:  None
         sx0:  None
       sx_ss:  None
           y:  [[0.1        0.4        0.7        0.2        3.4        0.1       ]
 [0.57995051 0.73365809 0.0951589  1.15990103 3.73365809 0.57995051]
 [0.55996496 0.71470091 0.0694127  1.11992992 3.71470091 0.55996496]
 [0.5462855  0.68030366 0.06349394 1.092571   3.68030366 0.5462855 ]
 [0.53561883 0.64937432 0.05923555 1.07123766 3.64937432 0.53561883]
 [0.52636487 0.62259567 0.05568686 1.05272975 3.62259567 0.52636487]
 [0.51822014 0.59943346 0.05268079 1.03644027 3.59943346 0.51822014]
 [0.51103767 0.57935661 0.05012037 1.02207534 3.57935661 0.51103767]
 [0.5047003  0.56191593 0.04793052 1.0094006  3.56191593 0.5047003 ]
 [0.49910666 0.54673518 0.0460508  0.99821332 3.54673518 0.49910666]
 [0.49416809 0.53349812 0.04443206 0.98833619 3.53349812 0.49416809]
 [0.48980688 0.52193768 0.043034   0.97961375 3.52193768 0.48980688]
 [0.48595477 0.51182733 0.04182339 0.97190953 3.51182733 0.48595477]
 [0.48255177 0.50297415 0.04077267 0.96510354 3.50297415 0.48255177]
 [0.47954512 0.49521321 0.03985882 0.95909024 3.49521321 0.47954512]
 [0.47688834 0.48840307 0.03906254 0.95377669 3.48840307 0.47688834]
 [0.4745405  0.48242201 0.03836756 0.94908099 3.48242201 0.4745405 ]
 [0.47246548 0.47716503 0.0377601  0.94493097 3.47716503 0.47246548]
 [0.47063147 0.4725413  0.03722844 0.94126294 3.4725413  0.47063147]
 [0.46901037 0.46847203 0.03676259 0.93802075 3.46847203 0.46901037]
 [0.4675774  0.46488882 0.03635397 0.93515479 3.46488882 0.4675774 ]
 [0.46631066 0.46173208 0.03599523 0.93262131 3.46173208 0.46631066]
 [0.46519082 0.45894988 0.03568002 0.93038165 3.45894988 0.46519082]
 [0.46420083 0.45649685 0.03540285 0.92840166 3.45649685 0.46420083]
 [0.4633256  0.45433333 0.03515899 0.9266512  3.45433333 0.4633256 ]
 [0.46255181 0.45242458 0.0349443  0.92510362 3.45242458 0.46255181]
 [0.46186769 0.45074017 0.0347552  0.92373537 3.45074017 0.46186769]
 [0.46126283 0.44925338 0.03458856 0.92252566 3.44925338 0.46126283]
 [0.46072804 0.44794076 0.03444166 0.92145609 3.44794076 0.46072804]
 [0.46025521 0.44678169 0.03431212 0.92051042 3.44678169 0.46025521]
 [0.45983714 0.44575805 0.03419785 0.91967428 3.44575805 0.45983714]
 [0.4594675  0.44485389 0.03409701 0.91893499 3.44485389 0.4594675 ]
 [0.45914066 0.44405515 0.03400802 0.91828132 3.44405515 0.45914066]
 [0.45885167 0.44334948 0.03392946 0.91770334 3.44334948 0.45885167]
 [0.45859615 0.44272595 0.03386009 0.91719229 3.44272595 0.45859615]
 [0.45837021 0.44217497 0.03379883 0.91674042 3.44217497 0.45837021]
 [0.45817044 0.44168806 0.03374473 0.91634087 3.44168806 0.45817044]
 [0.45799379 0.44125773 0.03369693 0.91598759 3.44125773 0.45799379]
 [0.4578376  0.44087739 0.03365471 0.9156752  3.44087739 0.4578376 ]
 [0.45769949 0.44054121 0.0336174  0.91539899 3.44054121 0.45769949]
 [0.45757737 0.44024405 0.03358444 0.91515475 3.44024405 0.45757737]
 [0.45746939 0.43998138 0.03355531 0.91493879 3.43998138 0.45746939]
 [0.45737391 0.43974917 0.03352956 0.91474783 3.43974917 0.45737391]
 [0.45728949 0.4395439  0.03350681 0.91457897 3.4395439  0.45728949]
 [0.45721483 0.43936243 0.0334867  0.91442967 3.43936243 0.45721483]
 [0.45714882 0.43920199 0.03346892 0.91429765 3.43920199 0.45714882]
 [0.45709046 0.43906015 0.03345321 0.91418091 3.43906015 0.45709046]
 [0.45703884 0.43893475 0.03343932 0.91407769 3.43893475 0.45703884]
 [0.45699321 0.43882387 0.03342704 0.91398641 3.43882387 0.45699321]
 [0.45695285 0.43872585 0.03341618 0.9139057  3.43872585 0.45695285]
 [0.45691717 0.43863918 0.03340658 0.91383434 3.43863918 0.45691717]
 [0.45688562 0.43856255 0.0333981  0.91377123 3.43856255 0.45688562]
 [0.45685772 0.43849479 0.0333906  0.91371543 3.43849479 0.45685772]
 [0.45683305 0.43843488 0.03338397 0.91366609 3.43843488 0.45683305]
 [0.45681123 0.43838191 0.0333781  0.91362246 3.43838191 0.45681123]
 [0.45679194 0.43833508 0.03337292 0.91358389 3.43833508 0.45679194]
 [0.45677489 0.43829366 0.03336833 0.91354977 3.43829366 0.45677489]
 [0.4567598  0.43825705 0.03336428 0.91351961 3.43825705 0.4567598 ]
 [0.45674647 0.43822467 0.0333607  0.91349294 3.43822467 0.45674647]
 [0.45673468 0.43819604 0.03335753 0.91346935 3.43819604 0.45673468]]
      sigmay:  [[1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]]
          sy:  None
     ssigmay:  None
           z:  None
          rz:  None
      sigmaz:  None
          sz:  None
         srz:  None
     ssigmaz:  None
        sllh:  None
       s2llh:  None
           J:  [[-2.04603672  0.69437133  0.21909802]
 [ 0.57163266 -0.62836734  0.22836734]
 [ 2.          2.         -3.        ]]
        xdot:  [-1.08973588e-05 -2.64550670e-05 -2.92781904e-06]
      status:  0.0
         llh:  nan
        chi2:  nan
         res:  [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. 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. 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. 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. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]
        sres:  None
         FIM:  None
wrms_steadystate:  nan
t_steadystate:  nan
newton_numlinsteps:  None
newton_numsteps:  [[0 0 0]]
    numsteps:  [  0 100 145 172 188 199 206 210 214 217 221 224 226 229 231 233 236 238
 240 243 245 248 250 252 255 256 258 259 261 262 264 265 267 269 270 272
 273 275 276 278 279 281 282 283 284 285 286 287 288 289 290 291 292 293
 294 295 296 297 298 299]
 numrhsevals:  [  0 115 161 189 206 218 226 232 236 239 243 247 250 253 255 257 260 262
 264 267 269 272 274 276 279 281 283 284 286 287 289 290 292 294 295 297
 298 300 301 303 304 307 308 309 310 311 312 313 314 315 316 317 318 319
 320 321 322 323 324 325]
numerrtestfails:  [0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]
numnonlinsolvconvfails:  [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]
       order:  [0 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5
 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5]
   numstepsB:  None
numrhsevalsB:  None
numerrtestfailsB:  None
numnonlinsolvconvfailsB:  None

A note on datatypes

You might have wondered about expression such as rdata['y'].flatten() or amici.ExpData(model.get()). This is currently required to convert between Python data types and internal C++ datatypes. The SWIG interface does not yet do this conversion automatically.


In [15]:
print(model.getParameters())


(1.0, 0.5, 0.4, 2.0, 0.1, 2.0, 3.0, 0.2)

Plotting tractories

The simulation results above did not look too appealing. Let's plot the trajectories of the model states and outputs them using matplotlib.pyplot:


In [16]:
import amici.plotting
amici.plotting.plotStateTrajectories(rdata)
amici.plotting.plotObservableTrajectories(rdata)


Computing likelihood

Often model parameters need to be inferred from experimental data. This is commonly done by maximizing the likelihood of of observing the data given to current model parameters. AMICI will compute this likelihood if experimental data is provided to amici.runAmiciSimulation as optional third argument. Measurements along with their standard deviations are provided through an amici.ExpData instance.


In [17]:
# Create model instance and set time points for simulation
model = model_module.getModel()
model.setTimepoints(np.linspace(0, 10, 11)) 

# Create solver instance, keep default options
solver = model.getSolver()

# Run simulation without experimental data
rdata = amici.runAmiciSimulation(model, solver)

# Create ExpData instance from simulation results
edata = amici.ExpData(rdata, 1.0, 0.0)

# Re-run simulation, this time passing "experimental data"
rdata = amici.runAmiciSimulation(model, solver, edata)

print('Log-likelihood %f' % rdata['llh'])


Log-likelihood -87.232753

Sensitivity analysis

AMICI can provide first- and second-order sensitivities using the forward- or adjoint-method. The respective options are set on the Model and Solver objects.

Forward sensitivity analysis


In [18]:
model = model_module.getModel()
model.setTimepoints(np.linspace(0, 10, 11)) 
model.requireSensitivitiesForAllParameters()              # sensitivities w.r.t. all parameters
# model.setParameterList([1, 2])                          # sensitivities 
# w.r.t. the specified parameters
model.setParameterScale(amici.ParameterScaling_none)         # parameters are used as-is (not log-transformed)

solver = model.getSolver()
solver.setSensitivityMethod(amici.SensitivityMethod_forward)        # forward sensitivity analysis
solver.setSensitivityOrder(amici.SensitivityOrder_first) # first-order sensitivities

rdata = amici.runAmiciSimulation(model, solver)

# print sensitivity-related results
for key, value in rdata.items():
    if key.startswith('s'):
        print('%12s: ' % key, value)


          sx:  [[[ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-2.00747251e-01  1.19873139e-01 -9.44167984e-03]
  [-1.02561396e-01 -1.88820454e-01  1.01855972e-01]
  [ 4.66193077e-01 -2.86365372e-01  2.39662449e-02]
  [ 4.52560292e-02  1.14631370e-01 -3.34067917e-02]
  [ 4.00672911e-01  1.92564093e-01  4.98877760e-02]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-2.23007240e-01  1.53979022e-01 -1.26885280e-02]
  [-1.33426939e-01 -3.15955239e-01  9.49575028e-02]
  [ 5.03470377e-01 -3.52731535e-01  2.81567411e-02]
  [ 3.93630714e-02  1.10770683e-01 -1.05673871e-02]
  [ 5.09580305e-01  4.65255489e-01  9.24843702e-02]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-2.14278105e-01  1.63465065e-01 -1.03268419e-02]
  [-1.60981967e-01 -4.00490453e-01  7.54810651e-02]
  [ 4.87746420e-01 -3.76014315e-01  2.30919334e-02]
  [ 4.28733678e-02  1.15473583e-01 -6.63571668e-03]
  [ 6.05168647e-01  7.07226039e-01  1.23870915e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-2.05888039e-01  1.69308690e-01 -7.93085689e-03]
  [-1.84663808e-01 -4.65451966e-01  5.95026120e-02]
  [ 4.66407067e-01 -3.87612082e-01  1.76410135e-02]
  [ 4.52451106e-02  1.19865711e-01 -4.73313045e-03]
  [ 6.90798450e-01  9.20396632e-01  1.49475828e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-1.98803165e-01  1.73327269e-01 -6.03008194e-03]
  [-2.04303740e-01 -5.16111389e-01  4.68785776e-02]
  [ 4.47070329e-01 -3.94304031e-01  1.32107443e-02]
  [ 4.69732052e-02  1.22961726e-01 -3.35899414e-03]
  [ 7.68998997e-01  1.10844286e+00  1.70889328e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-1.92789113e-01  1.75978657e-01 -4.54517629e-03]
  [-2.20500139e-01 -5.55540706e-01  3.68776525e-02]
  [ 4.30424856e-01 -3.97907707e-01  9.75257132e-03]
  [ 4.82793653e-02  1.24952071e-01 -2.30991631e-03]
  [ 8.40805132e-01  1.27504628e+00  1.89020151e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-1.87672774e-01  1.77588334e-01 -3.38318223e-03]
  [-2.33807210e-01 -5.86081383e-01  2.89236334e-02]
  [ 4.16201398e-01 -3.99295277e-01  7.06598589e-03]
  [ 4.92546647e-02  1.26089711e-01 -1.50412008e-03]
  [ 9.06806543e-01  1.42334018e+00  2.04522708e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-1.83320440e-01  1.78410041e-01 -2.47240697e-03]
  [-2.44690164e-01 -6.09568484e-01  2.25774268e-02]
  [ 4.04061655e-01 -3.99063011e-01  4.97908397e-03]
  [ 4.99612483e-02  1.26581014e-01 -8.85891392e-04]
  [ 9.67473969e-01  1.55589415e+00  2.17895305e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-1.79620591e-01  1.78640115e-01 -1.75822432e-03]
  [-2.53540124e-01 -6.27448860e-01  1.75019835e-02]
  [ 3.93704969e-01 -3.97656642e-01  3.35895468e-03]
  [ 5.04492283e-02  1.26586733e-01 -4.13401193e-04]
  [ 1.02322336e+00  1.67481439e+00  2.29524046e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]

 [[-1.76478441e-01  1.78430281e-01 -1.19867655e-03]
  [-2.60686972e-01 -6.40868688e-01  1.34365064e-02]
  [ 3.84873835e-01 -3.95414932e-01  2.10369506e-03]
  [ 5.07601806e-02  1.26231632e-01 -5.46464877e-05]
  [ 1.07443160e+00  1.78183962e+00  2.39710937e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00]]]
         sx0:  [[0. 0. 0.]
 [0. 0. 0.]
 [0. 0. 0.]
 [0. 0. 0.]
 [0. 0. 0.]
 [0. 0. 0.]
 [0. 0. 0.]
 [0. 0. 0.]]
       sx_ss:  [[nan nan nan]
 [nan nan nan]
 [nan nan nan]
 [nan nan nan]
 [nan nan nan]
 [nan nan nan]
 [nan nan nan]
 [nan nan nan]]
      sigmay:  [[1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]
 [1.  1.  1.  1.  1.  0.2]]
          sy:  [[[ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  1.00000000e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-2.00747251e-01  1.19873139e-01 -9.44167984e-03 -4.01494501e-01
    1.19873139e-01 -2.00747251e-01]
  [-1.02561396e-01 -1.88820454e-01  1.01855972e-01 -2.05122791e-01
   -1.88820454e-01 -1.02561396e-01]
  [ 4.66193077e-01 -2.86365372e-01  2.39662449e-02  9.32386155e-01
   -2.86365372e-01  4.66193077e-01]
  [ 4.52560292e-02  1.14631370e-01 -3.34067917e-02  9.05120584e-02
    1.14631370e-01  4.52560292e-02]
  [ 4.00672911e-01  1.92564093e-01  4.98877760e-02  8.01345822e-01
    1.92564093e-01  4.00672911e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.80072436e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-2.23007240e-01  1.53979022e-01 -1.26885280e-02 -4.46014480e-01
    1.53979022e-01 -2.23007240e-01]
  [-1.33426939e-01 -3.15955239e-01  9.49575028e-02 -2.66853878e-01
   -3.15955239e-01 -1.33426939e-01]
  [ 5.03470377e-01 -3.52731535e-01  2.81567411e-02  1.00694075e+00
   -3.52731535e-01  5.03470377e-01]
  [ 3.93630714e-02  1.10770683e-01 -1.05673871e-02  7.87261428e-02
    1.10770683e-01  3.93630714e-02]
  [ 5.09580305e-01  4.65255489e-01  9.24843702e-02  1.01916061e+00
    4.65255489e-01  5.09580305e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.60534515e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-2.14278105e-01  1.63465065e-01 -1.03268419e-02 -4.28556210e-01
    1.63465065e-01 -2.14278105e-01]
  [-1.60981967e-01 -4.00490453e-01  7.54810651e-02 -3.21963935e-01
   -4.00490453e-01 -1.60981967e-01]
  [ 4.87746420e-01 -3.76014315e-01  2.30919334e-02  9.75492840e-01
   -3.76014315e-01  4.87746420e-01]
  [ 4.28733678e-02  1.15473583e-01 -6.63571668e-03  8.57467356e-02
    1.15473583e-01  4.28733678e-02]
  [ 6.05168647e-01  7.07226039e-01  1.23870915e-01  1.21033729e+00
    7.07226039e-01  6.05168647e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.46870655e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-2.05888039e-01  1.69308690e-01 -7.93085689e-03 -4.11776078e-01
    1.69308690e-01 -2.05888039e-01]
  [-1.84663808e-01 -4.65451966e-01  5.95026120e-02 -3.69327617e-01
   -4.65451966e-01 -1.84663808e-01]
  [ 4.66407067e-01 -3.87612082e-01  1.76410135e-02  9.32814134e-01
   -3.87612082e-01  4.66407067e-01]
  [ 4.52451106e-02  1.19865711e-01 -4.73313045e-03  9.04902212e-02
    1.19865711e-01  4.52451106e-02]
  [ 6.90798450e-01  9.20396632e-01  1.49475828e-01  1.38159690e+00
    9.20396632e-01  6.90798450e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.36280366e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-1.98803165e-01  1.73327269e-01 -6.03008194e-03 -3.97606330e-01
    1.73327269e-01 -1.98803165e-01]
  [-2.04303740e-01 -5.16111389e-01  4.68785776e-02 -4.08607480e-01
   -5.16111389e-01 -2.04303740e-01]
  [ 4.47070329e-01 -3.94304031e-01  1.32107443e-02  8.94140657e-01
   -3.94304031e-01  4.47070329e-01]
  [ 4.69732052e-02  1.22961726e-01 -3.35899414e-03  9.39464104e-02
    1.22961726e-01  4.69732052e-02]
  [ 7.68998997e-01  1.10844286e+00  1.70889328e-01  1.53799799e+00
    1.10844286e+00  7.68998997e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.27091251e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-1.92789113e-01  1.75978657e-01 -4.54517629e-03 -3.85578227e-01
    1.75978657e-01 -1.92789113e-01]
  [-2.20500139e-01 -5.55540706e-01  3.68776525e-02 -4.41000277e-01
   -5.55540706e-01 -2.20500139e-01]
  [ 4.30424856e-01 -3.97907707e-01  9.75257132e-03  8.60849711e-01
   -3.97907707e-01  4.30424856e-01]
  [ 4.82793653e-02  1.24952071e-01 -2.30991631e-03  9.65587306e-02
    1.24952071e-01  4.82793653e-02]
  [ 8.40805132e-01  1.27504628e+00  1.89020151e-01  1.68161026e+00
    1.27504628e+00  8.40805132e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.18989205e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-1.87672774e-01  1.77588334e-01 -3.38318223e-03 -3.75345548e-01
    1.77588334e-01 -1.87672774e-01]
  [-2.33807210e-01 -5.86081383e-01  2.89236334e-02 -4.67614420e-01
   -5.86081383e-01 -2.33807210e-01]
  [ 4.16201398e-01 -3.99295277e-01  7.06598589e-03  8.32402797e-01
   -3.99295277e-01  4.16201398e-01]
  [ 4.92546647e-02  1.26089711e-01 -1.50412008e-03  9.85093295e-02
    1.26089711e-01  4.92546647e-02]
  [ 9.06806543e-01  1.42334018e+00  2.04522708e-01  1.81361309e+00
    1.42334018e+00  9.06806543e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.11829985e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-1.83320440e-01  1.78410041e-01 -2.47240697e-03 -3.66640879e-01
    1.78410041e-01 -1.83320440e-01]
  [-2.44690164e-01 -6.09568484e-01  2.25774268e-02 -4.89380328e-01
   -6.09568484e-01 -2.44690164e-01]
  [ 4.04061655e-01 -3.99063011e-01  4.97908397e-03  8.08123310e-01
   -3.99063011e-01  4.04061655e-01]
  [ 4.99612483e-02  1.26581014e-01 -8.85891392e-04  9.99224967e-02
    1.26581014e-01  4.99612483e-02]
  [ 9.67473969e-01  1.55589415e+00  2.17895305e-01  1.93494794e+00
    1.55589415e+00  9.67473969e-01]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  5.05500234e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-1.79620591e-01  1.78640115e-01 -1.75822432e-03 -3.59241182e-01
    1.78640115e-01 -1.79620591e-01]
  [-2.53540124e-01 -6.27448860e-01  1.75019835e-02 -5.07080249e-01
   -6.27448860e-01 -2.53540124e-01]
  [ 3.93704969e-01 -3.97656642e-01  3.35895468e-03  7.87409939e-01
   -3.97656642e-01  3.93704969e-01]
  [ 5.04492283e-02  1.26586733e-01 -4.13401193e-04  1.00898457e-01
    1.26586733e-01  5.04492283e-02]
  [ 1.02322336e+00  1.67481439e+00  2.29524046e-01  2.04644672e+00
    1.67481439e+00  1.02322336e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  4.99901907e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]

 [[-1.76478441e-01  1.78430281e-01 -1.19867655e-03 -3.52956881e-01
    1.78430281e-01 -1.76478441e-01]
  [-2.60686972e-01 -6.40868688e-01  1.34365064e-02 -5.21373944e-01
   -6.40868688e-01 -2.60686972e-01]
  [ 3.84873835e-01 -3.95414932e-01  2.10369506e-03  7.69747669e-01
   -3.95414932e-01  3.84873835e-01]
  [ 5.07601806e-02  1.26231632e-01 -5.46464877e-05  1.01520361e-01
    1.26231632e-01  5.07601806e-02]
  [ 1.07443160e+00  1.78183962e+00  2.39710937e-01  2.14886320e+00
    1.78183962e+00  1.07443160e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  4.94949118e-01
    0.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    1.00000000e+00  0.00000000e+00]
  [ 0.00000000e+00  0.00000000e+00  0.00000000e+00  0.00000000e+00
    0.00000000e+00  0.00000000e+00]]]
     ssigmay:  [[[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. 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. 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. 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. 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. 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. 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. 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. 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. 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. 1.]]]
      sigmaz:  None
          sz:  None
         srz:  None
     ssigmaz:  None
        sllh:  [nan nan nan nan nan nan nan nan]
       s2llh:  None
      status:  0.0
        sres:  [[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. 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. 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. 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. 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. 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. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0.]
 [0. 0. 0. 0. 0. 0. 0. 0.]]

Adjoint sensitivity analysis


In [19]:
# Set model options
model = model_module.getModel()
p_orig = np.array(model.getParameters())
p_orig[list(model.getParameterIds()).index('observable_x1withsigma_sigma')] = 0.1 # Change default parameter
model.setParameters(p_orig.flatten())
model.setParameterScale(amici.ParameterScaling_none)
model.setTimepoints(np.linspace(0, 10, 21)) 

solver = model.getSolver()
solver.setMaxSteps(10**4) # Set maximum number of steps for the solver

# simulate time-course to get artificial data
rdata = amici.runAmiciSimulation(model, solver)
edata = amici.ExpData(rdata, 1.0, 0)
edata.fixedParameters = model.getFixedParameters()
# set sigma to 1.0 except for observable 5, so that p[7] is used instead
# (if we have sigma parameterized, the corresponding ExpData entries must NaN, otherwise they will override the parameter)
edata.setObservedDataStdDev(rdata['t']*0+np.nan,
                            list(model.getObservableIds()).index('observable_x1withsigma'))

# enable sensitivities
solver.setSensitivityOrder(amici.SensitivityOrder_first) # First-order ...
solver.setSensitivityMethod(amici.SensitivityMethod_adjoint)        # ... adjoint sensitivities
model.requireSensitivitiesForAllParameters()              # ... w.r.t. all parameters

# compute adjoint sensitivities
rdata = amici.runAmiciSimulation(model, solver, edata)
#print(rdata['sigmay'])
print('Log-likelihood: %f\nGradient: %s' % (rdata['llh'], rdata['sllh']))


Log-likelihood: -862.683106
Gradient: [-1.91102072e+01 -7.10905588e+00  3.96506777e+01  9.94087022e-01
  1.27425627e+01  5.63237660e+00  5.14250302e+00  1.43782082e+04]

Finite differences gradient check

Compare AMICI-computed gradient with finite differences


In [20]:
from scipy.optimize import check_grad

def func(x0, symbol='llh', x0full=None, plist=[], verbose=False):
    p = x0[:]
    if len(plist):
        p = x0full[:]
        p[plist] = x0
    verbose and print('f: p=%s' % p)
    
    old_parameters = model.getParameters()
    solver.setSensitivityOrder(amici.SensitivityOrder_none)
    model.setParameters(p)
    rdata = amici.runAmiciSimulation(model, solver, edata)
    
    model.setParameters(old_parameters)
    
    res = np.sum(rdata[symbol])
    verbose and print(res)
    return res

def grad(x0, symbol='llh', x0full=None, plist=[], verbose=False):
    p = x0[:]
    if len(plist):
        model.setParameterList(plist)
        p = x0full[:]
        p[plist] = x0
    else:
        model.requireSensitivitiesForAllParameters()
    verbose and print('g: p=%s' % p)
    
    old_parameters = model.getParameters()
    solver.setSensitivityMethod(amici.SensitivityMethod_forward)
    solver.setSensitivityOrder(amici.SensitivityOrder_first)
    model.setParameters(p)
    rdata = amici.runAmiciSimulation(model, solver, edata)
    
    model.setParameters(old_parameters)

    res = rdata['s%s' % symbol]
    if not isinstance(res, float):
        if len(res.shape) == 3:
            res = np.sum(res, axis=(0, 2))
    verbose and print(res)
    return res

err_norm = check_grad(func, grad, p_orig, 'llh')
print('sllh: |error|_2: %f' % err_norm)
# assert err_norm < 1e-6
print()

for ip in range(model.np()):
    plist = [ip]
    p = p_orig.copy()
    err_norm = check_grad(func, grad, p[plist], 'llh', p, [ip])
    print('sllh: p[%d]: |error|_2: %f' % (ip, err_norm))

print()
for ip in range(model.np()):
    plist = [ip]
    p = p_orig.copy()
    err_norm = check_grad(func, grad, p[plist], 'y', p, [ip])
    print('sy: p[%d]: |error|_2: %f' % (ip, err_norm))

print()
for ip in range(model.np()):
    plist = [ip]
    p = p_orig.copy()
    err_norm = check_grad(func, grad, p[plist], 'x', p, [ip])
    print('sx: p[%d]: |error|_2: %f' % (ip, err_norm))

print()
for ip in range(model.np()):
    plist = [ip]
    p = p_orig.copy()
    err_norm = check_grad(func, grad, p[plist], 'sigmay', p, [ip])
    print('ssigmay: p[%d]: |error|_2: %f' % (ip, err_norm))


sllh: |error|_2: 0.003233

sllh: p[0]: |error|_2: 0.000035
sllh: p[1]: |error|_2: 0.000048
sllh: p[2]: |error|_2: 0.000063
sllh: p[3]: |error|_2: 0.000031
sllh: p[4]: |error|_2: 0.000013
sllh: p[5]: |error|_2: 0.000009
sllh: p[6]: |error|_2: 0.000022
sllh: p[7]: |error|_2: 0.003258

sy: p[0]: |error|_2: 0.000001
sy: p[1]: |error|_2: 0.000000
sy: p[2]: |error|_2: 0.000001
sy: p[3]: |error|_2: 0.000006
sy: p[4]: |error|_2: 0.000003
sy: p[5]: |error|_2: 0.000001
sy: p[6]: |error|_2: 0.000000
sy: p[7]: |error|_2: 0.000000

sx: p[0]: |error|_2: 0.000001
sx: p[1]: |error|_2: 0.000001
sx: p[2]: |error|_2: 0.000001
sx: p[3]: |error|_2: 0.000000
sx: p[4]: |error|_2: 0.000001
sx: p[5]: |error|_2: 0.000000
sx: p[6]: |error|_2: 0.000000
sx: p[7]: |error|_2: 0.000000

ssigmay: p[0]: |error|_2: 0.000000
ssigmay: p[1]: |error|_2: 0.000000
ssigmay: p[2]: |error|_2: 0.000000
ssigmay: p[3]: |error|_2: 0.000000
ssigmay: p[4]: |error|_2: 0.000000
ssigmay: p[5]: |error|_2: 0.000000
ssigmay: p[6]: |error|_2: 0.000000
ssigmay: p[7]: |error|_2: 0.000000

Export as DataFrame

Experimental data and simulation results can both be exported as pandas Dataframe to allow for an easier inspection of numeric values


In [21]:
# run the simulation
rdata = amici.runAmiciSimulation(model, solver, edata)

In [22]:
# look at the ExpData as DataFrame 
df = amici.getDataObservablesAsDataFrame(model, [edata])
df


Out[22]:
time datatype t_presim k0 k0_preeq k0_presim observable_x1 observable_x2 observable_x3 observable_x1_scaled observable_x2_offsetted observable_x1withsigma observable_x1_std observable_x2_std observable_x3_std observable_x1_scaled_std observable_x2_offsetted_std observable_x1withsigma_std
0 0.0 data 0.0 1.0 NaN NaN 0.024431 0.169439 -0.112275 -0.156135 4.591891 -0.461117 1.0 1.0 1.0 1.0 1.0 NaN
1 0.5 data 0.0 1.0 NaN NaN -0.048291 -0.672858 -0.166120 -0.246694 4.145442 -0.182429 1.0 1.0 1.0 1.0 1.0 NaN
2 1.0 data 0.0 1.0 NaN NaN 1.539596 0.619225 1.307691 0.821159 4.584188 2.428329 1.0 1.0 1.0 1.0 1.0 NaN
3 1.5 data 0.0 1.0 NaN NaN 2.930428 1.082986 0.713558 2.934550 3.943787 -0.612084 1.0 1.0 1.0 1.0 1.0 NaN
4 2.0 data 0.0 1.0 NaN NaN 1.734583 1.858265 1.333083 1.228891 3.776846 0.853197 1.0 1.0 1.0 1.0 1.0 NaN
5 2.5 data 0.0 1.0 NaN NaN 1.491161 1.330449 -0.187078 3.317399 4.241527 0.384136 1.0 1.0 1.0 1.0 1.0 NaN
6 3.0 data 0.0 1.0 NaN NaN 1.267422 2.563478 -0.815568 0.569572 2.035637 -0.196271 1.0 1.0 1.0 1.0 1.0 NaN
7 3.5 data 0.0 1.0 NaN NaN 3.507608 1.227280 0.009140 2.158644 2.811899 0.192859 1.0 1.0 1.0 1.0 1.0 NaN
8 4.0 data 0.0 1.0 NaN NaN 0.889391 2.717031 -0.966199 3.871869 3.722095 1.507652 1.0 1.0 1.0 1.0 1.0 NaN
9 4.5 data 0.0 1.0 NaN NaN -0.087687 0.513524 -1.007304 1.053572 3.759065 0.752053 1.0 1.0 1.0 1.0 1.0 NaN
10 5.0 data 0.0 1.0 NaN NaN 3.526411 -1.430472 1.923501 0.418906 2.919740 2.236638 1.0 1.0 1.0 1.0 1.0 NaN
11 5.5 data 0.0 1.0 NaN NaN 0.365648 -0.399467 0.777448 3.405444 3.593012 0.931441 1.0 1.0 1.0 1.0 1.0 NaN
12 6.0 data 0.0 1.0 NaN NaN 1.588553 0.994930 1.690150 1.003723 2.909611 1.193227 1.0 1.0 1.0 1.0 1.0 NaN
13 6.5 data 0.0 1.0 NaN NaN 1.163786 0.735783 -0.667029 0.966362 5.729384 0.230227 1.0 1.0 1.0 1.0 1.0 NaN
14 7.0 data 0.0 1.0 NaN NaN 0.851443 -0.544395 0.058794 1.945118 3.200816 0.901044 1.0 1.0 1.0 1.0 1.0 NaN
15 7.5 data 0.0 1.0 NaN NaN 0.344363 1.609499 0.934498 1.108364 4.041441 -0.804747 1.0 1.0 1.0 1.0 1.0 NaN
16 8.0 data 0.0 1.0 NaN NaN 1.023536 0.326682 0.665340 1.914369 5.218419 -0.466391 1.0 1.0 1.0 1.0 1.0 NaN
17 8.5 data 0.0 1.0 NaN NaN -0.270099 -1.417993 1.169516 1.239814 3.182998 0.209709 1.0 1.0 1.0 1.0 1.0 NaN
18 9.0 data 0.0 1.0 NaN NaN 0.363280 0.719389 1.861780 -0.647324 4.869684 -0.206207 1.0 1.0 1.0 1.0 1.0 NaN
19 9.5 data 0.0 1.0 NaN NaN 1.723910 1.294062 -0.641339 2.746924 3.686706 0.677104 1.0 1.0 1.0 1.0 1.0 NaN
20 10.0 data 0.0 1.0 NaN NaN 0.216199 -0.937043 -1.225093 1.984324 4.109546 0.803328 1.0 1.0 1.0 1.0 1.0 NaN

In [23]:
# from the exported dataframe, we can actually reconstruct a copy of the ExpData instace
reconstructed_edata = amici.getEdataFromDataFrame(model, df)

In [24]:
# look at the States in rdata as DataFrame 
amici.getResidualsAsDataFrame(model, [edata], [rdata])


Out[24]:
time t_presim k0 k0_preeq k0_presim observable_x1 observable_x2 observable_x3 observable_x1_scaled observable_x2_offsetted observable_x1withsigma
0 0.0 0.0 1.0 NaN NaN 0.075569 0.230561 0.812275 0.356135 1.191891 5.611173
1 0.5 0.0 1.0 NaN NaN 0.587658 1.357536 0.357611 1.325429 0.460763 7.217960
2 1.0 0.0 1.0 NaN NaN 0.959523 0.114062 1.211267 0.338986 0.850900 18.482569
3 1.5 0.0 1.0 NaN NaN 2.360028 0.352334 0.637482 1.793751 0.213135 11.824829
4 2.0 0.0 1.0 NaN NaN 1.174049 1.142429 1.263389 0.107822 0.061010 2.926629
5 2.5 0.0 1.0 NaN NaN 0.938105 0.631698 0.253379 2.211288 0.542776 1.689200
6 3.0 0.0 1.0 NaN NaN 0.720551 1.881516 0.879301 0.524169 1.646325 7.431420
7 3.5 0.0 1.0 NaN NaN 2.966248 0.561171 0.052367 1.075924 0.854209 3.485010
8 4.0 0.0 1.0 NaN NaN 0.353110 2.065728 1.025694 2.799308 0.070792 9.713711
9 4.5 0.0 1.0 NaN NaN 0.619225 0.123990 1.064957 0.009505 0.121550 2.205146
10 5.0 0.0 1.0 NaN NaN 2.999320 2.055154 1.867541 0.635276 0.704941 17.095472
11 5.5 0.0 1.0 NaN NaN 0.157266 1.012200 0.723048 2.359615 0.019721 4.085265
12 6.0 0.0 1.0 NaN NaN 1.069564 0.393327 1.637190 0.034255 0.691992 6.742378
13 6.5 0.0 1.0 NaN NaN 0.648487 0.144554 0.718658 0.064236 2.138155 2.850718
14 7.0 0.0 1.0 NaN NaN 0.339613 1.125950 0.008395 0.921458 0.380739 3.892141
15 7.5 0.0 1.0 NaN NaN 0.164205 1.036970 0.885239 0.091229 0.468912 13.133151
16 8.0 0.0 1.0 NaN NaN 0.518036 0.237421 0.617136 0.903369 1.654316 9.718916
17 8.5 0.0 1.0 NaN NaN 0.772714 1.974227 1.122292 0.234584 0.373236 2.929067
18 9.0 0.0 1.0 NaN NaN 0.136622 0.170508 1.815465 1.647128 1.320802 7.061092
19 9.5 0.0 1.0 NaN NaN 1.226561 0.752054 0.686809 1.752224 0.144698 1.797547
20 10.0 0.0 1.0 NaN NaN 0.278751 1.472624 1.269778 0.994425 0.573965 3.083792

In [25]:
# look at the Observables in rdata as DataFrame 
amici.getSimulationObservablesAsDataFrame(model, [edata], [rdata])


Out[25]:
time datatype t_presim k0 k0_preeq k0_presim observable_x1 observable_x2 observable_x3 observable_x1_scaled observable_x2_offsetted observable_x1withsigma observable_x1_std observable_x2_std observable_x3_std observable_x1_scaled_std observable_x2_offsetted_std observable_x1withsigma_std
0 0.0 simulation 0.0 1.0 NaN NaN 0.100000 0.400000 0.700000 0.200000 3.400000 0.100000 1.0 1.0 1.0 1.0 1.0 0.1
1 0.5 simulation 0.0 1.0 NaN NaN 0.539367 0.684679 0.191491 1.078734 3.684679 0.539367 1.0 1.0 1.0 1.0 1.0 0.1
2 1.0 simulation 0.0 1.0 NaN NaN 0.580072 0.733287 0.096424 1.160145 3.733287 0.580072 1.0 1.0 1.0 1.0 1.0 0.1
3 1.5 simulation 0.0 1.0 NaN NaN 0.570399 0.730652 0.076076 1.140799 3.730652 0.570399 1.0 1.0 1.0 1.0 1.0 0.1
4 2.0 simulation 0.0 1.0 NaN NaN 0.560535 0.715836 0.069694 1.121069 3.715836 0.560535 1.0 1.0 1.0 1.0 1.0 0.1
5 2.5 simulation 0.0 1.0 NaN NaN 0.553056 0.698751 0.066301 1.106112 3.698751 0.553056 1.0 1.0 1.0 1.0 1.0 0.1
6 3.0 simulation 0.0 1.0 NaN NaN 0.546871 0.681963 0.063733 1.093741 3.681963 0.546871 1.0 1.0 1.0 1.0 1.0 0.1
7 3.5 simulation 0.0 1.0 NaN NaN 0.541360 0.666109 0.061506 1.082720 3.666109 0.541360 1.0 1.0 1.0 1.0 1.0 0.1
8 4.0 simulation 0.0 1.0 NaN NaN 0.536280 0.651302 0.059495 1.072561 3.651302 0.536280 1.0 1.0 1.0 1.0 1.0 0.1
9 4.5 simulation 0.0 1.0 NaN NaN 0.531538 0.637515 0.057653 1.063076 3.637515 0.531538 1.0 1.0 1.0 1.0 1.0 0.1
10 5.0 simulation 0.0 1.0 NaN NaN 0.527091 0.624681 0.055960 1.054183 3.624681 0.527091 1.0 1.0 1.0 1.0 1.0 0.1
11 5.5 simulation 0.0 1.0 NaN NaN 0.522915 0.612733 0.054400 1.045829 3.612733 0.522915 1.0 1.0 1.0 1.0 1.0 0.1
12 6.0 simulation 0.0 1.0 NaN NaN 0.518989 0.601603 0.052960 1.037978 3.601603 0.518989 1.0 1.0 1.0 1.0 1.0 0.1
13 6.5 simulation 0.0 1.0 NaN NaN 0.515299 0.591229 0.051629 1.030598 3.591229 0.515299 1.0 1.0 1.0 1.0 1.0 0.1
14 7.0 simulation 0.0 1.0 NaN NaN 0.511830 0.581555 0.050399 1.023660 3.581555 0.511830 1.0 1.0 1.0 1.0 1.0 0.1
15 7.5 simulation 0.0 1.0 NaN NaN 0.508568 0.572529 0.049259 1.017136 3.572529 0.508568 1.0 1.0 1.0 1.0 1.0 0.1
16 8.0 simulation 0.0 1.0 NaN NaN 0.505500 0.564103 0.048203 1.011000 3.564103 0.505500 1.0 1.0 1.0 1.0 1.0 0.1
17 8.5 simulation 0.0 1.0 NaN NaN 0.502615 0.556234 0.047224 1.005231 3.556234 0.502615 1.0 1.0 1.0 1.0 1.0 0.1
18 9.0 simulation 0.0 1.0 NaN NaN 0.499902 0.548881 0.046315 0.999804 3.548881 0.499902 1.0 1.0 1.0 1.0 1.0 0.1
19 9.5 simulation 0.0 1.0 NaN NaN 0.497350 0.542008 0.045471 0.994700 3.542008 0.497350 1.0 1.0 1.0 1.0 1.0 0.1
20 10.0 simulation 0.0 1.0 NaN NaN 0.494949 0.535581 0.044686 0.989898 3.535581 0.494949 1.0 1.0 1.0 1.0 1.0 0.1

In [26]:
# look at the States in rdata as DataFrame 
amici.getSimulationStatesAsDataFrame(model, [edata], [rdata])


Out[26]:
time t_presim k0 k0_preeq k0_presim x1 x2 x3
0 0.0 0.0 1.0 NaN NaN 0.100000 0.400000 0.700000
1 0.5 0.0 1.0 NaN NaN 0.539367 0.684679 0.191491
2 1.0 0.0 1.0 NaN NaN 0.580072 0.733287 0.096424
3 1.5 0.0 1.0 NaN NaN 0.570399 0.730652 0.076076
4 2.0 0.0 1.0 NaN NaN 0.560535 0.715836 0.069694
5 2.5 0.0 1.0 NaN NaN 0.553056 0.698751 0.066301
6 3.0 0.0 1.0 NaN NaN 0.546871 0.681963 0.063733
7 3.5 0.0 1.0 NaN NaN 0.541360 0.666109 0.061506
8 4.0 0.0 1.0 NaN NaN 0.536280 0.651302 0.059495
9 4.5 0.0 1.0 NaN NaN 0.531538 0.637515 0.057653
10 5.0 0.0 1.0 NaN NaN 0.527091 0.624681 0.055960
11 5.5 0.0 1.0 NaN NaN 0.522915 0.612733 0.054400
12 6.0 0.0 1.0 NaN NaN 0.518989 0.601603 0.052960
13 6.5 0.0 1.0 NaN NaN 0.515299 0.591229 0.051629
14 7.0 0.0 1.0 NaN NaN 0.511830 0.581555 0.050399
15 7.5 0.0 1.0 NaN NaN 0.508568 0.572529 0.049259
16 8.0 0.0 1.0 NaN NaN 0.505500 0.564103 0.048203
17 8.5 0.0 1.0 NaN NaN 0.502615 0.556234 0.047224
18 9.0 0.0 1.0 NaN NaN 0.499902 0.548881 0.046315
19 9.5 0.0 1.0 NaN NaN 0.497350 0.542008 0.045471
20 10.0 0.0 1.0 NaN NaN 0.494949 0.535581 0.044686

In [ ]:


In [ ]: