Considering discrepancy when calibrating a mechanistic electrophysiology model

Figure 1.. A comparison of the ten Tusscher (Model T [28], blue) and Fink (Model F [29], green) kinetics. These currents are voltage-clamp simulations under the same action potential clamp (shown in the top row panels). Only those currents with different kinetics are shown; the kinetics of INa, INaCa and INaK are identical in both models. Two of the gates in ICaL are identical in the two models, one gate has a different formulation, and Model F has one extra gate compared to Model T. The two models use different formulations for IKr (IKr activates during depolarization in Model T but not Model F), different parametrizations of the kinetics for IKs and Ito, and different equations for IK1 steady state. Currents are normalized in this plot by minimizing the squared-difference between the two models’ currents such that we emphasize the differences in kinetics rather than the conductances (which are rescaled during the calibration). Only ICaL shows what we would typically consider to be a large difference in repolarization kinetics, with the rest of the currents apparently being close matches between Model T and Model F. (Online version in colour.)

Figure 2.. Model F fitting and validation results. (a) Model F is fitted to the synthetic data (generated from Model T), using five action potentials recorded under a 1 Hz pacing protocol. The calibrated Model F (blue dashed line) shows an excellent fit to the calibration data (grey solid line). (b,c) Model F predctions for validation and context of use (CoU) data. (b) The calibrated Model F predictions closely matches the validation data (2 Hz pacing), giving a (false) confidence in the model performance. (c) Notably, Model F gives catastrophic predictions for the IKr block (CoU) experiments (suggesting the validation data are not an appropriate test given the intended model use). The posterior predictions are model predictions made using parameter values sampled from the posterior distribution (figure 3); here, 200 samples/predictions are shown, but they overlay and are not distinguishable by eye. (Online version in colour.)

Figure 3.. Marginals of the posterior distribution of the Model F parameters, in terms of scaling factors for the conductances in Model T (si = gModel~Fi/gModel~Ti). Values of 1 would represent the parameters of Model T that generated the data; note that none of the inferred parameters for Model F are close 1. The red dashed lines indicate the result of the global optimization routine. Two of these parameters, SKs and SNaK, have distributions hitting the lower bound that was imposed by the prior, indicating that the calibration process is attempting to make them smaller than 10% of the original Model F parameter values. (Online version in colour.)

Figure 4.. A cartoon to illustrate the effect of model discrepancy on parameter fits in different models. Each cloud represents a range of possible outputs from each model, which they can reach with different parameter values. The true data generating process (DGP) lies outside either of our imperfect model classes 1 and 2, and neither can fit the data perfectly due to model discrepancy. When we attempt to infer parameters, we will converge upon models that generate outputs closest to the true DGP under the constraint of being in each model. Adding more data just increases the confidence in being constrained to model parameterizations on the boundary of the particular model, i.e. we become certain about θ∗, the pseudo-true parameter value for each model. Note that different models will have different pseudo-true parameter values. (Online version in colour.)

Figure 5.. Markov model representation of Models A, B and C used in the ion channel model tutorial where Model C is taken as ground truth and used to generate synthetic data, while Models A and B are candidate models that we attempt to fit and use for predictions, demonstrating the challenge of both model discrepancy and model selection. (Online version in colour.)

Figure 6.. (a) The fitted model predictions for Models A (blue) and B (orange) for the ion channel example (i.e. model predictions for the data they were trained with). Both models have been fitted to synthetic calibration data (grey) generated using Model C using the sinusoidal voltage-clamp protocol [43]. (b) Models A (blue) and B (orange) predictions for the validation data (grey) generated from Model C under the staircase protocol [26] (not used in training). Note that there are significant discrepancies around 12 000 ms. (Online version in colour.)

Figure 7.. Model A inferred marginal posterior distributions for the conductance, g in equation (3.3), and kinetic parameters p1, …, p8 (a list of parameters referring to Ai,j and Bi,j in equation (3.5)) with different discrepancy models: i.i.d. noise (blue), GP(t) (orange), GP(O, V) (green) and ARMA(2, 2) (red). (Online version in colour.)

Figure 8.. Model A fitted to the sinusoidal calibration protocol using the different discrepancy models: i.i.d. noise, GP(t), GP(O, V) and ARMA(2, 2). The plots show the mean (solid lines) and 95% credible intervals (shaded) of the posterior prediction for each model. (Online version in colour.)

Figure 9.. Model A’s prediction using the discrepancy models (i.i.d. noise, GP(t), GP(O, V) and ARMA(2, 2)), trained using the staircase voltage-clamp protocol [26]. We plot the posterior predictive mean (solid lines) with 95% credible intervals (shaded). The red arrows point to the tail current after the two activation steps, and mark an area of visible model mismatch: note the different performance of the four discrepancy models in this region. The blue arrow points to an obvious artefact at approximately 7000 ms induced by the GP(t) prediction which was trained on the sinusoidal protocol, and which does not take into account that we are now predicting for the staircase protocol. (Online version in colour.)