Advancement of subject-specific in silico medicine requires new imaging protocols tailored to specific anatomical features, paired with new constitutive model development based on structure/function relationships. In this study, we develop a new dual-velocity encoding coefficient (VENC) 4D flow MRI protocol that provides unprecedented spatial and temporal resolution of in vivo aortic deformation. All previous dual-VENC 4D flow MRI studies in the literature focus on an isolated segment of the aorta, which fail to capture the full spectrum of aortic heterogeneity that exists along the vessel length. The imaging protocol developed provides high sensitivity to all blood flow velocities throughout the entire cardiac cycle, overcoming the challenge of accurately measuring the highly unsteady nonuniform flow field in the aorta. Cross-sectional area change, volumetric flow rate, and compliance are observed to decrease with distance from the heart, while pulse wave velocity (PWV) is observed to increase. A nonlinear aortic lumen pressure–area relationship is observed throughout the aorta such that a high vessel compliance occurs during diastole, and a low vessel compliance occurs during systole. This suggests that a single value of compliance may not accurately represent vessel behavior during a cardiac cycle in vivo. This high-resolution MRI data provide key information on the spatial variation in nonlinear aortic compliance, which can significantly advance the state-of-the-art of in-silico diagnostic techniques for the human aorta.
Of the 97,000 km of blood vessels in the human body, the near-meter long segment connecting the left ventricle to the periphery, known as the aorta, is the most important. Diseases affecting the aorta such as aneurysm and dissection have long been documented, but still today remain difficult to treat. According to the Centre for Disease Control and Prevention, an average of 47,000 deaths each year in the United States are attributed to diseases of the aorta and its branches (excluding carotid and coronary disease). This exceeds the number of annual deaths due to breast cancer, pancreatic cancer, colon cancer, and prostate cancer .
Patients undergoing surgery of the aorta have two main treatment options: open surgical repair (OSR) or endovascular aortic repair (EVAR). The introduction of EVAR in the early 1990s was fueled by the need for a less invasive treatment option for comorbid patients and poor outcomes following OSR. In the quarter-century since its introduction, EVAR has shown superiority over OSR in the short-term, where studies continue to report mortality rates from 14% to 45% in the first 30 days post-OSR [2,3], but no significant benefits are apparent for EVAR patients in the long term .
As the first thoracic endovascular aortic repair (TEVAR) graft only received Food and Drug Administration approval in 2005 , long-term results are only now coming to light. A number of studies have reported high levels of cardiac complications following TEVAR, where Conrad  reports 34% mortality due to cardiac events in thoracic aortic aneurysms, while a study by Bischoff  reports 30% cardiac mortality for a larger thoracic aortic aneurysms cohort. A recent study by Concannon et al.  reports that, from a cohort of 151 patients with thoracoabdominal aortic aneurysms, 39% of total deaths were due to cardiac failure. Notably, all deaths due to new onset cardiac complications were in patients who underwent stenting of the supradiaphragmatic aorta. Altogether, these results suggest a dependence of postoperative cardiac outcomes on the location of stent deployment in the aorta.
A detailed biomechanical investigation of the influence of stent deployment on aortic deformation, hemodynamics, and pulse wave velocity (PWV) is required to uncover the mechanisms that cause cardiac complications post-TEVAR. As an important first step in this process, we propose a noninvasive 4D flow MRI protocol to accurately characterize spatial variations in biomechanical behavior throughout the entire aorta, in addition to dynamic variations throughout a cardiac cycle. The ability to characterize spatially dependent vessel geometry and deformation, blood flow patterns, and PWV will potentially guide the selection of stent-graft design and position in EVAR procedures in order to minimize the risk of cardiac complications postintervention. An increased PWV has been established as a strong risk factor for cardiac events, independent of traditional risk factors such as smoking, hypertension, and diabetes mellitus . The ability to accurately determine the spatially nonuniform PWV throughout the entire aorta, both pre- and postintervention, could potentially provide new insights.
The increase in clinical acceptance of EVAR has resulted in a reduction in the number of primary OSR cases , and a consequent reduction in the availability of tissue samples for in vitro biomechanical testing. Moreover, surgically excised tissue often consists of a small portion (approximately 1 cm2) of the aorta, presenting significant challenges in terms of bi-axial mechanical testing . Therefore, in vitro testing of excised tissue does not present a viable methodology to accurately determine the detailed spatial variations in compliance and PWV in a patient-specific aorta. Alternative approaches of combined medical imaging and computational analysis (finite element and computational fluid dynamics (CFD) modeling) to determine biomechanical properties noninvasively are highly promising, particularly in light of recent advances in medical imaging technology and computational capability.
Of the few studies that attempt to investigate the biomechanics of the aorta, its heterogeneity has been reasonably well established in animals through ex vivo testing of the excised vessel [12,13]. Previous in vivo analyses of the human aorta have focused on limited isolated segments, such as the thoracic  or abdominal aorta , which fail to provide the necessary anatomical coverage to capture the true heterogeneity and therefore cannot be taken to represent the entire vessel. Due to the lack of reliable and detailed information on the heterogeneity of the aorta, computational models have typically assumed that the wall stiffness is spatially uniform throughout the vessel [16–19]. An improved robust methodology to noninvasively characterize patient-specific spatial variation in aortic PWV and compliance throughout the cardiac cycle has the potential to provide accurate heterogeneous material properties for computational models, leading to significant improvements in EVAR device design, and subsequently, postoperative outcomes.
In silico tools are being considered as possible replacements for animal and human experimentation and the preclinical assessment . Advancement of subject-specific in silico medicine requires new imaging protocols tailored to specific anatomical features, paired with new constitutive model development based on structure/function relationships. In this study, a dual-velocity encoding coefficient (VENC) 4D flow MRI protocol is developed to achieve accurate measurement of the dynamically changing flow velocity field and lumen area throughout the entire cardiac cycle and the entire aorta. To the best of our knowledge, no previous medical imaging paper has reported such detailed spatial and temporal characterization of the human aorta. To date, 12 aortic dual-VENC 4D Flow MRI studies exist in the literature, four of which pertain to phantom geometries [21–24], while the remainder are focused on a single isolated segment of the aorta such as the ascending thoracic [25–31]. A nonlinear relationship between lumen area and pressure is observed in vivo over the duration of a cardiac cycle throughout the entire aorta, suggesting that aortic biomechanics may not be accurately characterized by a single value compliance coefficient, as commonly assumed [32–34]. Furthermore, our detailed in vivo measurements reveal that the lumen pressure–area relationship, and PWV are highly heterogeneous along the aortic length.
In this paper, a protocol is proposed to evaluate patient-specific hemodynamics and lumen deformation along the entire human aorta, and throughout the entire cardiac cycle, using phase-contrast magnetic resonance (PC-MRI) principles (specifically, 4D flow MRI). Further details of the applications and potential uses of 4D flow MRI can be found in: Refs. [35–40]. Generally, with the aim of assessing anatomical structures, it is the magnitude of the local spin magnetization vector that is used in the creation of typical magentic resonance images. However, important information regarding the movement of hydrogen protons is encoded in the phase of this vector. In the field of PC-MRI, such information is exploited to determine the flow velocity of targeted protons. A brief summary of the theoretical background to PC-MRI is presented in Sec. 2.1 to motivate the protocol proposed in this paper. Further details pertaining to spin dynamics and velocity encoding sensitization can be found in Appendix A.
2.1 Theoretical Background.
In this section, we provide a brief overview of the key theory and equations that motivate the dual-VENC protocol proposed in Sec. 2.2. The theoretical physics underlying MRI is extensively outlined in literature, e.g., Refs. [41–43]. In summary, MRI is a phase-sensitive modality that encodes information regarding the velocity of the targeted protons into the detected signal. The velocity is proportional to the phase of the local transverse magnetization vector. In the remainder of this paper, the term spins is used to refer to a finite group of protons within a given volume. In Fig. 1, cranial flow (in the positive z-direction) is indicated by positive grayscale values on the foot-head (FH) image in the ascending thoracic aorta. The flow direction is in the negative z-direction in the descending aorta, as indicated negative grayscale values in the FH image. In the anteroposterior (AP) image, posterior flow can be seen traversing the apex of the aortic arch while anterior flow is indicated by negative grayscale values as blood leaves the left ventricle into the ascending thoracic aorta. Similarly, flow sensitization is seen with the right-left (RL) image although velocity encoding is less obvious in the RL direction upon viewing a sagittal plane.
Equation (4) dictates how the scanner can sensitize to specific fluid velocities. For example, reducing the VENC from 200 cm/s to 50 cm/s requires a four-fold increase in the strength of , or an increase in the time over which it is activated. Thus, it is preferable that the strength of the magnetic field gradient be increased to achieve a reduction in velocity sensitization, instead of increasing T necessitating unfeasibly long scan times.
2.2 Proposed Dual-Velocity Encoding Coefficient Protocol for Complete Characterization of Aortic Flow.
Maximal sensitivity is obtained for spins moving at a velocity equal to the specified VENC value. This presents a particular challenge for determination of blood flow patterns in the aorta where flow is highly unsteady (temporally varying) and nonuniform (spatially varying). For example, a VENC of 200 cm/s may provide a suitable level of sensitivity to determine the high velocity blood flow patterns in the aortic arch during systole. However, such a VENC value is not suitable during the diastolic phase, where the fluid velocity is considerably lower. In fact, in using a VENC of 200 cm/s, low velocity blood flow during diastole cannot be distinguished from static tissue and the lumen of the aorta cannot be reliably identified. A reduced VENC is required to achieve sufficient resolution of the flow field during diastole.
Of course, such a low VENC is not suitable for systolic flow velocities; any fluid velocity greater than VENC will be misrepresented and aliased, as described elsewhere [44–46]. In an attempt to overcome this issue, previous studies have proposed phase unwrapping algorithms to estimate velocities higher than VENC. However, significant errors have been reported for such techniques, in addition to increased postprocessing time [47–49].
The dual-VENC protocol proposed in this study generates a composite dataset, with a high-VENC of 200 cm/s targeted to systole and a low-VENC of 50 cm/s targeted to diastole. As the only difference between our two datasets is the velocity sensitization, accurate velocity field measurement and lumen boundary isolation can be performed for each phase and plane throughout the entire cardiac cycle, all the while keeping acquisition parameters within the bounds specified in the most recent 4D flow MRI expert consensus statement . If the velocity of any pixel in an arbitrary plane of interest is greater than our low-VENC value (50 cm/s), we use the high-VENC matrix to calculate the cross-sectional area and volumetric flow rate; otherwise, we use the corresponding low-VENC matrix. This approach negates the need for phase unwrapping techniques and provides greater accuracy in flow quantification in areas where the fluid velocity is low than single high-VENC acquisitions.
2.3 Imaging Parameters.
The current study was approved by the institutional review board (Research Assessment Group (RAGp), Galway Clinic) and was conducted on a healthy 25-year-old male with a normotensive blood pressure measurement of 117/73 mmHg and a heart rate of 60 bpm. The subject was placed in a Philips Ingenia 3 T MRI scanner (Philips Medical Systems, Best, Netherlands) and a four lead electrocardiogram system was placed on the chest with retrospective synchronization to the scanner to image according to specific phases of the subjects' cardiac cycle. A noncontrast RF-spoiled gradient echo pulse sequence was employed in order to capture a sufficient number of heart phases under free-breathing conditions. The field of view was set to encompass the entire aorta. The longitudinal (FH) boundaries spanned from above the level of the aortic arch to distal to the common iliac bifurcation, while the lateral (AP) and (RL) bounds enclosed the breadth and width of the subject, respectively. The frequency encoding direction was set to AP to reduce artifact from respiratory motion. Important scan parameters are as follows: repetition time = 3.1 ms, echo time = 1.9 ms, flip angle = 8 deg, cardiac phases = 20, temporal resolution = 50 ms, isotropic in-plane resolution = 1 mm, slice thickness = 4 mm, VENC = 200 cm/s and 50 cm/s. VENC scouts were ran to obtain the minimum high-VENC value to prevent aliasing and optimize signal to noise ratio, while a 4 mm slice thickness was used to limit scanning time. The scan time for a VENC of 200 cm/s was 4 min, and 8 min for a VENC of 50 cm/s. In the case of the latter, the repetition time was increased to 10 ms to allow sufficient down-time for the gradient coils to prevent excessive overheating. A balanced four-point encoding scheme was used, further details of which can be found in Ref. .
All data were processed using in-house developed C++, python, and matlab code. Data processing was performed on an Intel Core i7 CPU with 16GB DDR3 RAM. Postprocessing time for the dual-VENC dataset was approximately 20 min. Raw MRI data files were sorted according to their encoding direction using radiantdicomviewer (v4.2.1, Medixant, Poznan, Poland) and subsequently organized according to the time-point in the cardiac cycle using a custom imagej plugin (Fig. 2(a)). paraview (5.4.1)2 visualization software served as the platform for reading the image data for each encoding direction, developing voxel associativity and subsequent calculation of local velocity magnitudes as illustrated in Fig. 2(b). It is necessary to ensure that observational planes are orthogonal to the mean direction of blood flow when attempting to characterize lumen deformations and volumetric flow rates due to the onset of a pressure pulse. A centerline detection algorithm was developed to ensure such requirements were fulfilled as shown in Fig. 2(c), where the centerline is defined as the centroid of the aortic flow domain.
Analyses were performed at 10 planes along the aorta (Fig. 2(d)), ranging from Plane 1 distal to the sinus of Valsalva to Plane 10 immediately proximal to the common iliac bifurcation, with an average section spacing along the centerline of 50 mm. Each time-point for each plane in both VENC datasets was then exported for all further postprocessing in matlab (R2017b, MathWorks Inc., Natick, MA). A bicubic interpolation algorithm is employed in order to attain further clarity for aortic lumen edge detection. A series of cubic splines were fit to the intensity values of individual pixels along both the x and y dimensions in a given plane (Fig. 2(e)) and the grid density was increased to 0.5 × 0.5 mm in-plane spatial resolution. Figures 2(f), (i) and (ii) highlight in-plane pixel data pre- and postinterpolation.
At this point, each velocity magnitude image is masked by the square of the corresponding magnitude (anatomical) image (Fig. 2(f) (iii)) in order to create a phasecontrast magnetic resonance angiography (PC-MRA) matrix according to methods described in Ref. . Using each PC-MRA image, the boundary of the fluid domain is isolated for each plane and phase of interest to determine the lumen area as a function of space and time based on a custom-built segmentation algorithm. Using the high-VENC velocity matrix, if the velocity of any pixel is greater than the low-VENC value (50 cm/s) we use this matrix to calculate the cross-sectional area and volumetric flow rate. Otherwise, the corresponding low-VENC matrix is used. The composite dataset generated by the dual-VENC protocol eliminates the need for any phase unwrapping techniques.
For each of the 10 planes analyzed, the percentage cross-sectional-area change () is defined according to , where subscripts “sys” and “dia” represent systole and diastole, respectively. After isolating the aortic lumen from surrounding structures, streamlines and flow vectors can be plotted as shown in Fig. 2(g). The integral of the velocity within the boundary of the aortic lumen provides the instantaneous volumetric flow rate (Q) as shown in Fig. 2(h). Finally, the local PWV at each plane is calculated according to the flow area method described in Vulliémoz  and shown in Fig. 2(i), where PWV is defined as the coefficient of proportionality between Q and bound by the systolic upstroke of the cardiac cycle. Additionally, spatial variance in PWV is calculated using the time-to-peak (TTP) method described in Ref. , where in this case the wave speed is defined as , where is the distance along the vessel centerline between regions of interest and is the time lag between flow peaks for the thoracic and abdominal aortic segments in this case.
To estimate local vessel compliance, the pressure must be estimated throughout the cardiac cycle at the location in question. The clinical definition of vessel compliance is typically given as during a cardiac cycle. Typically, vessel compliance is reported as a single value [34,55] as only the systolic (SBP) and diastolic pressures (DBP) are recorded. However, it is trivial to demonstrate that, even for the simplistic thin-walled linear-elastic cylindrical vessel undergoing infinitesimal deformation, a linear relationship does not exist between and and therefore a single value of local compliance cannot be identified. Moreover, the well-established nonlinear material behavior of arterial tissue, e.g., Refs. [56–58], further invalidates the concept of a single value of compliance. By considering the entire blood pressure waveform, we investigate the time-dependence in local compliance along the length of the aorta.
where the first term on the left in Eq. (7) contains the integral of the local acceleration of a fluid particle along a streamline between points 1 and 2. is the fluid density, is the pressure, and is the fluid velocity. We neglect the last term on the right-hand side as the subject is in the supine position in the scanner and hence the change in elevation along the vessel can be taken as zero.
Finally, the systolic peak is approximated by a second-order polynomial, which satisfies continuity and produces the prescribed mean arterial pressure.
Finally, we perform a patient-specific CFD simulation that solves for the pressure gradient at the ascending and abdominal aortic levels. A finite element mesh was generated by sweeping surface elements along the aortic centerline through each of the 10 planes analyzed previously, using matlab (R2017b, MathWorks Inc., Natick, MA) and GIBBON  (Fig. 4). The inlet and outlet faces were closed off and the fluid domain was then filled with tetrahedral elements. A parabolic velocity profile derived directly from the MRI data at the aortic root was prescribed at the inlet. The outlet boundary condition consisted of a resistance (R) of 1.5 × 108Pa·s m−3, representing downstream vasculature. The no-slip condition was prescribed at the lumen boundary of the aorta, while fluid back-flow and tangential stabilization (β = 1) was defined at the outlet to deal with flow reversal and increase stability of the solution (). Blood was modeled as a Non-Newtonian Carreau fluid (), with a density of 1060 kg/m3. The resultant pressure waveforms for both the ascending thoracic and distal abdominal aorta are shown in Fig. 3(d). The corresponding area versus pressure graphs are shown for each method in Figs. 5(a)–5(h), respectively.
3.1 Dual-Velocity Encoding Coefficient Protocol for Complete Characterization of Aortic Flow.
We employ PC-MRI principles to capture both the deformation and hemodynamics of the entire aorta. The proposed dual-VENC protocol provides high sensitivity to all blood flow velocities throughout the entire cardiac cycle, overcoming the challenge of accurately measuring the highly unsteady nonuniform flow field in the aorta. A single high-VENC approach, while providing accurate measurements of high velocities during systole, was found to have insufficient resolution at low velocities to differentiate blood flow during diastole from the surrounding static tissue; this observation has been previously reported [25,39,63]. Consequently, the lumen geometry cannot be accurately determined in any region of the aorta during diastole, as clearly illustrated in Fig. 6(a) (only the high velocity flow in thoracic aorta at a time-point of 200 ms (systole) is accurately measured). An inability to accurately determine the lumen geometry and velocity field in the entire aorta for the entire cardiac cycle prohibits the determination of clinically relevant quantities such as cross-sectional area, aortic compliance, volumetric flow rate, and PWV. As discussed in Sec. 2.2, a single acquisition low-VENC will not provide accurate measurement of high velocities during systole due to phase wrapping. This is evident in Fig. 6(a), where the high velocities at 200 ms are significantly underpredicted by the low-VENC acquisition, compared to the high-VENC that is specifically sensitized for accurate measurement during systole. However, flow velocities and the flow domain are determined with greater accuracy at all other time-points (0, 500, 700, and 900 ms) using a low-VENC, in contrast to the high-VENC measurements where flow is indistinguishable from the noise associated with surrounding static tissue. As the velocity to noise ratio is proportional to the velocity and inversely proportional to the VENC, comparably lower velocities (such as those in diastole) are measured with reduced accuracy [25,64,65]. Both phantom  and in vivo  studies have shown that a dual VENC approach results in a more accurate flow quantification than single VENC acquisitions. Figure 6(b) further highlights this motivation for a dual-VENC approach. At 200 ms (top row), high-VENC accurately represents the fluid domain for the thoracic plane, whereas velocity aliasing is evident in low-VENC. In fact, some velocity vectors over 50 cm/s are misrepresented as negative velocities traveling toward the heart for low-VENC at this thoracic plane during systole. At 900 ms (bottom row), high-VENC is incapable of distinguishing the fluid domain from static tissue in the abdominal plane, while the low-VENC accurately represents the flow field and aortic lumen boundary.
3.2 Spatial Deformation.
Figure 7 shows the spatial and temporal change in lumen cross-sectional area throughout a cardiac cycle. Clearly the lumen cross- sectional area (CSA) decreases with increasing distance from the heart at any given time-point in the cardiac cycle. For example, at time t = 250 ms, the CSA at Plane 2 in the ascending aorta is 644 mm2, compared to 295 mm2 at Plane 5 and 158 mm2 at Plane 10.
Figure 8(a) shows the lumen area at the end of diastole, Adia, for all 10 planes. The well-known tapering of the aorta is also evident, with a decrease in Adia with increasing distance from the heart. The percentage change in cross-sectional area, , due to the onset of the pressure pulse is presented for each plane in Fig. 8(b), where . First, it should be noted that ranges from 15% for plane 10 up to 65% for plane 1, providing an indication of the extremely large deformation of the aortic wall during a cardiac cycle. Indeed, it should be noted that the circumferential strains in the aortic tissue will be significantly larger than the values of reported, given that the undeformed reference area (at zero pressure) is significantly lower than (clearly the true undeformed reference area cannot be determined in a “live” aorta, so the measure is instead presented here to demonstrate the high aortic deformations during a cardiac cycle). Categorizing the planes into two subgroups, namely, “thoracic” and “abdominal,” a statistically significant difference in is observed between the two groups (p < 0.005).
3.3 Spatial Hemodynamics.
The integral of the velocity matrix within the boundary of the fluid domain yields the volumetric flow rate, (Fig. 9). The reduction of with increasing distance from the heart can be attributed, in part, to out-flow to visceral arteries including the supra-aortic, mesenteric and renal vessels. For example, the large drop in flow between planes 1 and 5 is associated with out-flow to the innominate, left common carotid, and left subclavian arteries supplying the head, neck, and upper body with a large volumetric blood flow. The opening of the aortic valve occurs at approximately 50 ms and closes at 500 ms, while the time lag between the flow peaks of each plane is related to the speed of the ejected pulse wave propagating through the aortic tree. Peak systolic blood flow ranges from 196 mL/s at plane 1 in the ascending aorta to 28 mL/s at plane 10 in the abdominal aorta, while diastolic flow at time-point 600 ms ranges from 53 mL/s at plane 1–2 mL/s in plane 10. The nonzero flow during diastole, illustrates the well-known Windkessel effect. The measurements presented here demonstrate that the diastolic flow due to the Windkessel effect is highest in the ascending aorta and reduces with increasing distance from the heart.
The coefficient of proportionality between and provides the wave propagation speed (i.e., the speed of a blood column as it travels through the aorta following ventricular ejection), formally known as the PWV. Figure 10(a) shows a higher wave velocity in the distal aorta than in the proximal aorta. Again, categorizing the planes into two subgroups, thoracic and abdominal, a statistically significant difference in PWV between the two groups is found (p < 0.005). The implementation of the TTP method to determine the PWV also provides a similar result, as shown in Fig. 10(b); the PWV in the abdominal aorta (8.2 m/s) is found to be approximately 28% higher than in the thoracic aorta (6.4 m/s). The increased PWV in the abdominal aorta is due, in part, to the tapered geometry, as shown in Fig. 8(a). However, spatial changes in vessel compliance also contribute to the observed increase in PWV.
3.4 Spatial and Temporal Compliance.
Spatial and temporal changes in vessel compliance are next investigated using the four blood pressure waveforms: ((i) uniform brachial pressure wave, (ii) spatially varying pressure wave computed using the unsteady Bernoulli approach, (iii) spatially varying pressure wave determined using the piecewise approach, and (iv) spatially varying pressure wave computed by solving a full patient specific CFD analysis) determined in Sec. 2.5. In Fig. 5, the instantaneous lumen cross-sectional area is plotted as a function of blood pressure. Results are presented for the four aforementioned pressure waveforms at the proximal ascending aorta and the distal abdominal aorta. The instantaneous compliance at a given lumen pressure is given by the slope of the pressure–area graph. In all cases, two distinct linear regions are observed such that a high vessel compliance occurs at low pressures, and a low vessel compliance occurs at high pressures. The decrease in compliance at higher lumen pressures is due to strain stiffening constitutive behavior of the aortic wall.
For each case presented in Fig. 5, the value of compliance is determined using linear regression fits for the two distinct regions of the pressure–area graphs (values are presented in Table 1). Clear evidence of strain stiffening is visible in each subplot of Fig. 5, where significantly higher aortic dilation for a given change in pressure are observed in the high compliance regime, compared to the low compliance regime. As an example, for Plane 2 (Fig. 5(a) (uniform blood pressure waveform)), the compliance at low pressure (C21 = 11.94 mm2/mmHg) is over eight times higher than the compliance at high pressure (C22 = 1.48 mm2/mmHg). For plane 2, the high compliance regime occurs for pressures below 85 mmHg. While a broadly similar bilinear behavior is also observed at plane 10 (abdominal aorta), compliance values are an order of magnitude lower than those at plane 2 (ascending aorta). As an example, for plane 10 (Fig. 5(e) (uniform blood pressure waveform)), a high compliance value of C101 = 0.67 mm2/mmHg is determined, with a low compliance value of C101 = 0.11 mm2/mmHg. Furthermore, at plane 10, the change in compliance regime is observed to occur at a pressure of –100 mmHg (compared to –85 mmHg at plane 2). These results highlight the dramatic differences in in vivo material behavior between the thoracic and abdominal aorta. Despite the fact that higher material strains occur in the thoracic aorta, as evident from Fig. 8, the instantaneous material stiffness is significantly higher in the abdominal aorta. The higher stiffness of the abdominal aorta explains, in part, the higher PWV in this region, as observed in Fig. 10.
Figure 11 shows an approximation of the circumferential stress versus strain, derived from the area-pressure curves for the proximal ascending and distal abdominal aorta, based on the law of Laplace where P is pressure, r is the radius, and t is wall thickness. In each case, significant strain stiffening is observed. It should be noted that these computed values do not consider the unloaded or stress-free reference configuration and are purely to demonstrate the significant stiffening observed throughout the cardiac cycle due to the straightening of collagen fibers. Follow-on work has been conducted by the authors incorporating the zero-pressure equilibrium configuration with a novel physiologically motivated constitutive law to capture the nonlinearity in a physical manner.
In this study, a dual-VENC 4D flow MRI protocol is developed to achieve accurate measurement of the dynamically changing flow velocity field and lumen area throughout the entire cardiac cycle and the entire aorta. To the best of our knowledge, no previous medical imaging paper has reported such detailed spatial and temporal characterization of the human aorta. A nonlinear relationship between lumen area and pressure is observed in vivo over the duration of a cardiac cycle throughout the entire aorta, suggesting that aortic biomechanics may not be accurately characterized by a single value compliance coefficient, as commonly assumed [32–34]. Furthermore, our detailed in vivo measurements reveal that the lumen pressure–area relationship and PWV are highly heterogeneous throughout the aorta.
4.1 Spatial Deformation.
We examine the deformation of the human aorta during the entire cardiac cycle at 10 planes, ranging from the sinus of Valsalva to immediately proximal to the common iliac bifurcation. The high levels of cross-sectional-area change, during a cardiac cycle, ranging from 15% in the abdominal aorta to 65% in the ascending aorta, highlight the extremely large deformations of the aortic wall. The levels of deformation observed over a cardiac cycle are similar to those reported in previous studies. Sonesson et al. report a 20% increase in diameter in the abdominal aorta in young adults using ultrasound  (approximately 44% area increase), while [68–70] all report an increase in area of >100% in thoracic mice aortae. Accurate characterization of such large deformations requires detailed imaging of the entire aorta throughout the entire cardiac cycle. While the observed trend that dilation reduces with distance from the heart is in broad agreement with previous noninvasive imaging studies by Refs.  and , the current study provides further insights by measuring dilation on a large number of planes spanning the entire aorta. A number of ex vivo studies also suggest that compliance decreases with distance from the heart [12,13]. A study by Tsamis et al.  reports that the ascending thoracic aorta contains 80 elastin lamellar units while the infrarenal abdominal aorta contains 32. The decrease in elastin and increase in collagen observed by Concannon et al.  provide a microstructural explanation for the decrease in compliance observed here with distance from the heart. Moreover, Tsamis et al.  also report a 50% decrease in elastin units between the descending thoracic and supraceliac aorta, possibly providing an explanation for locally varying cross-sectional-area change observed in the current study. A review paper by Sherif  reports that the aorta, from a developmental point of view, is neither a homogeneous structure nor one contiguous anatomical entity. Rather, it is suggested that the vessel can be split into discrete segments, each of which develops and differentiates under a distinct set of genetic and transcriptional factors. It is hypothesized that the regional differences in biomechanical behavior may be due to the development of the ascending thoracic from neural crest cells and descending thoracic aorta from the mesoderm. With distinct connections or “weld points” between such segments, this may be the cause for local differences in cross-sectional-area change and PWV measured between adjacent planes in the current study.
4.2 Spatial Hemodynamics.
A notable outcome of this study is the spatial variance in PWV along the aorta. Results show that the PWV increases with increasing distance from the heart. This finding is reinforced using the TTP method, uncovering a 28% increase in PWV between the thoracic and abdominal aorta (6.4 m/s versus 8.2 m/s). Generally, PWV is defined in the literature as a single value for the aorta [9,75–78]. The assumption of a uniform single valued PWV is primarily due to the method of clinical measurement, where the pressure pulse between two distinct sites, most commonly the carotid and femoral arteries (cfPWV) is recorded. The reference values for arterial stiffness' collaboration , report a mean PWV value of 6.2 m/s for a cohort of 1455 normal subjects < 30 year of age. However, the pathway over which cfPWV is defined does not include the highly compliant ascending aorta. The utilization of MRI techniques to quantify aortic PWV has the ability to quantify changes at a local level, producing an accurate patient-specific spatial map of PWV. A study by Quinaglia et al.  reported PWV readings targeted to the ascending aorta and found velocities of between 4 and 5.8 m/s, while  investigated the brachio-femoral pathway in 152 young adults and found mean PWV values of 8.7 m/s. Such measurements are comparable with the data presented in the current study for the thoracic and abdominal aorta, respectively.
4.3 Spatial and Temporal Compliance.
Compliance is generally presented as a single value, by taking the difference in area between diastole and systole and dividing this by patient's change in blood pressure. Aortic tissue is not a simple linear elastic material. Rather, it exhibits a significant increase in stiffness when it is stretched to a high level of deformation . Such mechanical behavior occurs due to the structural contribution of collagen fibers. At low arterial strains, collagen fibers are wavy, and an incremental increase in applied force will result in a significant increase in the length of the fiber, i.e., the fiber exhibits a low structural stiffness at low levels of deformation. A further incremental increase in force applied to a straightened collagen fiber will not result in a large increase in the length of the fiber. This is because the straightening of the fiber at high levels of deformation results in an increase of the structural stiffness .
The structural contribution of collagen results in the well-established nonlinear stress–strain relationship for arterial tissue, whereby the material exhibits low stiffness at low strains and high stiffness at high strains. The transition from the low stiffness regime to the high stiffness regime is commonly modeled using exponential strain stiffening material laws [83,84]. To date, such models have been motivated and calibrated using in vitro tests of excised arterial tissue. Our study provides evidence that significant strain stiffening of the aorta occurs in vivo over the deformation range of a cardiac cycle. This suggests that clinical compliance (defined as a change in lumen area with respect to a change in pressure) should not be characterized by a single value. Rather, a high compliance regime is observed for low pressures during diastole, followed by a transition to a low compliance regime for high pressures during systole. This in vivo observation is consistent with strain stiffening observed in in vitro testing, and it calls into question the accuracy of the common assumption that in vivo lumen area increases linearly with lumen pressure during a cardiac cycle (inherent in the description of compliance by a single coefficient, e.g., Refs. [33,85–90]).
Previous in vivo analyses of the human aorta include an investigation of compliance in the abdominal segment using ultrasound , where the authors report a decrease in compliance with age; however, such an imaging modality is impractical in portions of the thoracic aorta due to blind spots from bronchial air . Mohiaddin and colleagues investigated aortic compliance using MRI in the thoracic segment in a large cohort of 70 volunteers , and found that compliance was highest in the ascending segment; however, images were only acquired at diastole and systole, a temporal resolution too low to capture any nonlinearities.
Our study quantifies the values of high and low compliance during a cardiac cycle, and demonstrates that these values, and the associated transition pressures, are spatially heterogeneous. Results show that the aorta exhibits a high compliance regime at low pressures and a low compliance regime (LCR) at higher pressures, within each cardiac cycle. This nonlinearity in compliance has been widely observed in vitro whereby the instantaneous stiffness of arterial tissue increases with increasing uniaxial or biaxial stretch [92–98]. It is important to note when attempting to characterize aortic mechanical properties noninvasively, that the diastolic configuration extracted from in vivo analyses is neither the zero-pressure nor the stress-free configuration. Therefore, the fitting of in vivo pressure–area data without consideration of the subphysiological regime, will yield unphysical results. Follow-on work has been conducted by the authors, identifying an equilibrium vessel configuration at zero applied lumen pressure, which is observed to be critical step required in order to predict the key features of the pressure–area relationship observed in vivo. The role of elastin prestretch on the lumen pressure at which the aorta transitions from a high compliance to a low compliance regime due to collagen strain stiffening, is also investigated using a novel physically based constitutive law. This modeling approach is also shown to capture the key features of elastin and smooth muscle cell knockout experiments. Such detailed insights into vessel compliance are critical for development of an enhanced understanding of the relationship between pressure, blood flow, and PWV in the aorta, and will potentially lead to improved interventional procedures and device designs.
A number of limitations should be noted for the current study, providing motivation for follow-on studies. The purpose of this study was to develop a dual-VENC imaging protocol to generate high resolution subject-specific data on heterogeneous nonlinear aortic compliance and pulse wave velocity in a clinically feasible timeframe. While the data generated in the current study are limited to a single subject, the demonstration of this capability of our methodology provides a platform for extensive high-resolution characterization of aortic biomechanics for populations of healthy and diseased subjects. It should be noted that increased temporal resolution, spatial resolution, coverage, and signal-to-noise ratio all incur the cost of higher scan time and gradient coil capabilities in every MRI system. Hence, in order to maintain clinical feasibility, temporal resolution was sacrificed in this study. In the ideal situation, each phase would span a segment shorter than 50 ms, which may lead to greater accuracy in the quantification of area, flow and hence PWV, and so, more work may be justified in this area to see if any further optimization of parameters is possible for imaging the aorta in its entirety, while maintaining a short scan time. Further improvements in spatial resolution may be possible with the clinical integration of 7 T scanners that may aid in the quantification of compliance in older stiffer aortae. It is expected that the compliance estimates obtained from the young healthy case would be higher than those of an older subject, due to the natural process of arteriosclerosis that occurs with age. Moreover, for older/sicker patients, faster more irregular cardiac cycles will increase scan time, which presents challenges in whether the MRI machine is capable of capturing, for example, 20 phases in a shortened irregular RR-interval. It should also be noted that this protocol can be readily applied to a 1.5 T scanner, with the drawback of a significant time increase.
The same level of accuracy in area quantification cannot be achieved along the length of the aorta due to varying levels of compliance and a fixed spatial resolution, where deformations smaller than the pixel size are not registered. A possible solution to this exists in running a series of 2D PCMRI scans, with increasing in-plane resolution further from the heart; however, such scans do not account for flow continuity between planes and extreme care should be taken to prescribe imaging planes orthogonal to the mean direction of flow at each location of interest.
In terms of determining aortic compliance, a challenge remains to accurately measure a continuous location-specific blood pressure waveform throughout the aorta without resorting to an invasive catheterization procedure. In the absence of a clearly defined best strategy to compute a continuous pressure waveform noninvasively [50,99], we implemented four separate waveform generation methods, namely, “uniform,” “unsteady Bernoulli,” “piecewise,” and “CFD.” In any case, this study demonstrates that the bi-linearity of the measured compliance is not strongly affected by the method of approximating the lumen pressure waveform. Although there are differences in the slopes between each method, the bilinear nature of the pressure–area relationship is apparent in each case, and we rely on the goodness of fit to the data to provide the compliance estimates.
The results of this study have a number of potential implications for the fields of aortic biomechanics and cardiovascular surgery. The study presents a protocol that can provide accurate spatial and temporal measurements of compliance and PWV in the aorta. This may provide an incremental step in understanding why cardiac events occur post-TEVAR, through a better understanding of the relationship between PWV, aortic stiffness, and cardiac function. Stenting may have a spatially varying effect on the biomechanics of the aorta by inducing a cascade analogous to “accelerated arteriosclerosis” on the system. This in turn effects cardiac function, as documented elsewhere for arteriosclerosis developed during the ageing process [72,100,101]. During EVAR, however, a significant reduction in compliance may occur instantaneously due to stent deployment, in contrast to arteriosclerosis, where compliance gradually reduces over a period of decades. A follow-on study by the authors demonstrates the importance of accurate characterization of nonlinear aortic compliance and its implications on Nitinol stent–artery interactions. Simulations reveal that Nitinol stent grafts stretch the artery wall so that collagen is stretched to a straightened high stiffness configuration. The high compliance regime associated with low diastolic lumen pressure is eliminated, and the artery operates in the LCR throughout the entire cardiac cycle. The slope of the lumen pressure–area curve for the LCR postimplantation is almost identical to that of the native vessel during systole. This negligible change from the native systole slope occurs because the stent graft increases its diameter from the crimped configuration during deployment so that it reaches a low stiffness unloading plateau (The effective radial stiffness of which is negligible compared to the stiffness of the artery wall). This highlights the need for accurate quantification of nonlinear compliance in order to provide a mechanistic foundation for the common assumption that stents decrease aortic compliance [102–105]. The current study suggests that aortic compliance cannot be captured by a single value, and that the vessel is significantly less compliant in systole than diastole. Incorporating such detailed information into the design of EVAR devices with the aim of replicating the natural nonlinear compliance of the vessel may reduce the prevalence of the aforementioned complications.
A dual-VENC 4D Flow MRI protocol is developed and implemented in a commercial scanner for characterizing the biomechanics of the entire human aorta. A composite dataset approach is employed to maximally attenuate fluid contrast throughout the unsteady velocity profile of the cardiac cycle, providing an alternative method to phase unwrapping techniques. Pulse wave velocity increases from proximal to distal aorta, while cross-sectional-area change, volumetric flow rate, and compliance all reduce with distance from the heart. Finally, compliance is shown to alter significantly during the cardiac cycle, with significantly higher compliance being observed during periods of low blood pressure.
The author would like to acknowledge Mr. Torben Shneider and Mr. Matthew Clemence of Philips (Philips Medical Systems, Best, the Netherlands) and Evelyn Smith and Marie McMullen (Department of Radiology, Galway Clinic) for their support and guidance during this study.
Irish Research Council for Science, Engineering and Technology (EPSPG_2016_194; Funder ID: 10.13039/501100001596).