Reproducible model development in the Cardiac Electrophysiology Web Lab

The modelling of the electrophysiology of cardiac cells is one of the most mature areas of systems biology. This extended concentration of research effort brings with it new challenges, foremost among which is that of choosing which of these models is most suitable for addressing a particular scientific question. In a previous paper, we presented our initial work in developing an online resource for the characterisation and comparison of electrophysiological cell models in a wide range of experimental scenarios. In that work, we described how we had developed a novel protocol language that allowed us to separate the details of the mathematical model (the majority of cardiac cell models take the form of ordinary differential equations) from the experimental protocol being simulated. We developed a fully-open online repository (which we termed the Cardiac Electrophysiology Web Lab) which allows users to store and compare the results of applying the same experimental protocol to competing models. In the current paper we describe the most recent and planned extensions of this work, focused on supporting the process of model building from experimental data. We outline the necessary work to develop a machine-readable language to describe the process of inferring parameters from wet lab datasets, and illustrate our approach through a detailed example of fitting a model of the hERG channel using experimental data. We conclude by discussing the future challenges in making further progress in this domain towards our goal of facilitating a fully reproducible approach to the development of cardiac cell models.


Introduction
: Replicability versus reproducibility in cardiac electrophysiology modelling studies. Here we list what provisions enable replication and reproduction of cardiac electrophysiology modelling studies, and also estimate how often they feature in published studies. Our plan is for the next version of the Web Lab to become the missing tool for making model development reproducible.
There are at least three aspects to replicability and reproducibility in computational models that 48 we would like to distinguish between. We have outlined these aspects in Table 1, and we discuss 49 the entries in each row below: proposed to explain similar cardiac phenomena. 120 Knowing the provenance of models is especially important in the field of cardiac models, which 121 are often chimeric, employing data from multiple sources due to the difficulty of gathering data 122 sufficient to constrain today's complex models from a single source. Cardiac modellers often 123 borrow parameters or entire subunits from other models that may have been constructed for  When there exist a large number of parameterisations that describe the data equally well, the 130 model is deemed to be unidentifiable, and the modeller may wish to consider either an alternative 131 simpler model formulation, alterations to experiments to provide more information on parame-132 ters, or both (Raue et al., 2011). Considering overly-simplistic objective functions in parameter 133 tuning, such as matching just biomarkers from a single action potential, may lead to unidentified 134 parameters that can cause models to yield unexpected and erroneous behaviour when tested un-135 der new contexts of use (Daly et al., 2017). While model identifiability is a recognised problem in 136 cardiac cell models, its assessment has yet to be adopted as a standard practice during modelling 137 studies, perhaps due in part to the competing methodologies for doing so (Milescu et al., 2005;  Finally, reproducible reporting standards are important for understanding the effects of uncer-143 tainty and variability in biological models. Biological data is invariably affected by many sources 144 of variation, such as measurement error, intrinsic variation (e.g. "beat-to-beat" variability be-145 tween recordings on a single cardiomyocyte), and extrinsic variation (variation between individ-146 uals in a sample, e.g. inter-cell variability in a sample of cells, or inter-individual variability in a 147 population). This leads to uncertainty about (or variation in) the optimal model parameters to 148 describe the experimental data, and characterisation and interpretation of this uncertainty can 149 give us insights into biological variation or even the suitability of a given model to explain data 150 . Cases of extreme variation in optimal model parameters may indicate 151 the unsuitability of the model to represent the system, as it reduces faith in a direct biological 152 interpretation of each parameter.

153
In previous work, we sought to address the first two kinds of reproducibility listed above in 154 cardiac modelling studies through the development of the Cardiac Electrophysiology Web Lab.

155
The Web Lab is an online resource for the specification, execution, comparison, and sharing of 156 simulation experiments and their results (Cooper et al., 2016). This Web Lab was built using 157 the functional curation paradigm (Cooper et al., 2011): models are specified using a standard 158 format such as CellML, protocols are specified in a custom language capable of expressing a vast 159 range of experimental procedures, and interactions between the two are mediated by a domain-160 specific ontology. This allows for modular interactions between models and protocols, allowing 161 for multiple models to be compared under a common protocol or vice-versa, extending the 162 capabilities of contemporary online tools for analysing/comparing cellular model predictions such as WholeCellSimDB (Karr et al., 2012(Karr et al., , 2014. The Web Lab additionally provides visualisation 164 tools to aid in these comparative studies, as well as visibility settings that allow models and 165 protocols to be developed in private before being shared with the community.

166
In this paper, we will discuss our plans for, and initial progress towards, integrating experimen-167 tal data and reproducible model development into the Cardiac Electrophysiology Web Lab. We analysis. In section 2 we describe the steps needed to create such a tool (which we will refer to as 172 WL2 ) from our original implementation (which we call WL1 ). Section 3 showcases a prototype  We now outline the steps required to establish an improved Web Lab, WL2. An overview of each 181 step and the new capabilities it facilitates is shown in Table 2.

182
Step 1 Adding annotated data Comparing arbitrary data sets Structured queries Step 2 Linking data to protocols Comparing experimental protocol results Documenting data provenance Step 3 Comparing data to predictions Checking model applicability Documenting model provenance Continuous testing of models Step  pre-processing has been performed, and (3) code to reproduce the pre-processing process, to be 205 run off-line. This set-up would allow pre-processing code to be inspected, reviewed, and re-used.

206
It is also worth noting that some types of filtering (e.g. fitting a straight line through noisy data) 207 presuppose a certain structure in the data, and so mix modelling with pre-processing. By having 208 the raw data available online such work could be accommodated. Step 2: Linking data to experimental protocols 210 A crucial step in using data sets on WL2 will be to link them to an experimental protocol.  post-processing steps required to re-obtain that data. This will also facilitate a more careful 229 comparison of different data sets describing the outcomes of similar protocols.   Figure 1: A schematic overview of WL2. Experimental protocols, applied to biological models (e.g. myocytes, expression systems) give rise to experimental results. The association with a protocol, in combination with additional metadata, provides users with a thorough overview of how the data was obtained. Applied to computational models, the same protocols provide predictions. As in WL1, protocols are written in such a manner that they can be applied to several models on the Web Lab, and their predictions can be compared. A new feature will be the ability to compare predictions to predictions, experimental results to results, and results to predictions. By comparing experimental results and predictions from the same protocol, a fitting process can be initiated, leading to a set of parameter values represented either as singular points (optimisation) or distributions (inference).
2016b). The application of these techniques is a first step towards untangling experimental error  3. Prototype implementation 331 We now discuss the implementation of a prototype WL2, focused on performing statistical infer-332 ence over parameters of single-cell EP models and sub-models (e.g. ionic currents), and demon- 338  The rates are voltage-dependent functions, each parameterised by two scalar values as follows:

352
A Bayesian inference method was then applied to find the parameter values that provided the

379
In addition to the simulation specification, the WL2 prototype requires two files to complete the 380 specification of a fitting problem: a file containing experimental data, and a fitting specification 381 that shows how this data is employed to constrain the model. The content of these files will be 382 discussed in the subsequent sections. where the first row specifies names for each variable, and associated columns specify the corre-393 sponding data. We note that this structure currently expects zero-or one-dimensional data for 394 each named variable (although higher-dimensional arrays may be specified in flattened form), as 395 this was sufficient for our test case, but that the exact data representation can easily be changed 396 at a later stage. In the hERG current fitting experiment, the data file contains two columns of 397 equal length representing a series of (time, current) pairs. Our current implementation requires  The final component of a parameter fitting experiment is the fitting specification, which makes 403 use of a custom language that we will introduce in this section by means of a working example.

404
The fitting specification takes the form of a JSON-formatted text file (http://www.json.org).  Table 3, and will discuss the interpretation of each 410 required entry below.

411
Fitting specification entity Value algorithm AdaptiveMCMC arguments cmaOpt=5, cmaMaxFevals=20000, burn=50000, numIters=100000 output exp IKr=IKr input exp times=t prior (see Table 4) Table 3: Entries in the fitting specification for the hERG ion channel model. The value associated with the "algorithm" entry is a string of characters, and is represented as is, while all other value entries are nested JSON objects, and are presented in "key=value" format for clarity. This is also true for the prior specification, which is represented separately in Table 4 due to its size.
The first entry in the fitting specification is the "algorithm" to be used for parameter fitting, and approximate Bayesian inference algorithms that we are considering for inclusion in WL2.

418
Once we refine the list of algorithms we support, we will lobby for their inclusion in KiSAO (or 419 another accepted ontology) and adapt to this new form of algorithmic specification in future 420 iterations of the Web Lab.

421
In our prototype implementation, the adaptive MCMC algorithm uses a Gaussian likelihood 422 function, which is commonly assumed for time-series data. However, the prototype could easily 423 be extended to allow users to specify different likelihood functions, by adding an "objective" 424 entry to the fitting specification language.

425
The next entry we consider is a dictionary of "arguments," specific to the chosen fitting algorithm.

426
In the example shown in Table 3, these include the standard arguments for MCMC -the total 427 number of iterations "numIters" and the number of iterations discarded as burn-in "burn" -428 as well as two new arguments "cmaOpt" and "cmaMaxFevals" which deal with the selection of a

455
Parameter Name Prior Range Uniform(1e-7, 0.1) G Kr Uniform(0.0612, 0.612) obj:std 0.00463 Table 4: Prior distribution specified within the fitting specification for the 9-parameter hERG model. This prior is adapted from Beattie et al. (2018), who employ a wider prior in their MCMC inference but define this region as most likely to contain the optimal parameters. Parameters respect the shortened naming conventions of Equations (1)-(9) for clarity. An additional parameter, "obj:std", controls the observation noise standard deviation, part of the Gaussian likelihood function, and is set to a fixed value in this example (although in general it could be learned too).

456
We now present the results of our test-case fitting experiment, implemented in a WL2 prototype.

457
As with the WL1 implementation, an experiment may be carried out, and its results viewed, by 458 matching a model to a protocol in the 'Experiments' matrix view. Within the prototype WL2, 459 the only change to this set-up is that a 'protocol' entry now encompasses a simulation protocol, 460 fitting specification, and data file (as described in Section 3). The files used to represent this 461 fitting experiment are described in Section 3 (particularly Tables 3 and 4). Further details, 462 including links to the relevant online resources, can be found in Daly (2018).

463
After the execution of a fitting experiment, the first thing that the WL2 prototype allows us to 464 do is to compare the data simulated using our inferred maximum likelihood parameters to the 465 experimental data we employed during fitting. In Fig. 2, we see the results of overlaying data 466 simulated using the maximum posterior density parameters returned by our MCMC inference 467 onto the experimental data used to obtain these fitted parameter values. As the MCMC algorithm 468 returns a sample of parameter sets approximating the true distribution over parameters given 469 data, this visual shows us how well the most likely parameter set captures the observed behaviour.

470
The close agreement between these traces suggests that the inference strategy has produced a 471 distribution over parameterised models that captures observed behaviour well, which mirrors the  Clearly, there has not been sufficient incentive for scientists to provide annotated models and/or 507 data sets. However, very few tools make use of model or data annotations, so that the updated 508 Web Lab's annotation use could serve as a prime example of the benefits of annotation. As soon 509 as a link is established between experimental data and the models that use it (section 2.2), and 510 refitting (parts of) models becomes a routine task, questions of model-data provenance will arise 511 more naturally, and the need for well-annotated data will be felt by a wider audience. Comparing complex data sets. Biological systems are irreducibly complex; even when only a single 531 ionic current is measured, the 'background' is a living cell in which thousands of dynamical 532 mechanisms interact to create and maintain homoeostasis. As a result, any two independent 533 investigations into the same phenomenon will almost invariably differ in some details, some of 534 which may later turn out to be important. For annotation, this means that even when using a 535 standard such as the MICEE draft -which specifies around 100 properties to be recorded as 'minimum information' -some details will go unrecorded. It also means that the question of 537 whether two data sets describe 'the same' experiment is not always easy to decide. -which is a variable named in our ontology -could be tagged with the property "is a parameter 582 for the fast sodium current". Further properties could provide a more detailed description, e.g.

583
"influence fast sodium current activation", or "appears in a reaction rate of the form ae bV ". A 584 fitting algorithm for current "x" could then gather all variables tagged as "parameter for x" and At present developing cardiac electrophysiology models is "something of an art" (anonymous 604 senior cardiac modeller). We would like to see it become a science, defined by an unambiguous 605 algorithmic procedure. Our hope is that in the future a resource such as the Web Lab will 606 provide researchers with everything they need to know to reproduce a model's development.

607
The Web Lab will list what experimental protocols need to be performed in a wet lab in order 608 to parameterise a model, and then receive data from these experiments and use them to produce 609 parameterised mathematical models. Further in the future we aim to automate the process of  In this article we have described our work to date in developing a community-based resource to 619 support the cardiac electrophysiology research community. Our goal is that this resource should 620 become a repository for all aspects of the research of this community: experimental data and protocols, the computational models that are derived from that data and aid in its physiological 622 interpretation, and the instantiations in software of statistical inference techniques that are used 623 to derive those computational models from the experimental data and protocols. Whilst our work 624 as described here focuses on cardiac electrophysiology, the need that we address and the approach 625 that we have used are applicable across a large swathe of scientific research endeavour. In order 626 to make our work as widely accessible as possible, and in the hope that this approach might be