Cell to whole organ global sensitivity analysis on a four-chamber heart electromechanics model using Gaussian processes emulators.

Cardiac pump function arises from a series of highly orchestrated events across multiple scales. Computational electromechanics can encode these events in physics-constrained models. However, the large number of parameters in these models has made the systematic study of the link between cellular, tissue, and organ scale parameters to whole heart physiology challenging. A patient-specific anatomical heart model, or digital twin, was created. Cellular ionic dynamics and contraction were simulated with the Courtemanche-Land and the ToR-ORd-Land models for the atria and the ventricles, respectively. Whole heart contraction was coupled with the circulatory system, simulated with CircAdapt, while accounting for the effect of the pericardium on cardiac motion. The four-chamber electromechanics framework resulted in 117 parameters of interest. The model was broken into five hierarchical sub-models: tissue electrophysiology, ToR-ORd-Land model, Courtemanche-Land model, passive mechanics and CircAdapt. For each sub-model, we trained Gaussian processes emulators (GPEs) that were then used to perform a global sensitivity analysis (GSA) to retain parameters explaining 90% of the total sensitivity for subsequent analysis. We identified 45 out of 117 parameters that were important for whole heart function. We performed a GSA over these 45 parameters and identified the systemic and pulmonary peripheral resistance as being critical parameters for a wide range of volumetric and hemodynamic cardiac indexes across all four chambers. We have shown that GPEs provide a robust method for mapping between cellular properties and clinical measurements. This could be applied to identify parameters that can be calibrated in patient-specific models or digital twins, and to link cellular function to clinical indexes.

Modelling the interaction between stem cells derived cardiomyocytes patches and host myocardium to aid non-arrhythmic engineered heart tissue design.

Application of epicardial patches constructed from human-induced pluripotent stem cell- derived cardiomyocytes (hiPSC-CMs) has been proposed as a long-term therapy to treat scarred hearts post myocardial infarction (MI). Understanding electrical interaction between engineered heart tissue patches (EHT) and host myocardium represents a key step toward a successful patch engraftment. EHT retain different electrical properties with respect to the host heart tissue due to the hiPSC-CMs immature phenotype, which may lead to increased arrhythmia risk. We developed a modelling framework to examine the influence of patch design on electrical activation at the engraftment site. We performed an in silico investigation of different patch design approaches to restore pre-MI activation properties and evaluated the associated arrhythmic risk. We developed an in silico cardiac electrophysiology model of a transmural cross section of host myocardium. The model featured an infarct region, an epicardial patch spanning the infarct region and a bath region. The patch is modelled as a layer of hiPSC-CM, combined with a layer of conductive polymer (CP). Tissue and patch geometrical dimensions and conductivities were incorporated through 10 modifiable model parameters. We validated our model against 4 independent experimental studies and showed that it can qualitatively reproduce their findings. We performed a global sensitivity analysis (GSA) to isolate the most important parameters, showing that the stimulus propagation is mainly governed by the scar depth, radius and conductivity when the scar is not transmural, and by the EHT patch conductivity when the scar is transmural. We assessed the relevance of small animal studies to humans by comparing simulations of rat, rabbit and human myocardium. We found that stimulus propagation paths and GSA sensitivity indices are consistent across species. We explored which EHT design variables have the potential to restore physiological propagation. Simulations predict that increasing EHT conductivity from 0.28 to 1-1.1 S/m recovered physiological activation in rat, rabbit and human. Finally, we assessed arrhythmia risk related to increasing EHT conductivity and tested increasing the EHT Na+ channel density as an alternative strategy to match healthy activation. Our results revealed a greater arrhythmia risk linked to increased EHT conductivity compared to increased Na+ channel density. We demonstrated that our modeling framework could capture the interaction between host and EHT patches observed in in vitro experiments. We showed that large (patch and tissue dimensions) and small (cardiac myocyte electrophysiology) scale differences between small animals and humans do not alter EHT patch effect on infarcted tissue. Our model revealed that only when the scar is transmural do EHT properties impact activation times and isolated the EHT conductivity as the main parameter influencing propagation. We predicted that restoring physiological activation by tuning EHT conductivity is possible but may promote arrhythmic behavior. Finally, our model suggests that acting on hiPSC-CMs low action potential upstroke velocity and lack of IK1 may restore pre-MI activation while not promoting arrhythmia.

Correction: Linking statistical shape models and simulated function in the healthy adult human heart.

[This corrects the article DOI: 10.1371/journal.pcbi.1008851.].

Quantitative mapping of force-pCa curves to whole-heart contraction and relaxation.

The force-pCa (F-pCa) curve is used to characterize steady-state contractile properties of cardiac muscle cells in different physiological, pathological and pharmacological conditions. This provides a reduced preparation in which to isolate sarcomere mechanisms. However, it is unclear how changes in the F-pCa curve impact emergent whole-heart mechanics quantitatively. We study the link between sarcomere and whole-heart function using a multiscale mathematical model of rat biventricular mechanics that describes sarcomere, tissue, anatomy, preload and afterload properties quantitatively. We first map individual cell-level changes in sarcomere-regulating parameters to organ-level changes in the left ventricular function described by pressure-volume loop characteristics (e.g. end-diastolic and end-systolic volumes, ejection fraction and isovolumetric relaxation time). We next map changes in the sarcomere-regulating parameters to changes in the F-pCa curve. We demonstrate that a change in the F-pCa curve can be caused by multiple different changes in sarcomere properties. We demonstrate that changes in sarcomere properties cause non-linear and, importantly, non-monotonic changes in left ventricular function. As a result, a change in sarcomere properties yielding changes in the F-pCa curve that improve contractility does not guarantee an improvement in whole-heart function. Likewise, a desired change in whole-heart function (i.e. ejection fraction or relaxation time) is not caused by a unique shift in the F-pCa curve. Changes in the F-pCa curve alone cannot be used to predict the impact of a compound on whole-heart function. KEY POINTS: The force-pCa (F-pCa) curve is used to assess myofilament calcium sensitivity after pharmacological modulation and to infer pharmacological effects on whole-heart function. We demonstrate that there is a non-unique mapping from changes in F-pCa curves to changes in left ventricular (LV) function. The effect of changes in F-pCa on LV function depend on the state of the heart and could be different for different pathological conditions. Screening of compounds to impact whole-heart function by F-pCa should be combined with active tension and calcium transient measurements to predict better how changes in muscle function will impact whole-heart physiology.

In silico identification of potential calcium dynamics and sarcomere targets for recovering left ventricular function in rat heart failure with preserved ejection fraction.

Heart failure with preserved ejection fraction (HFpEF) is a complex disease associated with multiple co-morbidities, where impaired cardiac mechanics are often the end effect. At the cellular level, cardiac mechanics can be pharmacologically manipulated by altering calcium signalling and the sarcomere. However, the link between cellular level modulations and whole organ pump function is incompletely understood. Our goal is to develop and use a multi-scale computational cardiac mechanics model of the obese ZSF1 HFpEF rat to identify important biomechanical mechanisms that underpin impaired cardiac function and to predict how whole-heart mechanical function can be recovered through altering cellular calcium dynamics and/or cellular contraction. The rat heart was modelled using a 3D biventricular biomechanics model. Biomechanics were described by 16 parameters, corresponding to intracellular calcium transient, sarcomere dynamics, cardiac tissue and hemodynamics properties. The model simulated left ventricular (LV) pressure-volume loops that were described by 14 scalar features. We trained a Gaussian process emulator to map the 16 input parameters to each of the 14 outputs. A global sensitivity analysis was performed, and identified calcium dynamics and thin and thick filament kinetics as key determinants of the organ scale pump function. We employed Bayesian history matching to build a model of the ZSF1 rat heart. Next, we recovered the LV function, described by ejection fraction, peak pressure, maximum rate of pressure rise and isovolumetric relaxation time constant. We found that by manipulating calcium, thin and thick filament properties we can recover 34%, 28% and 24% of the LV function in the ZSF1 rat heart, respectively, and 39% if we manipulate all of them together. We demonstrated how a combination of biophysically based models and their derived emulators can be used to identify potential pharmacological targets. We predicted that cardiac function can be best recovered in ZSF1 rats by desensitising the myofilament and reducing the affinity to intracellular calcium concentration and overall prolonging the sarcomere staying in the active force generating state.

Linking statistical shape models and simulated function in the healthy adult human heart.

Cardiac anatomy plays a crucial role in determining cardiac function. However, there is a poor understanding of how specific and localised anatomical changes affect different cardiac functional outputs. In this work, we test the hypothesis that in a statistical shape model (SSM), the modes that are most relevant for describing anatomy are also most important for determining the output of cardiac electromechanics simulations. We made patient-specific four-chamber heart meshes (n = 20) from cardiac CT images in asymptomatic subjects and created a SSM from 19 cases. Nine modes captured 90% of the anatomical variation in the SSM. Functional simulation outputs correlated best with modes 2, 3 and 9 on average (R = 0.49 ± 0.17, 0.37 ± 0.23 and 0.34 ± 0.17 respectively). We performed a global sensitivity analysis to identify the different modes responsible for different simulated electrical and mechanical measures of cardiac function. Modes 2 and 9 were the most important for determining simulated left ventricular mechanics and pressure-derived phenotypes. Mode 2 explained 28.56 ± 16.48% and 25.5 ± 20.85, and mode 9 explained 12.1 ± 8.74% and 13.54 ± 16.91% of the variances of mechanics and pressure-derived phenotypes, respectively. Electrophysiological biomarkers were explained by the interaction of 3 ± 1 modes. In the healthy adult human heart, shape modes that explain large portions of anatomical variance do not explain equivalent levels of electromechanical functional variation. As a result, in cardiac models, representing patient anatomy using a limited number of modes of anatomical variation can cause a loss in accuracy of simulated electromechanical function.

Tracking the motion of intracardiac structures aids the development of future leadless pacing systems.

BACKGROUND: Leadless pacemakers preclude the need for permanent leads to pace endocardium. However, it is yet to be determined whether a leadless pacemaker of a similar design to those manufactured for the right ventricle (RV) fits within the left ventricle (LV), without interfering with intracardiac structures. METHODS: Cardiac computed tomography scans were obtained from 30 patients indicated for cardiac resynchronisation therapy upgrade. The mitral valve annulus, chordae tendineae, papillary muscles and LV endocardial wall were marked in the end-diastolic frame. Intracardiac structures motions were tracked through the cardiac cycle. Two pacemaker designs similar to commercially manufactured leadless systems (Abbott's Nanostim LCP and Medtronic's Micra TPS) as well as theoretical designs with calculated optimal dimensions were evaluated. Pacemakers were virtually placed across the LV endocardial surface and collisions between them and intracardiac structures were detected throughout the cycle. RESULTS: Probability maps of LV intracardiac structures collisions on a 16-segment AHA model indicated possible placement for the Nanostim LCP, Micra TPS, and theoretical designs. Thresholding these maps at a 20% chance of collision revealed only about 36% of the endocardial surface remained collision-free with the deployment of Micra TPS design. The same threshold left no collision-free surface in the case of the Nanostim LCP. To reach at least half of the LV endocardium, the volume of Micra TPS, which is the smaller design, needed to be decreased by 41%. CONCLUSION: Due to the presence of intracardiac structures, placement of leadless pacemakers with dimensions similar to commercially manufactured RV systems would be limited to apical regions.

Calibration of Cohorts of Virtual Patient Heart Models Using Bayesian History Matching.

Previous patient-specific model calibration techniques have treated each patient independently, making the methods expensive for large-scale clinical adoption. In this work, we show how we can reuse simulations to accelerate the patient-specific model calibration pipeline. To represent anatomy, we used a Statistical Shape Model and to represent function, we ran electrophysiological simulations. We study the use of 14 biomarkers to calibrate the model, training one Gaussian Process Emulator (GPE) per biomarker. To fit the models, we followed a Bayesian History Matching (BHM) strategy, wherein each iteration a region of the parameter space is ruled out if the emulation with that set of parameter values produces is "implausible". We found that without running any extra simulations we can find 87.41% of the non-implausible parameter combinations. Moreover, we showed how reducing the uncertainty of the measurements from 10 to 5% can reduce the final parameter space by 6 orders of magnitude. This innovation allows for a model fitting technique, therefore reducing the computational load of future biomedical studies.

Global Sensitivity Analysis of Four Chamber Heart Hemodynamics Using Surrogate Models.

Computational Fluid Dynamics (CFD) is used to assist in designing artificial valves and planning procedures, focusing on local flow features. However, assessing the impact on overall cardiovascular function or predicting longer-term outcomes may requires more comprehensive whole heart CFD models. Fitting such models to patient data requires numerous computationally expensive simulations, and depends on specific clinical measurements to constrain model parameters, hampering clinical adoption. Surrogate models can help to accelerate the fitting process while accounting for the added uncertainty. We create a validated patient-specific four-chamber heart CFD model based on the Navier-Stokes-Brinkman (NSB) equations and test Gaussian Process Emulators (GPEs) as a surrogate model for performing a variance-based global sensitivity analysis (GSA). GSA identified preload as the dominant driver of flow in both the right and left side of the heart, respectively. Left-right differences were seen in terms of vascular outflow resistances, with pulmonary artery resistance having a much larger impact on flow than aortic resistance. Our results suggest that GPEs can be used to identify parameters in personalized whole heart CFD models, and highlight the importance of accurate preload measurements.

Regional diastolic dysfunction in post-infarction heart failure: role of local mechanical load and SERCA expression.

AIMS: Regional heterogeneities in contraction contribute to heart failure with reduced ejection fraction (HFrEF). We aimed to determine whether regional changes in myocardial relaxation similarly contribute to diastolic dysfunction in post-infarction HFrEF, and to elucidate the underlying mechanisms. METHODS AND RESULTS: Using the magnetic resonance imaging phase-contrast technique, we examined local diastolic function in a rat model of post-infarction HFrEF. In comparison with sham-operated animals, post-infarction HFrEF rats exhibited reduced diastolic strain rate adjacent to the scar, but not in remote regions of the myocardium. Removal of Ca2+ within cardiomyocytes governs relaxation, and we indeed found that Ca2+ transients declined more slowly in cells isolated from the adjacent region. Resting Ca2+ levels in adjacent zone myocytes were also markedly elevated at high pacing rates. Impaired Ca2+ removal was attributed to a reduced rate of Ca2+ sequestration into the sarcoplasmic reticulum (SR), due to decreased local expression of the SR Ca2+ ATPase (SERCA). Wall stress was elevated in the adjacent region. Using ex vivo experiments with loaded papillary muscles, we demonstrated that high mechanical stress is directly linked to SERCA down-regulation and slowing of relaxation. Finally, we confirmed that regional diastolic dysfunction is also present in human HFrEF patients. Using echocardiographic speckle-tracking of patients enrolled in the LEAF trial, we found that in comparison with controls, post-infarction HFrEF subjects exhibited reduced diastolic train rate adjacent to the scar, but not in remote regions of the myocardium. CONCLUSION: Our data indicate that relaxation varies across the heart in post-infarction HFrEF. Regional diastolic dysfunction in this condition is linked to elevated wall stress adjacent to the infarction, resulting in down-regulation of SERCA, disrupted diastolic Ca2+ handling, and local slowing of relaxation.