Using CellML models of cardiac cells in Chaste
- Using CellML models of cardiac cells in Chaste
This page contains some notes on generating Chaste code for cardiac cell models from CellML files, using PyCml, the Chaste CellML toolkit. PyCml is included with the Chaste source release, in the folder python/pycml.
The process is mostly automatic, although some of the options require some human intervention, and the CellML file itself may need to be annotated in order to allow successful conversion (for example, to indicate which variables represent the transmembrane potential and stimulus current). Many annotated CellML files may be found in the Chaste repository, either in the heart/src/odes/cellml folder, or the cellml project.
There are two main ways of using CellML within Chaste. If you are a cardiac executable user, then providing you compiled the executable from source yourself, CellML files may be loaded on the fly, as described below. If however you are using the source release of Chaste directly, the sections on using dynamically loaded CellML models and using CellML files as sources in Chaste will be relevant.
Use of CellML in the cardiac executable
As of release 2.0 of Chaste, the cardiac executable has gained the ability to automatically load cell models encoded as CellML files at run-time, rather than needing them to be incorporated within Chaste when it is compiled. In order to take advantage of this, you need (at present) to have built the executable from source yourself, as it uses your Chaste source tree to convert the CellML file into runnable code.
This change means that the definition of ionic models within the parameters file has changed. To specify one of the models included within Chaste, you now need to wrap it in a <Hardcoded> element, e.g.
<Hardcoded>FaberRudy2000</Hardcoded>
Specify a dynamically loaded model using the <Dynamic> element instead, e.g.
<Dynamic> <Path relative_to="chaste_source_root">heart/dynamic/libDynamicallyLoadableLr91.so</Path> </Dynamic>
The path may either point to a pre-compiled shared library (as in the example), or a CellML file. If the latter, a shared library will be created in the same folder as the CellML file, and loaded by the executable at run-time. See source:trunk/heart/test/data/xml/ChasteParametersFullFormat.xml for a full example parameters file.
In release 2.1 of Chaste, support for CellML has been extended to allow you to modify certain parameters of the models within user-defined regions of the tissue. This is incorporated within the cell heterogeneity support in the parameters file:
<CellHeterogeneity> <Location unit="cm"> <Cuboid> <LowerCoordinates x="-2.0" y="-0.1" z="-1.0"/> <UpperCoordinates x="-0.5" y="0.1" z="1.0"/> </Cuboid> </Location> <SetParameter name="example" value="0.0"/> <SetParameter name="example2" value="2.0"/> </CellHeterogeneity>
Parameters to be modified in this way must be annotated specially in the CellML file - you cannot just modify any variable in the model. Each variable must have a cmeta:id specified (this gives the name used in the SetParameter element) and an RDF annotation added, e.g.
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#example"> <modifiable-parameter xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">yes</modifiable-parameter> </rdf:Description> </rdf:RDF>
Such a block may be placed inside any CellML element, but placing it within or next to the annotated variable is recommended. See below for more about RDF annotations.
The older scale factor support, which was hardcoded for a few of the models shipped with Chaste, can also be implemented for your own models using named parameters. If you have a parameter named ScaleFactorIto this will be given the value specified in the ScaleFactorIto element, for example. Using SetParameter instead is encouraged, however, as the explicit scale factor elements will be removed in the future.
Variable annotation as a parameter or derived quantity is also important to support the OutputVariables functionality in the executable. Any variable thus annotated, or any state variable, may be specified to be included in the output data along with V and phi_e.
Options controlling the code generation process for the CellML file may be specified using a file-specific configuration file, as described below. This may be used to apply optimisations to generate faster code, for example.
Using dynamically loaded CellML models in Chaste tests or projects
Your own code built on Chaste may also use dynamically compiled and loaded CellML models. There is now a helper class CellMLLoader to simplify this process, but all the functionality it uses was already in release 3.0. See that class' documentation for details of usage, but note that it allows you to pass options for the helper script described below to the model loading process.
The cellml project contains many annotated CellML files which you may use in this fashion.
Using CellML files as sources directly in Chaste
The Chaste build system uses the helper script described below to generate source code from CellML files 'on the fly'. CellML files located within src or test folders will generate code that can be used like normal C++ classes. It can be very helpful to have access to a range of models, by using the Chaste cellml project. This project contains many of the common models, pre-annotated with metadata for automatic conversion, and tagging of many of the common parameters you may wish to vary/investigate. See the README.txt at the previous link for instructions on including this project within your own project.
By default, all possible variants will be supplied to the helper script, i.e. it will generate CVODE code if your build is using CVODE, and Backward Euler code if a .out file is available. This can be fine-tuned by providing a file-specific configuration file. For a CellML file called "MyModel.cellml", this is a file "MyModel-conf.xml" located in the same folder. It may contain any content allowable in the main PyCml configuration file (see below, especially the section on command-line arguments), and will override any options specified therein. See source:trunk/heart/src/odes/cellml/TenTusscher2006Epi-conf.xml for an example. You can now also provide project-wide arguments, by adding them to your project SConscript file. For example, insert:
# Change some flags just for this project env = SConsTools.CloneEnv(env) env['PYCML_EXTRA_ARGS'] = ['--expose-annotated-variables']
before the DoProjectSConscript() call, to generate cell models that have parameters available for all metadata tagged variables.
If the CellML file is placed within a dynamic folder (e.g. source:trunk/heart/dynamic), however, it will be compiled to a shared library suitable for loading dynamically by the executable. In this case, only a single output variant (.hpp & .cpp pair) can be generated, and so an error will be produced if you supply arguments to the helper script (via a model-specific config file) that imply generating multiple variants. This functionality is used by the executable itself as described above.
See the CellML files in source:trunk/heart/src/odes/cellml for examples.
Using the helper script
The script source:trunk/python/ConvertCellModel.py should ease the process of using PyCml. Run the script supplying CellML files as arguments. Extra options for PyCml can also be supplied. It will generate both normal and 'fully' optimised versions of models (i.e. with both partial evaluation and lookup tables done).
Examples:
./python/ConvertCellModel.py --assume-valid heart/src/odes/cellml/LuoRudy1991.cellml ./python/ConvertCellModel.py heart/src/odes/cellml/AnotherCellModel.cellml --no-member-vars
Some options are specified implicitly. At the time of writing these are --conf=config.xml --use-chaste-stimulus --convert-interfaces --row-lookup-method -c -o. You can use
./python/pycml/translate.py --help
to see a full list of options to PyCml itself (--Wu is possibly the most useful - it demotes dimensional inconsistency errors to warnings). The helper script also defines some options, primarily --cvode, --normal, and --opt (run ./python/ConvertCellModel.py -h to see the full list). These control whether to create code for solving with CVODE or Euler's method, and whether to use PyCml's optimisation techniques. It defaults to --normal --opt, and hence generates two sets of source code: one with optimisation, and one without.
Generating Analytic Jacobians for Backward Euler and CVODE solvers
In release 2.1 support was added for generating Backward Euler style cell models using the helper script. In release 3.2 support is to be added for generating 'Analytic Jacobians' for CVODE cells in the same way (CVODE will run fine with numerical Jacobians if this step is not done, but slightly faster with analytic). This is a multi-stage process, as the analysis uses Maple to perform symbolic differentiation. First run PyCml directly to generate a Maple script, then run Maple on this script (again you may need to add the --Wu option to PyCml):
./python/pycml/translate.py -J --conf=python/pycml/config.xml <cellml_file.cellml> maple -i <cellml_file.mpl> > <cellml_file.out>
Finally, place the resulting .out file adjacent to the CellML file. Supplying the --backward-euler option to the helper script will then generate a backward Euler class. As of r17984 this process should now be fully automatic.
The helper script will name output files and classes based on the input CellML file name. For example, with an input file 'MyModel.cellml' and all variations generated (i.e. --normal --opt --cvode --backward-euler) the following files will be produced, which define the relevant classes:
- MyModel.hpp, MyModel.cpp: class CellMyModelFromCellML with no optimisation
- MyModelOpt.hpp, MyModelOpt.cpp: class CellMyModelFromCellMLOpt PyCml optimisations
- MyModelCvode.hpp, MyModelCvode.cpp: class CellMyModelFromCellMLCvode which can be solved with CVODE
- MyModelCvodeOpt.hpp, MyModelCvodeOpt.cpp: class CellMyModelFromCellMLCvodeOpt which can be solved with CVODE and uses PyCml optimisations
- MyModelBackwardEuler.hpp, MyModelBackwardEuler.cpp: class CellMyModelFromCellMLBackwardEuler which is solved using backward Euler, and uses PyCml optimisations.
By default the helper script will use the PyCml shipped with Chaste. If, however, you have the environment variable PYCML_DIR set, the script will assume this points to an installation of PyCml that you want to use instead.
Model annotation with RDF
CellML files may include metadata through the use of RDF, the Resource Description Framework (see the CellML metadata specification for more information). PyCml makes use of several different annotations when generating C++ source code for Chaste.
Some annotations may be required to enable successful code generation. In particular, PyCml needs to know which variables represent the transmembrane potential and stimulus current, in order to link the models into the mono/bi-domain equations. (If the model does not have a stimulus current because it represents a self-excitatory cell, it must be annotated to indicate this.) Many common naming conventions are supported by the default configuration file (source:trunk/python/pycml/config.xml), but name annotations can be used for unrecognised cases. You may also need to annotate the membrane capacitance if its value is required for automatic units conversions (if the model uses amps as the dimensions of transmembrane currents).
All metadata annotations must occur within an <rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> element. Such an element may occur within any CellML element, so you may place annotations next to the variable being annotated, or in a single rdf:RDF at the top level of the model, as you prefer. Variables to be annotated must have a cmeta:id attribute, which is used to identify the variable in the annotation.
Each annotation takes the form of an rdf:Description element:
<rdf:Description rdf:about="#cmeta_id"> <!-- Annotation element goes here --> </rdf:Description>
where "cmeta_id" is replaced by the value of the cmeta:id attribute on the variable being annotated.
The following annotations are understood by PyCml. See the CellML files in source:trunk/heart/src/odes/cellml for examples of their use.
Standardised names
<bqbiol:is xmlns:bqbiol="http://biomodels.net/biology-qualifiers/" rdf:resource="https://chaste.comlab.ox.ac.uk/cellml/ns/oxford-metadata#standard-name"/>
This can be used to specify the 'standard' (to PyCml) name of this variable. Replace 'standard-name' with the name in question.
The full list of allowed names can be found within the cellml_metadata.py file in the PyCml source code.
Currently only a few of these names, shown below, have special meaning to the Chaste CellML translation process, but this will be expanded in the future.
- membrane_voltage - the transmembrane potential
- membrane_capacitance - this is useful if converting the model current units to those expected by Chaste requires knowledge of the cell's membrane capacitance
- cytosolic_calcium_concentration - if this variable is present, the generated class will define the GetIntracellularCalciumConcentration method
- membrane_stimulus_current - specifies which variable represents the intracellular stimulus current
- membrane_stimulus_current_amplitude, membrane_stimulus_current_duration, membrane_stimulus_current_period, membrane_stimulus_current_offset - if the stimulus coded in the model is a regular square wave, the relevant parameters may be annotated to allow PyCml to define the UseCellMLDefaultStimulus method. The presence of membrane_stimulus_current_offset is optional.
Supplying the --expose-annotated-variables options to PyCml will also annotate all non-state-variables, annotated in this way, as parameters or derived quantities, as appropriate. This will thus make them accessible via the GetAnyVariable and OutputVariables functionality.
Full example:
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#voltage"> <bqbiol:is xmlns:bqbiol="http://biomodels.net/biology-qualifiers/" rdf:resource="https://chaste.comlab.ox.ac.uk/cellml/ns/oxford-metadata#membrane_voltage"/> </rdf:Description> </rdf:RDF>
We also use these standardised names to interface between models and protocol definitions in the FunctionalCuration system.
Modifiable parameters
<modifiable-parameter xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">yes</modifiable-parameter>
This specifies that the variable should become a named parameter within the C++ class, with name given by its cmeta:id (unless it has a standardised name, in which case that is used). The variable can then be altered at run-time, via calls to the SetParameter method. Such parameters can also be set via the ChasteParameters XML file, as described above, and selected for output using the OutputVariables functionality.
Full example:
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#max_conductance"> <modifiable-parameter xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">yes</modifiable-parameter> </rdf:Description> </rdf:RDF>
Derived quantities
<derived-quantity xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">yes</derived-quantity>
This specifies that the variable should become a derived quantity within the C++ class, with name given by its cmeta:id (unless it has a standardised name, in which case that is used). Derived quantities may be computed using the ComputeDerivedQuantities and ComputeDerivedQuantitiesFromCurrentState methods, and selected for output using the OutputVariables functionality. This is useful for quantities which are not necessary for computing the time evolution of the model, but still of interest.
Full example:
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#ionic_current"> <derived-quantity xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">yes</derived-quantity> </rdf:Description> </rdf:RDF>
Variable range checking
Chaste has the ability to check at each time step that cell model variables have not gone out of an expected range. This is especially useful for gating variables, which must lie between 0 and 1, being probabilities. In order to take advantage of this, you must annotate your model to specify the expected range, using the pycml:range-low and pycml:range-high annotations, e.g.
<range-low xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">0</range-low> <range-high xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">1</range-high>
The range specified is treated as a closed interval.
Full example:
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#gating_variable"> <range-low xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">0</range-low> <range-high xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">1</range-high> </rdf:Description> </rdf:RDF>
Self-excitatory models
This stops PyCml looking for a stimulus current in the model, so that sino-atrial node models can be used in Chaste. In this case, you need to annotate the model itself, i.e. the cmeta:id used should be that on the model element. For example,
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#demir_model_1994"> <is-self-excitatory xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#">yes</is-self-excitatory> </rdf:Description> </rdf:RDF>
Named attributes
After release 2.1 Chaste has gained the ability to use real-valued named attributes associated with ODE classes. PyCml can use this functionality to transfer named attributes from a CellML into Chaste, so they are available for querying by other Chaste code. This will not typically be of use to executable users, but developers working with Chaste may find it useful. For example, such annotations could specify various "typical values" for the cell model, such as normal cycle length, or suggested forward Euler time step.
These are model annotations, i.e. the cmeta:id used should be that on the model element. For example,
<rdf:RDF xmlns:rdf="http://www.w3.org/1999/02/22-rdf-syntax-ns#"> <rdf:Description rdf:about="#luo_rudy_1991"> <named-attribute xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#"> <rdf:Description> <name rdf:datatype="http://www.w3.org/2000/10/XMLSchema#string">SuggestedCycleLength</name> <value rdf:datatype="http://www.w3.org/2000/10/XMLSchema#double">750</value> </rdf:Description> </named-attribute> <named-attribute xmlns="https://chaste.comlab.ox.ac.uk/cellml/ns/pycml#"> <rdf:Description> <name rdf:datatype="http://www.w3.org/2000/10/XMLSchema#string">SuggestedForwardEulerTimestep</name> <value rdf:datatype="http://www.w3.org/2000/10/XMLSchema#double">0.005</value> </rdf:Description> </named-attribute> </rdf:Description> </rdf:RDF>
The PyCml configuration file
Various parts of PyCml's functionality may not work properly without a configuration file being specified. Notably this is used to define typical names for the variables in the model which represent the transmembrane potential, stimulus current, and other transmembrane ionic currents. This removes the need to annotate these variables (as described above) in common cases. The configuration file also specifies which variable(s) are used to index lookup tables, and the ranges over which these vary. Model-specific configuration files can also be used to change the command-line options supplied to PyCml.
The configuration file is XML, and contains both general and model-specific settings. A sample file is provided with PyCml (source:trunk/python/pycml/config.xml). This may be modified to add new model-specific configuration, or change the global settings.
Further information on the configuration file is given in the documentation for ConfigurationStore.read_configuration_file in translate.py. The key settings are as follows.
Global options are specified in a global element, for example:
<global> <lookup_tables> <lookup_table> <!-- Only one var allowed here --> <var type='config-name'>transmembrane_potential</var> <min>-100.0001</min> <max>59.9999</max> <step>0.01</step> </lookup_table> </lookup_tables> <currents> <stimulus> <var type='oxmeta'>membrane_stimulus_current</var> <var type='name'>membrane,i_Stim</var> <var type='name'>membrane,I_stim</var> <var type='name'>cell,i_Stim</var> <var type='name'>membrane,i_stim</var> <var type='name'>membrane,I_st</var> <var type='name'>membrane,i_st</var> <var type='name'>membrane,i_pulse</var> </stimulus> <ionic_match> <!-- regexp on var names --> <var type='name'>membrane,i_.*</var> <var type='name'>cell,i_.*</var> <var type='name'>cell,xi.*</var> </ionic_match> <!-- Note that the stimulus current is automatically excluded from being an ionic current. Also note that there are implicit ^ and $ around the regexp. --> </currents> <transmembrane_potential> <var type='oxmeta'>membrane_voltage</var> <var type='name'>membrane,V</var> <var type='name'>cell,V</var> <var type='name'>membrane,E</var> </transmembrane_potential> <membrane_capacitance> <!-- TODO: say how this is used --> <var type='oxmeta'>membrane_capacitance</var> <var type='name'>membrane,Cm</var> <var type='name'>membrane,C</var> <var type='name'>membrane,C_m</var> <var type='name'>membrane,C_sc</var> </membrane_capacitance> </global>
There are 3 ways of specifying variables:
- By name (var type='name'). Variable names are given in full form, i.e. "component,variable". The component given should be that from which the variable is exported.
- By standardised name (var type='oxmeta'). Use the name from the oxmeta annotations. This is still experimental.
- By reference to a section of this config file (for use when defining lookup table keys) e.g. <var type='config-name'>transmembrane_potential</var>
Where multiple options are given, the first one that matches something in the input CellML model is used. Ionic currents are specified using a regular expression to match variable names.
Model-specific settings are given in for_model or for_models elements. They can contain the same options as in the global element. Models can be specified by name (in the case of for_model) or by id (in both cases). For example:
<for_model id="mahajan_shiferaw_2008_version01"> ... </for_model> <for_model name="fox_model_2002"> ... </for_model> <for_models> <ids> <id>ten_tusscher_model_2004_endo</id> <id>ten_tusscher_model_2004_M</id> <id>ten_tusscher_model_2004_epi</id> </ids> ... </for_models>
Supplying command-line arguments via the configuration file
The configuration file may also specify additional command line arguments for PyCml, by including a command_line_args element, e.g.
<command_line_args> <arg>--Wu</arg> <arg>--cvode</arg> <arg>--opt</arg> </command_line_args>
The arguments provided override those specified on the command line if there are any conflicts. Note that arguments for the helper script ConvertCellModel.py may also be provided, although an error will occur if such a configuration file is used when running translate.py directly.
Invoking PyCml directly
With recent advances in PyCml integration with Chaste, as described above, it should now rarely be necessary to run PyCml directly. However, if you do want to, this section outlines some common use cases.
A general point is that, with any of the transformations, you may need to include the --Wu flag if the model does not pass units validation. This demotes the error messages to warnings, allowing transformation to proceed (assuming no other errors were found).
Generating a standard AbstractCardiacCell subclass
This one is easy, and the most automatic.
./python/pycml/translate.py --conf=python/pycml/config.xml --use-chaste-stimulus --convert-interfaces <cellml_file>
Generating optimised AbstractCardiacCell subclasses
To apply both partial evaluation and lookup tables, use
./python/pycml/translate.py --conf=python/pycml/config.xml --use-chaste-stimulus --convert-interfaces -a -p -l <cellml_file>
Generating an AbstractBackwardEulerCardiacCell subclass
This is a multi-stage process, as the analysis uses Maple to perform symbolic differentiation. Maple can be found on linux.cs.ox.ac.uk.
./python/pycml/translate.py -J --conf=python/pycml/config.xml <cellml_file.cellml> maple -i <cellml_file.mpl> > <cellml_file.out> ./python/pycml/translate.py -j <cellml_file.out> --conf=python/pycml/config.xml --use-chaste-stimulus --convert-interfaces -a <cellml_file.cellml>
Generating an optimised AbstractBackwardEulerCardiacCell subclass
The following commands will apply all optimisations, including the use of backward Euler.
./python/pycml/translate.py -J --conf=python/pycml/config.xml <cellml_file.cellml> maple -i <cellml_file.mpl> > <cellml_file.out> ./python/pycml/translate.py -j <cellml_file.out> --conf=python/pycml/config.xml --use-chaste-stimulus --convert-interfaces -a -p -l <cellml_file.cellml>
Installing PyCml
See InstallPyCml.