Haemodynamic effects of hypertension and type 2 diabetes: Insights from a 4D flow MRI‐based personalized cardiovascular mathematical model

Type 2 diabetes (T2D) and hypertension increase the risk of cardiovascular diseases mediated by whole‐body changes to metabolism, cardiovascular structure and haemodynamics. The haemodynamic changes related to hypertension and T2D are complex and subject‐specific, however, and not fully understood. We aimed to investigate the haemodynamic mechanisms in T2D and hypertension by comparing the haemodynamics between healthy controls and subjects with T2D, hypertension, or both. For all subjects, we combined 4D flow magnetic resonance imaging data, brachial blood pressure and a cardiovascular mathematical model to create a comprehensive subject‐specific analysis of central haemodynamics. When comparing the subject‐specific haemodynamic parameters between the four groups, the predominant haemodynamic difference is impaired left ventricular relaxation in subjects with both T2D and hypertension compared to subjects with only T2D, only hypertension and controls. The impaired relaxation indicates that, in this cohort, the long‐term changes in haemodynamic load of co‐existing T2D and hypertension cause diastolic dysfunction demonstrable at rest, whereas either disease on its own does not. However, through subject‐specific predictions of impaired relaxation, we show that altered relaxation alone is not enough to explain the subject‐specific and group‐related differences; instead, a combination of parameters is affected in T2D and hypertension. These results confirm previous studies that reported more adverse effects from the combination of T2D and hypertension compared to either disease on its own. Furthermore, this shows the potential of personalized cardiovascular models in providing haemodynamic mechanistic insights and subject‐specific predictions that could aid in the understanding and treatment planning of patients with T2D and hypertension.


Introduction
Hypertension and type 2 diabetes (T2D) are powerful risk factors for cardiovascular diseases such as coronary artery disease, atrial fibrillation, stroke and renal failure (Forouzanfar et al., 2017;Sarwar et al., 2010). The diseases frequently co-exist: hypertension is twice as common in subjects with T2D than in normoglycaemic subjects (af Geijerstam et al., 2022;Naseri et al., 2022) and co-existing hypertension and T2D can increase the risk of cardiovascular disease by almost three-fold (Verdecchia et al., 2004). The interaction between T2D, hypertension and cardiovascular disease can be understood through subject-specific effects on the regulation of haemodynamics such as blood flow and blood pressure. These mechanisms are difficult to assess non-invasively and thus not fully understood.
In both T2D and hypertension, these haemodynamic changes are complex and impact the entire cardiovascular system. In hypertension, one can see changes such as vascular stiffening, ventricular hypertrophy, increased peripheral resistance and diastolic dysfunction (Nadar & Lip, 2020;Petrie et al., 2018). In T2D, similar changes can be seen, with impaired left ventricular function and decreased compliance leading to prolonged relaxation, increased filling pressure and eventually increased atrial pressure as the disease progresses (Tan et al., 2020). Especially in T2D, metabolic and renal dysfunction contribute to the haemodynamic changes, such as insulin resistance leading to increased insulin-dependent sodium retention and risk of renal damage (Brands & Manhiani, 2012;Petrie et al., 2018), both causing increased plasma volume resulting in increased cardiac output and increased blood pressure (Kawasoe et al., 2017). Over time, these haemodynamic changes can lead to adverse remodelling of the heart and cardiovascular system, including left ventricular hypertrophy and atherosclerosis. These diseases in turn can cause more severe cardiovascular diseases such as heart failure, ischemic heart disease or stroke (Koren et al., 1991;Nadar & Lip, 2020;Tan et al., 2020).
Although all of these haemodynamic changes are seen on a group level, there are large patient-specific variations in the haemodynamic mechanisms related to both hypertension and T2D. In hypertension, these patient-specific factors often result in a trial-and error approach to anti-hypertensive treatment (Williams et al., 2018). Similarly, the understanding of subject-specific haemodynamic variations in T2D is not yet fully incorporated into clinical practice (Chung et al., 2020). Most studies have focused on long-term outcomes such as the risk of cardiovascular diseases, and not on understanding the complex and patient-specific haemodynamic 0 Kajsa Tunedal is a PhD student at Linköping University, in the Department of Biomedical Engineering and the Centre for Medical Imaging and Visualization (CMIV). She received her MSc in Engineering Biology in 2021 after completing an additional 1-year research internship in systems biology. Currently, she splits her time between one research group within systems biology, one group within cardiovascular imaging and modelling, and collaborations with medical staff. In her PhD project, she focuses on understanding the mechanisms of short-and long-term blood pressure regulation in health and disease by creating patient-specific mechanistic cardiovascular models based on non-invasive measurements. mechanisms behind those risks (Chung et al., 2020;Forouzanfar et al., 2017;Sarwar et al., 2010;Williams et al., 2018). New tools and a deeper understanding of the associated haemodynamic mechanisms are needed to personalize the diagnosis and treatment of these two common diseases. Subject-specific haemodynamics in hypertension and T2D can be measured using various techniques. To measure pressure changes, cuff-based methods are the main non-invasive methods and can measure pressure in the limbs, whereas invasive methods such as catheter-based pressure sensors are needed to measure pressure in the heart and larger vessels (McEniery et al., 2014). Blood flow can be measured through non-invasive imaging techniques such as cardiac ultrasound and magnetic resonance imaging (MRI). Ultrasound is widely used in the clinic to assess velocity and geometry but has limited precision as a result of the beam angle and assumptions about haemodynamics (Kadappu & Thomas, 2015). By contrast, MRI provides more standardized measurements of cardiac anatomy and function with large spatial coverage. Four-dimensional phase-contrast MRI (4D flow MRI) quantifies the time-resolved 3D velocity field in the heart and thorax . 4D flow MRI velocity data can be used to calculate haemodynamic information such as blood flow over time, wall shear stress or pulse wave velocity in the aorta . However, 4D flow MRI cannot directly measure absolute pressures or ventricular and atrial contractility, and it cannot be used on its own to comprehensively analyse all the interrelated components of the complex haemodynamics at play.
To be able to comprehensively analyse haemodynamics in T2D and hypertension, one can use cardiovascular mathematical models. Cardiovascular models can be trained on all subject-specific data and then used to predict how the haemodynamics add up and also to provide new model-derived haemodynamic measurements (Casas et al., 2017(Casas et al., , 2018. The combination of cardiovascular models and subject-specific data in so-called personalized cardiovascular models is increasingly used to understand haemodynamics. Personalized models range from larger 3D electrophysiological and fluid dynamics models (Lantz et al., 2016;Sermesant et al., 2012), which are computationally expensive, to smaller J Physiol 601.17 and less computationally heavy lumped parameter (0D) or one-dimensional (1D) models, which provide fast patient-specific simulations of both central and global haemodynamics suitable in a clinical setting (Shi et al., 2011). Personalized 0D-1D models can estimate subject-specific aortic pressure in healthy subjects (Gallo et al., 2021;Kondiboyina et al., 2022;Mariscal-Harana et al., 2021) and describe baroreflex regulation (Randall et al., 2019), but do not focus on assessing other haemodynamic mechanisms important in hypertension and diabetes. Several non-patient-specific models have been used to investigate haemodynamic aspects of hypertension such as sodium balance (Guyton et al., 1972) and the contribution of both cardiac remodelling and increased peripheral resistance to increased blood pressure (Segers et al., 2000). Additionally, Ouyang et al. (2020) used a simple lumped parameter model to classify T2D patients and controls. We have previously developed a personalized lumped-parameter model that uses 4D flow data: the cardiovascular avatar (Casas et al., 2017(Casas et al., , 2018. However, the potential of personalized cardiovascular mathematical models has not yet been used to investigate the patient-specific haemodynamic mechanisms in hypertension and T2D.
In the present study, we aim to elucidate the haemodynamic mechanisms in hypertension and T2D by combining 4D flow MRI and a personalized cardiovascular mathematical model. To do this, we first aim to compare model-derived haemodynamic parameters between controls, patients with T2D, patients with hypertension, and patients with both T2D and hypertension. To further explore the group-level and patient-specific mechanisms behind the differences, we perform clustering, correlation with markers of chronic disease and concurrent blood pressure and heart rate, and use the cardiovascular model to create patient-specific predictions of haemodynamics.

Cohort and study design
In a sub-study of the SCAPIS study (Bergström et al., 2015), 46 patients with T2D and 46 age-sex-and smoking-matched controls were recruited for an MRI examination at rest including brachial pressure measurement (Edin et al., 2022). Additionally, cohort characteristics from the core SCAPIS study were retrieved (Bergström et al., 2015;Edin et al., 2022). Informed consent was obtained from all participants and the study was approved by the Regional Ethical Review Boards in Linköping (Dnr 2016/229-31 and2018/478-31) and Umeå (Dnr 2010-228-31M) and adheres to the standards in the Declaration of Helsinki. T2D patients were identified as either self-reported T2D in health forms, fP-glucose ≥ 7.0 mmol L -1 on two occasions or HbA1c ≥ 48 mmol L -1 . Controls were included according to no self-reported T2D, fP-glucose < 6.0 mmol L -1 and HbA1c < 42 mmol L -1 . Exclusion criteria for all subjects were significantly irregular rhythm and contraindications to MRI. Additionally, all subjects without brachial pressure measurement (n = 8, 6 T2D and two controls), two subjects with self-reported valvular heart disease (n = 2, T2D) and two subjects with bad MRI data quality including excessive phase wraps (n = 1, T2D) or patient movement (n = 1, T2D) were retrospectively excluded, resulting in 44 controls and 36 subjects with T2D (i.e. 80 subjects in total). In retrospective grouped analysis, hypertensive subjects were defined as subjects with a mean brachial home blood pressure of systolic blood pressure (SBP) ≥ 135 or diastolic blood pressure (DBP) ≥ 85 mmHg using additional data from the SCAPIS study (Johansson et al., 2021). Home blood pressure was measured in the morning and evening for 7 days, three times per occasion, in a seated position with the semi-automatic oscillometric Omron m10-IT (Omron Healthcare, Kyoto, Japan), resulting in a mean value of 13 × 3 measurements, as described in detail in Johansson et al. (2021). The 80 subjects were divided into four groups: controls without T2D or hypertension, subjects with only T2D, subjects with only hypertension, and subjects with both T2D and hypertension.

MRI data acquisition
During the MRI session, 4D flow data as well as cardiac 3D short-axis cine and 2D 3-chamber cine balanced steady state-free precession (bSSFP) images were acquired. The 4D flow MRI acquisitions were retrospectively cardiac gated and respiratory navigator gated echo-planar imaging sequence with read-out factor 7 on a 1.5-T Philips Ingenia (Philips Healthcare, Best, The Netherlands). Other scan parameters included sagittal-oblique slab covering the whole heart and the thoracic aorta, velocity encoding 120 cm s -1 , flip angle 5°, echo time (TE) 5 ms, repetition time (TR) 9.1 ms, parallel imaging (sensitivity encoding) speed up factor 2 in the AP and the RL direction, acquired spatial resolution 2.9 × 3.1 × 2.9 mm 3 , and acquired temporal resolution 36 ms. The 4D flow data were reconstructed to 40 timeframes and a spatial resolution of 2.7 × 2.7 × 2.9 mm 3 . Typical scan time was ∼7-8 min, navigator excluded. ECG was used for all subjects except for 14 of the included subjects where peripheral pulse measurement (Ppu) was used instead.
The cine 3-chamber images were acquired with TE 1.6 ms, TR 3.2 ms, flip angle 60°, acquired spatial resolution 2.5 × 2.5 mm 2 , slice thickness 6 mm, acquired temporal resolution 49-69 ms, breath hold duration 16 s, reconstructed spatial resolution 1.0 × 1.0 mm, and 30 reconstructed heart phases. Finally, the cine 3D sequence was acquired in a single breath-hold with TR 2.7-2.8 ms, TE 1.4 ms, flip angle 50°, acquired spatial resolution 2.5 × 2.5 × 8.0 mm 3 , reconstructed spatial resolution 1.0 × 1.0 × 8.0 mm 3 , acquired temporal resolution 49-69 ms, and breath-hold duration 16 s. Directly before the MRI acquisition, the systolic and diastolic brachial blood pressure was measured in a supine position using an automatic blood pressure monitor.

MRI post-processing
Each 4D flow MRI dataset was corrected for phase wraps with Laplacian unwrapping (Loecher et al., 2016) and for background phase errors using a weighted second-order polynomial fit to static tissue (Ebbers et al., 2008). Atlas-based time-resolved segmentations of the aorta and left atrium were automatically created (Bustamante et al., 2015). Blood flow through the mitral and aortic valves was calculated using valve-tracking accounting for valve motion (Viola et al., 2020;Westenberg et al., 2005), as in Casas et al., 2017 andCasas et al., 2018. The flow in the ascending aorta was calculated in a plane automatically placed in the mid-ascending aorta. The pulse wave velocity (PWV) in the whole aorta was calculated by using the cross-correlation method on flow rate curves extracted from planes automatically positioned along the centerline of the aorta (Dyverfeldt et al., 2014). The ascending aortic volume was calculated by taking each plane area times the distance to the next plane and summing the resulting volume segments along the aorta. The end-diastolic volume (EDV) and end-systolic volume (ESV) in the left ventricle were derived from short-axis MRI from the cardiac 3D cine bSSFP (Edin et al., 2022), and the left ventricular stroke volume was calculated as EDV -ESV.

Subject-specific cardiovascular models and parameter estimation
For all 80 subjects, personalized models were created by fitting the model parameters to patient-specific data extracted from 4D flow MRI, cine images and blood pressure (Fig. 1). The cardiovascular mathematical model used was a previously published cardiovascular lumped parameter model (Casas et al., 2017) translated into differential algebraic equations and applied in the Matlab AMCI toolbox (Fröhlich et al., 2021). The model describes the left side of the central and systemic circulation with seven compartments modelled by combinations of resistances R, compliances C and inductances L, representing frictional losses, vessel wall elasticity and mass flow inertia, respectively (Fig. 1C). The first compartment, the pulmonary venous system, has a constant pulmonary capillary pressure source Ppu and two vessel compartments. The contraction and relaxation of the atrial and ventricular compartments are described as time-varying elastance functions Ela(t) and Elv(t), driving the time-dependent changes in pressure and volume throughout the whole system (see Supporting information, Figure S1). The mitral and aortic valves are modelled as pressure difference-driven diodes, where the resistance Rav is a function of the energy loss coefficient of the aortic valve (ELCo) to account for pressure recovery. Finally, the ascending aorta and peripheral arteries are described as vessel segments, where the resistance Rpr = Rtot -Rao and the compliance Cpc = Ctot -Caa. A more detailed description of the model is provided in the first section in the Supporting information and Casas et al. (2017Casas et al. ( , 2018, and all equations are listed in the model file (Tunedal, 2023). The estimated model parameters are described in Table 1 and the constant parameters are described in the Supporting information (Table S1). The model was simulated until steady state (an absolute difference of all model states < 1 between two consecutive heartbeats or 20 heartbeats) and the last simulated steady-state heartbeat was compared to data.
The data used to personalize the models included blood flow over one cardiac cycle at four locations (pulmonary veins, mitral valve, aortic valve and ascending aorta), systolic and diastolic brachial pressure before MRI, and data on five model parameters calculated based on brachial pressure and MRI data. The simulated blood flow curves for the aortic valve and ascending aorta were only fitted to the data during systole, and the simulated blood flow curves for the mitral valve were only fitted to the data during diastole. To ensure that the first time point in the flow curve corresponded to isovolumetric contraction, the isovolumetric time point was set to the average of three rounds of manually chosen points based on the volumetric flow at the mitral valve, aortic valve and ascending aorta. In the same way, the time point for the start of diastole was chosen. The brachial pressures were estimated in the model as a simple function of the simulated aortic pressure as in Marazzi et al. (2022), with brachial DBP = aortic DBP and brachial SBP = aortic SBP + SBPdiff, where SBPdiff is a subject-specific parameter ranging between 0.1 and 18.9, which is the variation found in the Anglo-Cardiff Collaborative Trial II study for subjects aged 50-70 years old (McEniery et al., 2008). During parameter estimation, SBPdiff was set to the difference between measured SBP and simulated aortic SBP, or to 0.1 if the difference was smaller than 0.1 and 18.9 if the difference was larger than 18.9. The difference between brachial and aortic diastolic pressure was assumed to be neglectable.
Out of the five data-based parameter values, the total compliance (Ctot), total resistance (Rtot), maximum elastance in the left ventricle (Emax LV ) and the energy loss For the data-based parameters ( * ), the range given is (minimum bound of all subjects -maximum bound of all subjects). The literature values of l given are the same as in Casas et al. (2018). The new parameter bounds herein were in general set to: parameters describing compliance [l/7, l * 7], resistance [l/5, l * 5], inertance [l/6, l * 6], and left atrium and ventricle [l/3, l * 3], except for the five parameters that were fitted directly to calculated data values, where the bounds were instead set to ±75% of the data value. coefficient at the aortic valve (ELCo) were all calculated as in Casas et al. (2017Casas et al. ( , 2018) (Supporting information, first section). Finally, the total compliance in the aorta, Caa, was calculated from the aortic PWV and ascending aortic volume V based on the Bramwell-Hill equation (Bramwell & Hill, 1922): Caa = V ρ * PWV 2 mL mmHg -1 . To personalize the model, 28 of the model parameters were estimated (Table 1), whereas 11 parameters were kept constant (see Supporting information, Table S1). This is based on the approach employed in Casas et al. (2018), with the addition that we include an estimation of parameters describing the pulmonary veins, as a result of the added flow measurement in the pulmonary veins. The parameter estimation was achieved by minimizing the squared difference between simulation and data, the cost V: where y n t is the measured value for variable n at timepoint t, y n is all time points of the measured variable n, T n is the total number of time points for the variable n, y n t is the corresponding simulated value, θ is the set Figure 1. Creation of subject-specific cardiovascular mathematical models to derive non-invasive haemodynamic parameters A, subject-specific measurements from MRI and supine brachial blood pressure were collected and used to calculate five data-based model parameters. Four flow curves from the pulmonary veins (PV), mitral valve (MV), aortic valve (AV) and ascending aorta (AA) were also derived. B, the subject-specific measurements were used to estimate subject-specific parameters by minimizing the difference between simulations and data, shown for one subject here. The dotted lines in B correspond to the parameter bounds used in the optimization. C, the structure of the cardiovascular model, where each colour-coded part in the cardiovascular model (right) corresponds to an anatomical part in the heart and cardiovascular system (left). D, the resulting subject-specific models are used to derive more haemodynamic variables, which are then compared between the different groups. E, the subject-specific models are used to make patient-specific predictions of haemodynamics. of estimated parameters and w n is any weight for the variable n. The costs for flow at the aortic valve and aorta were weighted with 0.5 to even out the weight of flow information in the systole and diastole. Additionally, the costs for blood pressure and data-based parameters were weighted with a factor of 3, and penalization Z(k) for unphysiological oscillations in the blood flow in the aortic valve and ascending aorta was added with a weight of 30 per oscillatory flow curve k. In the parameter estimation, the parameter bounds (Table 1) were increased compared to the bounds used in Casas et al., 2017 to include variations in the larger and non-healthy population at the same time as remaining within physiological values by basing the bounds on the previously described literature values (Casas et al., 2017(Casas et al., , 2018. Parameter estimation was performed for each of the 80 subjects using the MEIGO toolbox (Egea et al., 2014), starting with a random guess within the parameter bounds and then running 40 × 5 consecutive enhanced scatter search optimizations. After the random start guesses, the optimization was repeated 10 times using the previous best parameter set as the start guess until the cost did not improve or until 20 consecutive optimizations. Finally, a last optimization with an increased maximal time of optimization was run to ensure that the solver converged. To estimate the parameter uncertainty, Markov chain Monte Carlo (MCMC) sampling with 10 5 iterations was performed for each subject using the PESTO toolbox (Stapor et al., 2018). Based on the best cost from the MCMC sampling, all parameter sets from both MCMC sampling and the found optima's from eSS optimizations with costs ≤ best MCMC cost + 10% of the best MCMC cost were included in the uncertainty analysis. In subject-specific simulations including uncertainty, all parameter sets or a maximum of 2000 parameter sets were simulated. In group comparisons, the median value of all included parameter values for each subject was used. All optimization problems were run in parallel at the Swedish National Supercomputer Centre (NSC) (Linköping University, Linköping, Sweden).
To predict the effect of reduced left ventricular relaxation, the dimensionless relaxation rate constant of the LA (m2LV) in the subject-specific optimal parameter set was lowered to 20, which corresponds to one of the lowest values across all subjects (Fig. 4M), and a lower value than the values of ∼30-50 found in previous studies of healthy volunteers (Casas et al., 2017(Casas et al., , 2018. In all other simulated parameter sets, m2LV was lowered with the same fraction between the optimized and predicted values to calculate the uncertainty of the prediction. All other parameters were kept the same. Finally, to compare the model-based and data-based stroke volumes, the stroke volume was calculated in both the data and simulation at the mitral valve, aortic valve and ascending aorta by integrating the blood flow volume at each location with respect to time.

Statistical analysis
To test whether the model-derived parameters and cohort characteristics differed between the four subject groups, pairwise comparisons were used. First, a Shapiro-Wilk test was applied for each variable to test if the variable was normally distributed for all groups. For normally distributed variables, pairwise t tests were then used, and, for non-parametric variables, a Wilcoxon rank sum test was used. For pairwise comparisons of categorical variables, a chi-squared test was used if all expected values ≥ 1 and at least 50% of the expected values ≥ 5. Otherwise, Fisher's exact test was used. To correct for the six multiple comparisons between groups, Benjamini-Hochberg correction was performed with a false discovery rate of 5%. When choosing which of the 28 model parameters to study further, all parameters where any of the six pairwise comparisons had P < 0.05 were included, but the parameters with comparisons that remained statistically significant after Benjamini-Hochberg correction were marked with an asterisk ( * ). Data are expressed as the mean ± SD if normally distributed, as median (interquartile range) if any group is not normally distributed, and as a percentage if categorical. A detailed summary of all statistical tests and their respective results is presented in the Statistical Summary Document.
To investigate the haemodynamic mechanisms behind the group differences, eight selected parameters were checked for correlation with indicators of chronic disease (i.e. home blood pressure, HbA1c and diabetes duration) and with markers of stress concurrent with data acquisition (i.e. heart rate and blood pressure during MRI). Pearson correlation was used if both variables were normally distributed, otherwise, Spearman's correlation was used. First, three of the eight parameters were selected as the parameters that were significantly different between any of the four groups. Second, because previous studies have shown changes in aortic stiffness and resistance in hypertension and diabetes (Nuamchit et al., 2020;Tedesco et al., 2004), the correlations with aortic compliance (Caa) and aortic resistance (Rao) were checked. Finally, in a previous modelling study of the effects of dobutamine, which include increases in heart rate and blood pressure, we showed changes in four haemodynamic parameters describing ventricular function: a decrease in the LV diastolic time constant and an increase in the maximum elastance of the left ventricle, ventricular relaxation rate and ventricular contraction rate with dobutamine (Casas et al., 2018). Thus, these four parameters were also correlated with the indicators of chronic disease and concurrent stress.

Principal component analysis and clustering
To better understand the haemodynamic differences between the subjects with both hypertension and T2D and the rest of the subjects, a principal component analysis (PCA) was performed on the subject-specific median values of the 28 optimized parameters. Using the unsupervised machine learning method from Jones et al. (2021), convex hulls of the four disease groups were created, and the subjects were clustered into two groups to account for one hypothetically more healthy and one hypothetically more diseased group. Before the PCA and clustering were applied, the parameters were standardized by centring and dividing each parameter value by the SD of that parameter. Both k-means and hierarchical clustering were used, and clustering groups were formed where both methods clustered the subjects into the same group (clusters 1 and 4), and where the two methods clustered the subject into different groups, so-called non-consistent clustered (NCC) (clusters 2 and 3). Finally, the cluster groups were combined with the HT+T2D convex hull, and all subjects within the same cluster group as well as the same HT+T2D hull ended up in the same final clustered group. Thus, the final clustered groups were: r HT+T2D group: all HT+T2D-subjects. r HT+T2D-like cluster group: non-HT+T2D-subjects within the same k-means and hierarchical clusters as the majority of HT+T2D subjects, and within the HT+T2D PCA hull.
r Non-HT+T2D-like cluster group: non-HT+T2Dsubjects within different k-means and hierarchical clusters from the majority of HT+T2D subjects, and outside the HT+T2D PCA hull.
r Non-consistently clustered group (NCC): • Any non-HT+T2D subject within different k-means and hierarchical clusters (clusters 2 and 3). • Any non-HT+T2D subject within the same k-means and hierarchical clusters as the majority of HT+T2D subjects, but outside the HT+T2D PCA hull. • Any non-HT+T2D subject within the HT+T2D PCA hull but in a different cluster than the majority of HT+T2D subjects.

Group differences in cohort characteristics
The cohort characteristics are presented in Table 2. Differences between any of the four groups were seen in weight, BMI, HR, DBP before MRI, systolic home BP, diastolic home BP, capillary P-glucose, HbA1c and low-density lipoprotein (LDL). Even though the groups were created based on home blood pressure, the SBP before MRI did not differ significantly between groups, and DBP before MRI was only higher in subjects with HT+T2D compared to controls, but not between any other groups. The heart rate during MRI was lower in controls compared to all other groups, whereas there was no significant difference between the three patient groups.
In the two groups of subjects with diabetes, there was no significant difference in cohort characteristics except for the difference in home blood pressure. In the two groups of subjects with hypertension on the other hand, except for the difference in capillary P-glucose and HbA1c, the BMI was higher and the LDL was lower in the HT+T2D group compared to the HT group. The LV E/e ratio was significantly higher in HT+T2D compared to controls, but there was no difference in the E/A ratio between the groups.

The subject-specific cardiovascular models fit to data
For each of the 80 subjects, a subject-specific cardiovascular mathematical model was created and fitted to subject-specific MRI and blood pressure data. All personalized models agree with the patient-specific estimation data of blood flows, blood pressure and five haemodynamic parameters (Fig. 2). The median of the subject-specific root mean square error (RMSE) of the simulated blood flows in the pulmonary veins, mitral valve, aortic valve, and ascending aorta compared to data was 24.1 mL, which is similar to our previously reported RMSE in healthy controls (Casas et al., 2018) and reasonable compared to the blood flow values with peaks of ∼100-500 ml. The median RMSE for the blood pressure is 0.9 mmHg, which was smaller than the expected SD of the brachial cuff measurements (Table 3). The model simulations compared to blood flow data of three representative subjects and the agreement with parameter data and brachial pressure data for all subjects are shown in Fig. 2. All subject-specific simulations compared to blood flow data are shown in the Supporting information ( Figure S2). Finally, the model prediction of left ventricular stroke volume and ejection fraction agreed well with the 3D cine MRI-based stroke volume data not used for model training (Fig. 3D and E).
The cine 3D-based stroke volume is slightly higher  than the flow-based model prediction as a result of the commonly found difference between flow-based and cine-based stroke volume data. Because of the lack of volume data in the model training, the large model uncertainty of the end-diastolic volume is propagated to the uncertainty in the ejection fraction calculated as stroke volume/end-diastolic volume, whereas the model uncertainty of stroke volume and other model predictions is smaller. As reference values, the model simulations compared to 4D flow-based measurements of stroke volume at the mitral valve, aortic valve and ascending aorta are seen in Fig. 3A-C. The model agreement with flow-based stroke volume data was good even though the maximum difference in stroke volume between the data in the three locations was 16 mL with respect to median (Table 3), making it impossible for the model to have a perfect fit to all three stroke volumes at the same time. The median of the subject-specific sum of the  absolute difference between data and simulation in the three locations was 21.4 mL.

Haemodynamic comparison between controls, subjects with hypertension and subjects with T2D
The subject-specific median parameter values were compared between the four subject groups and differences (P < 0.05) in parameter values were found between any of the groups in 11 of 28 parameters (Fig. 4). After Benjamini-Hochberg correction for several comparisons, three significant differences remained, one in the pulmonary veins and two in the left ventricle. The remaining differences were between controls and T2D in the parameter Cpvc (P = 0.006), corresponding to the capacitance of pulmonary capillaries and veins, and, in ksystLV (P = 0.016), the systolic time constant of the LV. The systolic time constant was also higher in T2D than in HT+T2D (P = 0.005). Finally, m2 LV, the relaxation rate constant of the LV, was lower in HT+T2D compared to controls (P = 0.002), HT (P = 0.002) and T2D (P = 0.013). No significant differences were found between HT and T2D, nor between controls and HT.   The subject-specific uncertainty for each parameter was large compared to the group uncertainty for the LV systolic time constant (127%, 109% and 132% for the C, T2D and HT+T2D groups), intermediate in the capacitance of pulmonary capillaries and veins (100% and 93% for the C and T2D groups), and smallest in the LV relaxation rate constant (59%, 43%, 126% and 69% for the groups C, T2D, HT and HT+T2D) [see Supporting information (Table S2); median values of the subject-specific interquartile ranges compared to the interquartile ranges of each group in %].
The only group difference in predicted ventricular and atrial pressure between the four groups was in the maximum left ventricular pressure, which only differed between the controls and all other groups (Fig. 5).

Clustering of HT+T2D-like and non-HT+T2D-like subjects
To analyse how the haemodynamic parameters differed in subjects with co-existing hypertension and T2D, clustering and PCA were performed. In the PCA, the first two principal components explained 35.6% of the variance in the data, where PCA1 explained 20.3%, PCA2 explained 15.2%, and all following components explained less than 15% (see Supporting information, Figure S3). The PCA coefficients were spread out among almost all of the 28 optimized parameters (see Supporting information, Table S3), where neither m2LV, ksystLV, nor Cpvc stands out from the rest. The k-means and hierarchical clusters did not match well with the PCA-based complex hull of any of the disease groups; neither visually ( Fig. 6A), nor computationally. When the HT+T2D group was compared to the clusters, 8 HT+T2D subjects were in cluster 1 (chosen as the HT+T2D-like cluster), 7 HT+T2D subjects were in cluster 4 (chosen as the non-HT+T2D-like cluster) and 2 HT+T2D-subjects were classified differently between the clustering methods and thus were classified as NCC. The NCC group including all non-consistently clustered non-HT+T2D-subjects contains 44% (n = 35) of all 80 subjects. Only 10 subjects ended up in the non-HT+T2D-like cluster-based group, and 18 ended up in the HT+T2D-like cluster-based group. All the resulting groups can be seen in Fig. 6C, where the non-HT+T2D-like subjects ended up close to the area covered only by the control hull, whereas the HT+T2D-like subjects ended up closer to the HT-hull and the T2D-hull. When plotting the groups against the left ventricular relaxation rate and the pulmonary venous compliance, which both had significant differences between the disease-based groups, the HT+T2D-hull stands out in the LV-relaxation rate axis (Fig. 6D). Additionally, many of the non-HT+T2D-like subjects ended up in the HT+T2D-hull, indicating that based on these two parameters, some of the non-HT+T2D-like subjects were instead similar to the HT+T2D subjects.

Haemodynamic parameter correlations with heart rate, blood pressure and diabetic status
To investigate the haemodynamic mechanisms behind the group differences, correlation tests between model parameters and indicators of chronic disease and markers of stress concurrent with acquisitions were performed (Fig. 7). Both the pulmonary venous compliance and the LV systolic time constant were slightly correlated with heart rate (P < 0.001, r 2 = 0.22 and P = 0.016, r 2 = 0.07). Additionally, pulmonary venous compliance was correlated with SBP during MRI (P = 0.017, r 2 = 0.05) and HbA1c (P = 0.047, r 2 = 0.01). The LV relaxation rate, m2 LV, was not correlated with any of the variables, and none of the three parameters were correlated with diabetes duration or home blood pressure. Furthermore, both aortic resistance and compliance were weakly correlated with heart rate (P = 0.001, r 2 = 0.11 and P = 0.003, r 2 = 0.07) and SBP during MRI (P < 0.001, r 2 = 0.24 and P < 0.001, r 2 = 0.20). Aortic compliance was correlated with HbA1c (P = 0.036, r 2 = 0.05) and systolic (P < 0.001, r 2 = 0.1) and diastolic home blood pressure (P = 0.003, r 2 = 0.09). Finally, the LV diastolic time constant (kdiastLV) had a clear correlation with heart rate during MR (P < 0.001, r 2 = 0.67). The maximum elastance (EmaxLV) was more weakly correlated with heart rate (P < 0.001, r 2 = 0.12) and SBP during MRI (P < 0.001, r 2 = 0.10). Neither the ventricular contraction rate (m1LV), nor relaxation rate (m2LV) was correlated with any of the variables.

Subject-specific predictions of haemodynamic function
As an example of the subject-specific variations and haemodynamic differences, the haemodynamics of four subjects from each of the four groups were examined (Fig. 8). The subjects were selected as the subjects with the best fit to data among those with matching hypertension levels at home and before MRI in each of the four groups. The two hypertensive subjects both had higher pressure in the aorta and ventricle compared to the normotensive groups, as expected. The time-varying elastance was also clearly higher in the hypertensive subjects, at least for the HT subject. By contrast to the group-level differences, no clear reduction in ventricular relaxation was seen in the subject with HT+T2D compared to the other subjects even though all four subjects had different blood flow patterns at the mitral valve. The pressure and time-varying elastance in the left atrium both had a larger model uncertainty than the other predictions, with no clear differences between the four subjects. The aortic and mitral blood flow volumes and the aortic and ventricular pressures had smaller model uncertainty, which was expected because they were either directly or indirectly fitted to the data. To investigate the effect of impaired left ventricular relaxation (the main haemodynamic difference between the four groups), a model prediction of decreased LV relaxation rate in the control subject was conducted. The effect of impaired relaxation was mainly seen in the mitral blood flow with a relatively shorter diastole with a shorter early filling, the early and late peaks closer in time, longer isovolumic relaxation time, and increased flow during late filling (Fig. 9). Additionally, the slower relaxation directly affected the time-varying elastance function of the left ventricle (Fig. 9C), where the slope was slower in the predicted 'diseased' state compared to the original healthy control subject. The same pattern appeared in the ventricular pressure (Fig. 9F), which took a longer time to reach its minimum value. The other model predictions (time-varying elastance in the LA, blood flow at the aortic valve, aortic pressure, and LA pressure) were not much affected by the reduced relaxation rate.

Discussion
Using subject-specific cardiovascular mathematical models, we show group differences in haemodynamics related to diastolic dysfunction in subjects with hypertension, subjects with T2D, and subjects with both hypertension and T2D compared to controls. However, the subject-specific cardiovascular models, created in 80 subjects, also reveal that the haemodynamic changes behind hypertension and T2D are complex, highly individual and cannot be explained only by diastolic dysfunction. Instead, a more comprehensive subject-specific approach such as the cardiovascular model used in the presented study is needed.

Left ventricular relaxation rate as a marker for early diastolic dysfunction in HT+T2D subjects
The significantly reduced left ventricular relaxation rate in the subjects with both hypertension and T2D (Fig. 4) causes impaired diastolic filling and is a sign of early diastolic dysfunction (Nishimura et al., 1997). The diastolic dysfunction in HT+T2D is further reflected in the LV E/e ratio, another marker for diastolic dysfunction, which is significantly higher in HT+T2D compared to controls (Table 2). Following a reduction in left ventricular relaxation rate, the next step in the typical progression of diastolic dysfunction is an increase in left atrial pressure and left ventricular filling pressure (Nishimura et al., 1997). However, no increase in atrial pressure can be seen when comparing the predicted mean left atrial pressures between the HT+T2D group and the other groups (Fig. 5). Additionally, there is no additional effect of HT+T2D on the ventricular or aortic pressure ( Fig. 5 and Table 2). The absence of an increase in atrial pressure in this HT+T2D group with a decreased relaxation rate underscores that these patients are probably at an early stage of diastolic dysfunction. This aligns with the clinical setting, as the subjects were recruited from a relatively healthy general population and their disease duration is probably relatively short given their age. The adverse effect of T2D and hypertension on ventricular function has been confirmed by previous studies. In the Strong Heart study, both T2D and hypertension impacted LV relaxation, and their combination had a larger impact on LV relaxation than either process in isolation (Liu et al., 2001). One explanation for this is more extensive fibrosis in the left ventricle that has been seen in subjects with combined disease (van Hoeven & Factor, 1990). Previous studies have indicated that pulse pressure, a marker of arterial stiffness, in hypertensive patients is a predictor of the development of diabetes (Yasuno et al., 2010), whereas a Mendelian randomization study suggests that T2D can cause hypertension, but not the other way around (Sun et al., 2019).
In the present study, the E/e values are larger in HT+T2D compared to the controls (Table 2). Additionally, an increased E/e ratio in T2D compared to non-T2D was previously found in this cohort, suggesting that diabetes affects diastolic dysfunction (Edin et al., 2022). Biological explanations for the haemodynamic effects of diabetes include endothelial dysfunction, insulin resistance and hyperglycaemia, which can lead to the stiffening of both large and small vessels as well as the left ventricle (Petrie et al., 2018;Tan et al., 2020). Diastolic dysfunction has also been associated with the amount of intramyocardial fat, which is increased in T2D (Rijzewijk et al., 2008). However, in the present study, the LV relaxation rate is not correlated with markers of chronic disease nor with concurrent blood pressure or heart rate (Fig. 7), which indicates that it is the combined burden of hypertension and T2D affecting ventricular relaxation. Either way, changes in ventricular function are related to cardiovascular disease in both hypertension and T2D (Ernande et al., 2017) and the CVD risk is increased in co-existing T2D and hypertension (Chen et al., 2011;Verdecchia et al., 2004).
The patient-specific haemodynamic effects of reduced ventricular relaxation are illustrated in Fig. 9, where a clear impaired relaxation is seen in the blood flow at the mitral valve The control subject already has a higher peak during late filling and a flow volume-based E/A ratio below 1 as is often present in patients at this age, although this pattern is further exaggerated as impaired relaxation worsens. This impaired relaxation pattern is typical in early diastolic dysfunction (Nishimura et al., 1997). As seen in the group comparisons, in this patient there is no increased atrial pressure in the model prediction, again indicating an early stage of diastolic dysfunction. Further studies would be needed to investigate how T2D and hypertension affect subject-specific ventricular function in different disease stages.

Reduced relaxation rate is not enough to explain the complex and subject-specific haemodynamics in T2D and hypertension
Both the correlations with markers of long-term disease, the PCA, the clustering and the subject-specific model predictions show that, even though reduced relaxation rate is one of the most prominent haemodynamic differences found in this group of subjects, the haemodynamic changes in T2D and hypertension are more complex and subject-specific than simply a decrease in relaxation rate. The clustering resulted in two new groups of subjects that do not exactly correspond to any of the disease-based groups ( Fig. 6A and B). Rather, it shows that there is a gradual change from healthy to the patterns found in T2D and hypertension and that some of the subjects from the most diseased HT+T2D group share haemodynamic patterns with some of the subjects from the control group. This can be a sign that not only some of the control subjects were closer to developing hypertension and/or T2D, but also some of the subjects with T2D and hypertension in some aspects had less derangement of their haemodynamic patterns, either as a result of an early stage of disease progression or because of compensatory mechanisms such as remodelling of the ventricle. Nonetheless, neither the PCA, nor the clustering methods could separate the four disease groups based on the parameter values. However, other clustering and visualization methods or other cardiovascular variables such as ventricular geometry or mass might be better at explaining the complex variation in this population.
The complex variation was demonstrated in the PCA analysis where ventricular relaxation was only one of multiple haemodynamic parameters that together explained the haemodynamic variation in these subjects (see Supporting information, Table S3). One additional important parameter is aortic compliance, which did not significantly differ between the groups. Aortic compliance was slightly correlated with HbA1C, which indicated that T2D affected aortic stiffness (Fig. 7C). Furthermore, aortic compliance is the only one of the eight selected parameters that correlated with home blood pressure, indicating increased arterial stiffness in hypertension ( Fig. 7E and F). Increased arterial stiffness with both T2D and hypertension is in line with several previous studies, where an increased effect of co-existing T2D and hypertension has been reported (Eren et al., 2004;Nuamchit et al., 2020;Tedesco et al., 2004). Additionally, Eren et al. (2004) found a relationship between aortic stiffness and left ventricular diastolic function in patients with T2D, hypertension, and both T2D and hypertension.
The complexity and subject-specific variation in haemodynamic mechanisms could, to some extent, be further illustrated when patient-specific haemodynamic predictions were analysed. The aortic pressure, ventricular pressure and time-varying elastance were, as expected, increased in the two selected hypertensive subjects compared to the control and the subject with T2D (Fig. 8C, D and F). However, when analysing the single subjects, not all group differences were seen. As an example, there was no clear diastolic dysfunction in the selected subject with both hypertension and T2D compared to other subjects (Fig. 8A). Compared to the group median, some subjects were maybe at different stages of hypertension and diabetes, more or less adequately treated, and/or with varying amounts of acute stress at the time of the MRI scan. Finally, the complexity of haemodynamic effects in T2D and hypertension was illustrated in the subject-specific model prediction (Fig. 9), where changing only the ventricular relaxation did not impact blood flows and pressures except for the flow at the mitral valve. This showed that, on a subject-specific basis, it was not only the relaxation rate that was affected, but also a combination of several haemodynamic parameters that together caused a modified overall haemodynamic pattern in HT and T2D.
Finally, the correlation in Fig. 7 showed the importance of concurrent haemodynamics, where for example the concurrent heart rate and blood pressure had a larger impact on model parameters than the markers of chronic disease in resting conditions. Additional chronic group differences might be blunted by subject-specific differences in stress during the acquisition resulting in increased heart rate and blood pressure. Even though group differences in heart rate and blood pressure are expected, especially between the HT group and the T2D group, almost no group differences in concurrent blood pressure were found, and the heart rate was similar between the disease groups (Table 2). This indicates that subject-specific increases in blood pressure and heart rate during the scan might obscure the group-level differences in parameters such as aortic stiffness or resistance, which both correlate slightly with heart rate and systolic blood pressure in the scanner (Fig. 7A and  B). Additionally, other heart rate-related parameters like the maximum ventricular elastance or the LV diastolic time constant could have chronic effects that are blunted by concurrent effects in the scanner. To further reveal those group differences that are affected by concurrent haemodynamics, longitudinal studies of patient-specific disease progression or stress interventions such as exercise might be needed.

Limitations
This work has some limitations. The model used in the present is a simple description of the haemodynamics in the heart and it lacks some of the more complex geometrical-and spatial-dependent mechanisms known to affect blood flow and blood pressure. To ensure that the model predictions are correct, validation with data not used in the model training can be used. The model predictions of pulmonary, atrial and ventricular pressure were not validated experimentally here because this would require invasive and undesired methods such as inserting a catheter into the heart or aorta, which is not ethically motivated in these relatively healthy subjects. Instead, the model was trained using the non-invasive brachial pressure measurement to ensure that the simulated aortic pressure was within a reasonable range. To further improve the estimate of the left atrial haemodynamics, we included measurements of blood flow in the pulmonary veins, which were not used in previous publications of this model (Casas et al., 2017(Casas et al., , 2018. Additionally, the model was trained on the blood flow curves, which contain indirect information about the stroke volume, but the stroke volume of the left ventricle was left out and compared to the volumes in the model. However, to further validate the model, more measurements are needed. In the present study, the uncertainties of MRI and blood pressure measurements were not considered. MRI suffers from a range of errors such as non-zero background phase offset and merging several heartbeats during the 8 min acquisition into one heartbeat . Moreover, the brachial cuff pressure is known to vary both in the presence of medical staff and with factors such as stress or time of the day, as well as body posture (Stergiou et al., 2021). Because only single measurements of blood flow and blood pressure were performed for each subject, the scan-to-scan variation for each subject was unknown. J Physiol 601.17 Furthermore, both the pulmonary venous flow and the aortic compliance data may have had larger measurement uncertainty compared to the other blood flow volumes and parameters, which could explain why the model was less good at explaining those data as well. As an example of the model output sensitivity to training data, the sensitivity of model parameters and predictions to measured brachial pressure is evaluated in the last section of the Supporting information, where the sensitivity to measured SBP of the ventricular relaxation rate constant is only −0.6 %/mmHg, resulting in a model overestimation of 9.6% or 3.2 (unitless) if the measured SBP is underestimated by 16 mmHg.
Except for the uncertainty in our model and the experimental data, there are additional factors that could have affected the haemodynamic results. The cross-sectional design of the study limits the ability to understand the subject-specific effects of T2D and hypertension. To further understand these effects, a longitudinal study following each patient would be ideal, although demanding, to perform. Another potentially important factor that can affect the haemodynamics of hypertensive subjects is anti-hypertensive medication, and the subjects with T2D may use glucose-and lipid-lowering medication. In the present study, medication information was not available, and therefore not included. Additionally, no sex differences were considered in the study. There is no significant difference with respect to sex between the groups, although the majority of the HT+T2D subjects (76%) were men, whereas only 67% of the controls, 47% of subjects with T2D and 64% of the subjects with HT were men. Furthermore, the subject-specific model uncertainty was not included in the comparative statistics. However, the relative interquartile range of individual subjects compared to the interquartile range of the groups was relatively low for the ventricular relaxation rate (see Supporting information, Table S2). Finally, even though the total number of 80 subjects is relatively large, the number of subjects within the four subgroups was smaller. The conclusions based on the group comparisons should therefore be interpreted as exploratory rather than generalizable conclusions.

Conclusions
The combination of a subject-specific cardiovascular model and 4D flow MRI data demonstrated that the haemodynamic mechanisms in hypertension and T2D are complex, subject-specific and dependent on concurrent haemodynamics. The main group-level haemodynamic effect visible in a resting state was an early diastolic dysfunction in subjects with co-existing hypertension and T2D compared to controls and subjects with T2D or hypertension alone. However, the subject-specific cardiovascular mathematical model showed that a more comprehensive and subject-specific analysis is needed to fully understand the patient-specific haemodynamic mechanisms in T2D and hypertension, including a gradual change of a combination of several subject-specific haemodynamic parameters such as aortic stiffness, heart rate, and ventricular contraction and relaxation. These new insights into group-level and subject-specific haemodynamics could, together with a personalized model approach and further studies, aid in the treatment planning of patients with both T2D and hypertension.