Reproducibility study of the Fabbri et al. 2017 model of the human sinus node action potential

ORIGINAL Abstract The sinoatrial node (SAN) is the natural pacemaker of the mammalian heart. It has been the subject of several mathematical studies, aimed at reproducing its electrical response under normal sinus rhythms, as well as under various conditions. Such studies were traditionally done using data from rabbit SAN cells. More recently, human SAN cell data have become available, resulting in the publication of a human SAN cell model (Fabbri et al., 2017), along with its CellML version. Here, we used the CellML ﬁle provided by the model authors, together with some SED-ML ﬁles and Python scripts that we created to reproduce the main results of the aforementioned modelling study.


Introduction
The sinoatrial node (SAN) plays an important role in cardiac function and although it has been extensively studied, some of its intrinsic mechanisms are still open for debate. Most of the experimental data used in SAN modelling have been carried out on animals and on rabbits in particular (Lakatta and DiFrancesco, 2009;DiFrancesco, 2010;Lakatta, 2010;Maltsev and Lakatta, 2010;Noble et al., 2010;Verkerk et al., 2007;Himeno et al., 2008;DiFrancesco and Noble, 2012;Lakatta and Maltsev, 2012;Rosen et al., 2012;Monfredi et al., 2013;Yaniv et al., 2013Yaniv et al., , 2015. This has resulted in the development of comprehensive SAN models (Wilders, 2007). Yet, this body of work can hardly be transposed to humans.
Human SAN action potentials were first recorded by Drouin (1997), followed by Verkerk et al. a decade later (Verkerk et al., 2007(Verkerk et al., , 2013. The first human SAN cell model was developed as a proof-of-concept by Seemann et al. (2006). This was followed by the model of Chandler et al. (2009). More recently, Verkerk and Wilders (2015) studied the effect of mutations on human SAN cells, highlighting the need for a human-specific SAN cellular electrophysiology model. Such a model was formulated by Pohl et al. (2016), but its action potential shape does not match that of experimental recordings. Fabbri et al. (2017) addressed this shortcoming by developing their human SAN cell model using available human electrophysiological data. Fabbri et al. (2017) published a CellML version (Cuellar et al., 2003) of their model on the Physiome Model Repository (Yu et al., 2011). However, the CellML file on its own is not sufficient to reproduce all predictions presented in the primary paper. Some SED-ML files (Waltemath et al., 2011) and Python scripts were therefore created and used with the aforementioned CellML file to reproduce the main results from Fabbri et al. (2017). No modifications were made to the CellML file mathematics or parameters and all the equations and parameters can be found in the original paper. Fabbri et al. (2017) developed a human SAN cell model, based on the rabbit SAN cell model of Severi et al. (2012) and on recent electrophysiological data from human SAN cells. The resulting action potential and calcium transient are in agreement with experimentally recorded values. Mutations associated with sinus node dysfunction were also modelled and their effects on pacing rate agree with clinical observations.

Model description
The model was developed in Simulink and simulations performed using MATLAB's ode15s solver (Shampine and Reichelt, 1997). Simulations were run until calcium dynamics reached steady-state. Custom MATLAB (2013a) code was used for automatic optimization and feature extraction. A CellML-encoded version of the model is available at https://models.physiomeproject.org/e/568.
The simulation results presented here were produced using the 2021-07-09 snapshot of OpenCOR (Garny and Hunter, 2015) together with various Python scripts that rely on a SED-ML file to configure (the solver to use, the duration of the simulation, the model parameters to track, etc.) and run a given simulation using the model encoded in CellML. Python scripts are also used to generate the figures using Matplotlib (Hunter, 2007).   Figure 2 of the primary paper with a cycle length of 814 ms which corresponds to a beating rate of 74 beats min -1 . In Figure 2, the time course of the simulated action potential (A) and its underlying components (B-K) are shown, which corresponds to Figure  3 in the primary paper.  Figure 3 shows the membrane potential and its associated currents during diastolic depolarization. The resulting behaviour corresponds to that presented in Figure 4 in the primary paper.

Model results
In Figure 4, the model is used to reproduce the results of a progressive block of the funny current, I f , to show its effect on the cycle length (CL) and, therefore, on the pacing rate. This corresponds to Figure 7 in the primary paper.   Figure 5 shows the effect of a shift in the activation curve of I f on the CL (A), the diastolic depolarization rate over the first 100 ms of diastolic depolarization (DDR 100 ; B), the maximum diastolic potential (MDP; C), and the action potential duration at 90% repolarization (APD 90 ; D). This corresponds to Figure 8 in the primary paper. As reported in Fabbri et al. (2017), little variation in APD 90 is observed, albeit with slight differences between the Python implementation presented here and the original MATLAB implementation used in the primary paper. The contribution of the sodium-calcium exchanger, I NaCa , to the action potential is illustrated in Figure 6, where a progressive block of I NaCa was performed. This corresponds to Figure 9 in the primary paper, although the results are slightly different.  Figure 7 shows the effect of acetylcholine (ACh) and isoprenaline (Iso) on the membrane potential, its net current, the respective target currents for ACh and Iso, and on the sarco-endoplasmic reticulum Ca 2+ -ATPase (SERCA) pump uptake rate. This corresponds to Figure 10 in the primary paper.

7/9
There are details/documentation on how the source code was compiled There are details on how to run the code in the provided documentation The initial conditions are provided for each of the simulations Details for creating reported graphical results from the simulation results Source code: a declarative language is used (e.g. SBML, CellML, NeuroML) The algorithms used are defined or cited in previous articles The algorithm parameters are defined Post-processing of the results are described in sufficient detail

Executable model provided:
The model is executable without source (e.g. desktop application, compiled code, online service) There are sufficient details to repeat the required simulation experiments