A biomechanics-based parametrized cardiac end-diastolic pressure–volume relationship for accurate patient-specific calibration and estimation

A simple power law has been proposed in the pioneering work of Klotz et al. (Am J Physiol Heart Circ Physiol 291(1):H403–H412, 2006) to approximate the end-diastolic pressure–volume relationship of the left cardiac ventricle, with limited inter-individual variability provided the volume is adequately normalized. Nevertheless, we use here a biomechanical model to investigate the sources of the remaining data dispersion observed in the normalized space, and we show that variations of the parameters of the biomechanical model realistically account for a substantial part of this dispersion. We therefore propose an alternative law based on the biomechanical model that embeds some intrinsic physical parameters, which directly enables personalization capabilities, and paves the way for related estimation approaches.

The passive behavior of the heart translates into the so-called end-diastolic pressure-volume relationship (EDPVR), which gives the volume of the left ventricle (LV) reached under a given filling pressure. The patientspecific determination of this curve is an important objective, both for clinical-e.g. for diagnosing pathological values of the myocardium compliance-and for biomechanical modeling purposes-for calibrating the constitutive parameters of the passive behavior. In the paper 1 a most interesting property of the EDPVR was established, i.e. that provided the volume is adequately normalized this curve displays limited inter-individual-or even inter-species-variability, and can be approximated in the single form with V ref and V 30 being the volumes associated with zero pressure and 30 mmHg, respectively. Throughout this paper we will refer to this as the "Klotz et al. relationship". This relationship is extremely valuable to estimate the complete EDPVR with very few measured pressure-volume data points, and also to calibrate the coefficients of a biomechanical model based on the full EDPVR curve estimated with this method, rather than with the scarce original data points, see [2][3][4][5] for a few examples. Conversely, the literature is very scarce concerning how a biomechanical model could be used to shed more light on the Klotz et al. relationship, and to provide a critical assessment tool for this relationship, whether it be for its shape itself, or for the dispersion of the data points around the curve. The paper 6 goes in this direction, by first proposing a modified curve with a non-zero slope at the origin, and concluding that the fitting of a biomechanical model with this modified curve leads to yet another different shape.
In the present work, we show by using the biomechanical model of 7 that the remaining dispersion of the data points around the Klotz et al. relationship in the normalized volume/pressure space may be accounted for by parametric variations in the model, in particular with a stiffness parameter. Based on this finding, we propose a new parametrized pressure-volume relationship that can be made patient-specific by adjusting the embedded parameters. Therefore, our proposed relationship has the built-in ability to be more accurate-for a given individual-than the "averaged" relationship of Klotz et al. www.nature.com/scientificreports/ In "Methods" section (Methods) we summarize the passive cardiac biomechanical model considered in this paper, and perform a calibration of this model based on actual data collected in a clinical setting of general anesthesia. In "Results" section (Results) we report on parameter sensitivity with respect to variations of a stiffness parameter and of the thickness of the ventricular wall. Moreover, we compare personalization capabilities of the biomechanical model and of the Klotz et al. relationship, when considering typical patient-specific configurations associated with different values of the stiffness parameter. We then discuss these results in "Discussion" and "Limitations" sections, and finally provide conclusions in "Conclusions" section.

Methods
Biomechanical model. We consider the simplified cardiac model proposed in 7 , itself based on the 3D model proposed and analysed in 8 and derived by approximating the left ventricular geometry as a hollow sphere, and assuming that fiber directions are isotropically distributed across the thickness in order to have complete spherical symmetry. The constitutive behavior, by contrast, is preserved in its full complexity as in 8 , and since we are concerned here with the EDPVR we focus on the static-or quasi-static-passive part, modeled based on a hyperelastic potential inspired from 9 and in the following exponential form with J 1 and J 4 being the so-called first and fourth reduced invariants of the right Cauchy-Green deformation tensor C , i.e., where e fib denotes the relative extension in the fiber direction, and [C 0 , C 1 , C 2 , C 3 ] is a set of constitutive parameters. We also assume that the cardiac tissue is perfectly incompressible, i.e., det C = 1 . Under the simplified geometric and kinematical assumptions considered here, the fiber extension is uniform throughout the whole cardiac tissue at any time, and related to the changes in spherical radius by where R and R ref respectively denote the radius in the current and reference configurations. Note that, although the spherical approximation may seem quite drastic, this simplified model has demonstrated an excellent adequacy-when confronted with actual clinical or experimental data-to represent the pressure-volume behavior of a single ventricle in health and disease 7,10,11 .
Denoting now by d and d ref the thickness of the wall in the current and reference configurations, respectively, tissue incompressibility dictates that 4π(R ref The internal volume of the ventricle in the current and reference configurations is given by and therefore where ǫ denotes the dimensionless thickness Finally, in the framework of our biomechanical model, it can be established that the EDPVR can be written in the form where P model is a function that interrogates the derivatives of the hyperelastic potential (2) with respect to J 1 and J 4 , see 7 (Eq. (19)), restricted to the static passive case. Note that in order to interpret (7) as a function of the volume V we can invert the relationship (6) in the form of an implicit function e fib V /V ref .

Model calibration.
Our purpose here is to calibrate the above biomechanical cardiac model for a typical human heart. To that aim we use some ultrasound data obtained for a collection of patients anesthetized for a neuroradiological procedure, and select one specific individual with a heart of "average" dimensions within the population considered in the study 12 in which 45 patients were included, see Table 1 for the specific characteristics of the patient that we selected. We emphasize that this selection of an individual dataset is not performed with a view to a patient-specific study, but to provide a reference EDPVR that can be compared with the data of 1 and from which whole families of EDPVRs can be generated by varying geometric and mechanical parameters.
Data description and determination of (P, V) point. The ultrasound data were obtained using trans-thoracic echocardiography (TTE), performed after the general anesthesia induction. This study was approved by the www.nature.com/scientificreports/ appropriate Institutional Review Board -ethical committee of the Société de Réanimation de Langue Française (CE-SRLF 14-34)-which waived the need for written informed consent. Consequently, oral informed consent was obtained after providing a protocol information letter, and the subject had the possibility to withdraw from the study at any time.
The apical 4-chambers view was recorded to estimate the left ventricular (LV) volume at end-diastole V ED , based on the recommendations given in 13 . In short, the LV cavity is discretized into N disks characterized by their heights (h i ) and diameters (D i ) . Therefore, V ED can be approximated as For our selected patient, using N = 550 disks with h i = 0.1546 mm, this rule led to V ED = 133 ml.
The associated end-diastolic pressure P ED is then estimated using a semi-quantitative method described in detail in 14 . In short, by analyzing the characteristics of the transmitral flow, a classification of the LV filling pressure as normal or high is determined. In our case, this procedure concluded to a normal LV filling pressure, i.e. P ED = 7 mmHg.
Calibration of geometric parameters. The main objective of 1 is to provide a method for estimating a patientspecific EDPVR based on very few measured pressure-volume points, including with one single such point. This has proven to be extremely valuable to calibrate the passive behavior of biomechanical cardiac models in a number of studies in which the model parameters have been calibrated by using the EDPVR obtained by the method of 1 applied with the available data, instead of directly calibrating the biomechanical model against the data. This strategy clearly makes the calibration-or estimation-step much better posed than when using the scarce data points. In our case, we only use the method of 1 to estimate the patient-specific volume parameters V ref and V 30 instead of the full EDPVR. We thus use the estimate of 1 (Eq. (8)), i.e.
which we then substitute in (1) to obtain In addition, we estimate the LV wall thickness in the reference configuration-corresponding to V ref . To that purpose, we first determine the external volume V ext by using the technique already mentioned above to obtain V ED in (8) albeit here with the diameters associated with the epicardium. Within the simplified spherical geometry assumption of the biomechanical model considered here-recall "Biomechanical model" section-the thickness d ref is directly inferred by subtracting the radii associated with V ext and V ref , i.e.
where we have used the fact that-due to incompressibility-the total volume of the ventricular wall is the same between the reference and end-diastolic configurations, hence, given by V ext − V ED . Finally, the radius of the midsurface in the spherical geometry is given by www.nature.com/scientificreports/ Calibration of mechanical parameters. Our objective here is to calibrate the mechanical parameters (C i ) for the selected patient. As already mentioned, one possible strategy would be to use the complete EDPVR provided by the methods of 1 based on the single point (P ED , V ED ) , in order to formulate an identification problem as wellposed as possible. However, this would introduce a bias towards the type of curve shape associated with the Klotz et al. relationship, and therefore we instead choose to use some of the data considered in 1 , namely, the ex vivo human heart data corresponding to 80 normal and diseased hearts of various etiologies. We will thus propose an identification criterion adapted to the data points provided in 1 (Fig. 3), i.e. with normalized volumes. Noting that the normalized volume V n can be written as we have a very simple linear relationship between V n and V /V ref that only depends on the parameter V 30 /V ref , which in our case is calibrated as described above. Therefore, for any choice of the set of mechanical parameters (C i ) we can plot the modeled EDPVR associated with (7) in the (V n , P) space, and compare it with the data points (V j n , P j data ) N data j=1 considered in 1 . In this paper we then calibrate the mechanical parameters by minimizing the following root mean square error (RMSE) where V j /V ref is directly inferred from V j n by inverting the linear relationship (11). In practice, we will carry out a minimization process in two steps. First, using the Klotz et al. relationship, we generate a synthetic set of P data to be used in the minimization problem (12) with a steepest gradient descent algorithm. This will provide preliminary mechanical parameters that can subsequently be taken as an initial guess for minimizing (12) with the actual data points of 1 , and using the quasi-Newton algorithm provided in the Matlab fminunc function.

Results
Biomechanical model calibration. The geometric parameters calibrated as explained in "Model calibration" section are listed in Table 1.
The mechanical parameters resulting from this minimization are found to be The resulting curve in the (V n , P) space is shown in Fig. 1 and compared with the relationship and data of Klotz et al. 1 .
Sensitivity with respect to stiffness parameter. We call "stiffness parameters" the material parameters that have a direct multiplicative effect in the mechanical stress-strain law, i.e. in our case C 0 and C 2 . Therefore we will investigate the impact of multiplicative stiffness variations by considering  www.nature.com/scientificreports/ with [C 0 ,C 1 ,C 2 ,C 3 ] being the optimized parameters obtained in "Biomechanical model calibration" section, and α a single coefficient driving stiffness variations. In this case the multiplicative effect of C 0 and C 2 in the stressstrain law directly translates for the pressure law (7)  [C 0 , C 1 , C 2 , C 3 ] = [αC 0 ,C 1 , αC 2 ,C 3 ],  www.nature.com/scientificreports/ Sensitivity with respect to thickness parameter. The biomechanical model embeds the thickness parameter ǫ as an explicit parameter, and it is interesting to study the effect of thickness variations on the EDPVR. By contrast, there is no such parameter in the Klotz et al. relationship, and therefore we will focus on the biomechanical model in this section. As the biomechanical model corresponds to a complex membrane, it is expected-from classical structural mechanics-that the associated stiffness is roughly proportional to the membrane thickness. To investigate this more precisely, we plot in Fig. 4 the scaled quantity P model (C i , ǫ, e fib )/ǫ as a function of e fib for various values of ǫ . Since we observe that the resulting curves are remarkably insensitive to variations in ǫ , we can directly conclude that-in an excellent approximation-the function P model (C i , ǫ, e fib ) is linear in ǫ , indeed.
Comparison of personalization capabilities. In this section we will illustrate and compare the personalization capabilities of the biomechanical model and of the Klotz et al. relationship, in the case of patientspecific stiffness variations.
We start from the calibration performed in "Biomechanical model calibration" section and consider the cases of two different hearts associated with stiffness multiplicative coefficients α = 0.5 and α = 2 , with all other parameters left unchanged. The biomechanical model is assumed to provide the reference for the resulting modified EDPVR curves, since stiffness variations are adequately represented by construction in such a model via the definition of the stress-strain law.
We compare in Fig. 5 the EDPVR curves obtained with the biomechanical model and the Klotz et al. relationship, the latter being adapted by applying V 30 given by the biomechanical model.

Discussion
We now discuss some specific aspects pertaining to the above results.
Biomechanical model calibration The RMSE value obtained with the calibrated biomechanical model is smaller by 6% than that associated with the Klotz et al. relationship (1.52 mmHg versus 1.61 mmHg). From a qualitative point of view, as observed in Fig. 1, the major difference lies in the initial part of the EDPVR where the biomechanical model provides a better fit on average, which was expected due to the vanishing slope of the Klotz et al. relationship at the origin-an unphysical feature because this corresponds to the stiffness near the reference configuration. It should also be pointed out that the optimization process used in this calibration has revealed a rather limited sensitivity of the mechanical parameters-in particular for C 0 and C 2 -albeit this is of little consequence in this study, because we only use the shape of the resulting EDPVR in the normalized (V n , P) space-which directly corresponds to the optimization criterion-and not the specific values of the mechanical parameters.
Sensitivity with respect to stiffness and thickness parameters As directly observed in Fig. 1, there exists some significant dispersion in the EDPVR data considered in the normalized (V n , P) space. The Klotz et al. relationship-based on a nonlinear regression of the data in this space-treats this dispersion as measurement noise. By contrast, the biomechanical model incorporates some physical parameters, and we see in Fig. 3 that stiffness variations induce deformations of the EDPVR that can instead provide some interpretation for the data dispersion. The same holds with thickness variations since, as seen in "Sensitivity with respect to thickness parameter" section such variations have a similar-linear-effect on the structural stiffness. This is both valuable from a qualitative-physical interpretation-and quantitative point of view, since by adjusting these parameters on a personalised basis we can expect to obtain more accurate results.
Comparison of personalization capabilities Upon recalibrating the EDPVR according to stiffness variations associated with multiplicative coefficients α = 0.5 and α = 2 , the RMSE with the synthetic data is larger by 63% and 113%, respectively, when using the Klotz et al. relationship instead of the biomechanical model, see Fig. 5 (left). All the RMSE values may seem rather small in absolute value, but they should be compared to the noise www.nature.com/scientificreports/ present in the synthetic data, i.e. 1 mmHg. Moreover, the comparison between the synthetic data range and the actual data points of 1 in the right-hand side of Fig. 5 further illustrates the fact that the dispersion in the actual data may be realistically accounted for-up to some reasonable amount of measurement noise-by parametric variations associated with the stiffness parameter, or equivalently with the thickness of the ventricular wall. This also provides an interpretation of the benefits provided by the biomechanical model that allows to fit the behavior with the relevant subregion of the Klotz et al. data, as shown by the respective positions of the two different derived EDPVRs in the same graphs. Of course, this personalization study is only preliminary because it is performed based on synthetic data-i.e. EDPVR data generated by the biomechanical model-and a definite validation of this would require actual data of a very invasive type, as also discussed in "Limitations" section below.
Effective computation of biomechanical model law In some applications it may be valuable-or even necessary-to be able to compute the function (7) in a very effective manner. This is the case in particular when inverse problems are considered e.g. when seeking to estimate patient-specific parameters based on actual data. To that purpose, we can benefit from the sensitivity explorations performed in "Sensitivity with respect to stiffness parameter" and "Sensitivity with respect to thickness parameter" sections and from the linearity of this function with respect to both the stiffness-when considering a single multiplicative coefficient α-and thickness parameters. As a consequence, indeed, we can consider the function with ǭ being the (fixed) thickness coefficient calibrated for our-supposedly "average"-case above, and perform a fit to an adequate-e.g. polynomial-function to any desired accuracy. As an example, a fit with a polynomial of degree 7 gives an RMSE of 0.03 mmHg for e fib between 0 and 0.25 (corresponding to ventricular pressures between 0 and 30 mmHg), which of course is more than sufficient in actual applications. We provide the details of this polynomial fit as supplementary material for reproducibility purposes. Denoting then by P fit (e fib ) this approximation of the function (16), for other values of the stiffness and thickness parameters we can approximate P model αC i , ǫ, e fib by α(ǫ/ǭ)P fit (e fib ) with excellent accuracy. Moreover, we emphasize that some very valuable (16) e fib → P model C i ,ǭ, e fib , Figure 5. Comparison of personalized models for stiffness coefficients α = 0.5 and α = 2 . The synthetic data are generated from the biomechanical model with stiffness coefficients α = [0.5, 2] , with random noise added according to a centered reduced normal law with standard deviation 1 mmHg. The RMSE calculated as in (12) between synthetic data and the biomechanical model or the Klotz

Limitations
This study has two main limitations. The first limitation pertains to the scarcity of the data considered. Even though we conjecture that using different datasets for other normal hearts of average dimensions would not significantly impact the results, since-again-we essentially use the fitted shape of the biomechanical model in the normalized volume vs. pressure space, this warrants some detailed verification in further studies. An extensive verification would in fact require some invasive data with measured multiple-points EDPVRs for several normal hearts, as well as for several pathological hearts for which the passive stiffness is suspected to be altered. The second limitation concerns the model used here, in particular due to the spherical approximation made for the geometry in order to avoid the complexity of a 3D model. Other types of reduced models-such as those proposed in 15,16 -could be considered to verify that our findings are not significantly impacted. A 3D computational model could also be employed and calibrated as in 6 against the data of Klotz et al., albeit it is well-known that 3D models are very sensitive to the specific design of their artificial boundary conditions where the geometric model is truncated, and therefore would not necessarily produce more trustworthy EDPVRs in practice. Other effects could also be included in the model to assess their impact, such as variations in tissue volume. Such variations have been evidenced, indeed-see e.g. 17 -and have often been attributed to the effects of blood perfusion 18 . A reduced model where such effects are embedded-see e.g. 19 -could thus be used for comparative assessment purposes, and the same holds for 3D poromechanical models 20 , although detailed validations of such models are still needed in the cardiac setting.

Conclusions
We have proposed a law based on a biomechanical model to parametrize a cardiac EDPVR, as an alternative to the phenomenological law proposed in the pioneering work of Klotz et al. 1 . Our proposed law is expressed via the natural variable associated with the fiber extension, and features two major advantages: (1) it incorporates some intrinsic, easily interpretable, physical parameters-namely, the volume in the reference configuration, the wall thickness, and a stiffness parameter-that directly enable personalization capabilities; (2) variations of these parameters may realistically account for the dispersion observed in actual data. An additional valuable finding is that thickness variations have a linear (multiplicative) effect in this law, in an excellent approximation. Finally, polynomial fits can be used as surrogate models of arbitrary accuracy to avoid having to deal with the intricacies of the biomechanical model.
A very natural perspective concerns the use of our proposed law for inverse problem purposes, in order to adapt the EDPVR to measured data points by estimating patient-specific values of the embedded parameters. Using physically interpretable parameters will enrich the output of the estimation process, indeed, and also make it more amenable to quantitative assessments (since such parameters can sometimes be independently measured by specific means).

Data availability
The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request.