Haemodynamics in the mouse aortic arch computed from MRI-derived velocities at the aortic root

Mice are widely used to investigate atherogenesis, which is known to be influenced by stresses related to blood flow. However, numerical characterization of the haemodynamic environment in the commonly studied aortic arch has hitherto been based on idealizations of inflow into the aorta. Our purpose in this work was to numerically characterize the haemodynamic environment in the mouse aortic arch using measured inflow velocities, and to relate the resulting shear stress patterns to known locations of high- and low-lesion prevalence. Blood flow velocities were measured in the aortic root of C57/BL6 mice using phase-contrast MRI. Arterial geometries were obtained by micro-CT of corrosion casts. These data were used to compute blood flow and wall shear stress (WSS) patterns in the arch. WSS profiles computed using realistic and idealized aortic root velocities differed significantly. An unexpected finding was that average WSS in the high-lesion-probability region on the inner wall was actually higher than the WSS in the low-probability region on the outer wall. Future studies of mouse aortic arch haemodynamics should avoid the use of idealized inflow velocity profiles. Lesion formation does not seem to uniquely associate with low or oscillating WSS in this segment, suggesting that other factors may also play a role in lesion localization.


INTRODUCTION
Mouse models of chronic diseases such as atherosclerosis are valuable tools for understanding human pathology. Indeed, the use of genetically modified strains has been a major breakthrough in investigations of the biological pathways involved in atherogenesis [1,2]. Lesions are readily induced in the larger arteries of susceptible strains, and these vessels are also accessible for imaging; the aortic arch has received particular attention [1 -4].
It is well known that atherogenesis is influenced by haemodynamic factors, the most prominent of which is wall shear stress (WSS). Several studies [5 -10] have examined WSS patterns in the mouse aortic arch by numerically simulating blood flow, for which both the accurate arch geometry and inflow boundary conditions (blood velocity profiles at the aortic root) must be known. Previous studies have used idealized (e.g. parabolic or Wormersley) velocity profiles at the root. Using a more physiologically realistic profile could have a significant impact on computed WSS patterns. More specifically, flow in a curved vessel such as the aortic arch typically results in a velocity profile that is skewed towards the outer curvature, reducing WSS on the inner curvature [11], although this effect is reversed near the inlet of a curved tube under certain circumstances [12]; a skewed inlet velocity profile could intensify or weaken this effect. The use of physiologically accurate profiles is therefore important to our understanding of the relationship between atherogenesis and biomechanical factors.
Our goals in this work were to measure temporally and spatially resolved mouse-specific blood velocity profiles at the aortic root and to determine their effect on aortic arch haemodynamics. We used MRI to measure velocities, corrosion casts to determine aortic arch geometries and numerical modelling to compute WSSs. We found that using realistic aortic root velocity profiles had an impact on computed aortic arch haemodynamics, and that our blood flow simulations predicted a surprising similarity between the haemodynamic environment in regions known to differ in their predilection for the formation of atherosclerotic-like lesions.

Animals
Animals were housed in a temperature-(21 + 28C), humidity-(50 + 10%), and light-(non-inverted 12 L : 12 D cycle) controlled environment and given ad libitum access to food and water. Ten adult C57BL6 male mice (Harlan, UK) were imaged, of which five animals had datasets of sufficient quality for subsequent use (table 1).
Anaesthesia was induced with 4 per cent isoflurane in 100 per cent oxygen and maintained with approximately 1.5 per cent isoflurane in 1.5 l min -1 O 2 via a facemask, while the mouse lay in the supine position on a platform. Subdermal ECG needle electrodes and a respiratory cuff (SA Instruments, Stony Brook, NY, USA) were placed for monitoring heart and respiratory rate, respectively. Body temperature was monitored via a rectal probe and maintained at 378C (36.7 + 0.28C) using a heated-air system (SA Instruments, Stony Brook, NY, USA). For the remainder of this work, all mice will be referenced according to the date of the MRI scan (e.g. 210909AM).

Magnetic resonance imaging
The MRI slice used for velocimetry was centred on the pulmonary artery (figure 1). Details of imaging protocols are given in the electronic supplementary material, supplementary methods. For one-dimensional blood flow velocity encoding, two sets of gradient echo fast low angle shot (FLASH) images were acquired with all imaging parameters kept constant, except for the first moment of the flow encoding gradient. Flow was measured separately along three orthogonal directions: through-plane, and in-plane along phase and frequency encoding directions. Based on preliminary images, the waveforms of the through-plane flow encoding gradient were calculated to fulfil the requirement that a peak velocity (VENC) of 140 cm s 21 corresponded to a phase shift of p. VENC for the smaller in-plane (secondary) velocities were iteratively decreased to 40 cm s 21 , except for mouse 070909AM which was scanned with an in-plane VENC of 50 cm s 21  Data acquisition was performed with dual cardiac and respiratory gating while preserving steady-state condition of longitudinal magnetization. The respiratory rate during acquisition was 36 + 10 breaths per minute (mean + s.d.) and the average heart rate was 520 + 45 bpm for the five mice used in this study. Measurements of blood flow were taken every 6 ms during the cardiac cycle at time points from 0 to 48 ms after the onset of systole, as defined by the R wave and a 4 ms trigger delay.

Processing of the MRI velocity data
Several post-processing steps were required (figure 2). First, a mask was created that identified voxels in the lumen of the ascending aorta. It was determined manually using a combination of the anatomical and velocity data. The size and shape of the mask changed over the cardiac cycle as the aorta moved. Second, the in-plane velocity fields contained vectors that differed substantially from their neighbours and were thus likely erroneous. Such spurious vectors are also a feature of particle image velocimetry datasets, and methods have been developed to detect and correct them. We used the dynamic mean operator [13] to detect spurious vectors and an inverse distance interpolation of nearest, non-spurious neighbours to correct them.
Blood velocity was measured in the magnet's reference frame but the subsequent blood flow modelling and intermouse comparison required velocities expressed relative to the moving vessel wall. In order to compensate for vessel movement, vessel wall velocity was calculated using the no-slip condition (fluid velocity next to the wall must equal the wall velocity). 'Border nodes' that were inside the artery but had a neighbour outside the artery were identified from the ascending aorta mask, and the average velocity over all border nodes was taken as vessel wall velocity at each time step. To compensate for vessel movement, this velocity was subtracted from the measured blood velocity.
To allow inter-mouse comparisons, datasets were mapped to a common circular domain and rotated so that the dorsoventral axis of each mouse was in a common orientation. A consistent circular domain was created by warping the aortic root to fit to a circular domain without changing any velocity values. The dorsoventral axis was defined by determining the geometric centres of the spine and sternum from the MRI magnitude images. Further details of post processing are given in the electronic supplementary material, supplementary methods.

Micro-CT post processing
Aortic geometries were obtained by micro-CT scanning of corrosion casts obtained post-mortem. Scans were segmented to isolate the geometry of interest and create a three-dimensional surface model that could be used as the basis for a computational mesh for blood flow modelling (see the electronic supplementary material, supplemental methods). The cast from 210909PM showed artefacts and was not used in computational blood flow modelling, although the MR blood velocity data were unaffected and were included when computing average velocities.

Blood flow modelling
Blood flow and WSS patterns were computed numerically as described in the electronic supplementary material, supplementary methods. On the basis of the four mice with good casts and MR velocimetry data from five mice, seven flow simulations were completed as follows. One simulation was performed to compute unsteady aortic blood flow in each of the four mice with good casts, using mouse-specific inlet velocities determined by MRI and the corresponding cast geometry. Three additional unsteady simulations were performed using the cast geometry of mouse 210909AM. In the first, only the through-plane MRI velocities were specified at the aortic root, the inplane components being set to zero. This simulation explored the importance of the in-plane velocity components on the WSS field. The second used a uniform aortic root velocity profile in the through-plane direction and no in-plane velocity components, with the flow rate matched to the full MRI measurements. This simulation was compared with the previous one to explore the effects of using idealized rather than realistic through-plane inlet velocities. In the third, mouseaveraged, unsteady MRI-determined velocities were applied in order to determine how a realistic but nonsubject-specific set of velocity conditions compares to subject-specific aortic conditions.

RESULTS
A plot of flow rate as a function of time for the five mice with good MRI scans is shown in figure 3. There was good agreement among mice with respect to the peak flow rate and the general shape of the flow waveform. Mouse 210909PM had notably more prolonged flow than the other mice, and also the lowest heart rate (table 1). Aortic root blood velocity data for these five mice showed complex flows with some common features. In order to better understand blood flow features at the aortic root that were common among mice, we computed a mouse-averaged velocity field and the standard error of the mean for this velocity field (figures 4 and 5). (Individual data are given in the electronic supplementary material, supplementary results.) The velocity patterns showed several differences to analogous measurements in humans [14]. The through-plane standard error (contours in figure 5) was relatively low for the first six time points (up until t ¼ 30 ms); for t . 30 ms, the error was of the same order as the velocity itself. The standard error in the in-plane directions (vectors in figure 5) was of the same order as the in-plane velocity at all time points.
At t ¼ 0 ms, the mouse-averaged velocity field exhibited only small velocities. At t ¼ 6 ms, there was through-plane acceleration of the blood while the inplane flow was noticeably disorganized, with few structures being evident. At t ¼ 12 ms, close to peak systole, the through-plane velocities reached a maximum near the centre of the vessel with the highest-speed fluid zone extending slightly in the dorsal and ventral directions. The in-plane flow displayed more organized behaviour, although distinct structures were difficult to discern. At 18 ms there was a deceleration of the through-plane flow compared with the 12 ms pattern, and the highest-speed fluid zone moved towards one side of the vessel. This behaviour was maintained in the through-plane flow pattern at 24 ms. At both 18 and 24 ms, the in-plane data displayed clockwise rotating vortices that slowly dissipated over the duration of the cardiac cycle. At 30 and 36 ms, the deceleration continued, as did the movement of the highest-speed fluid zone towards the right-ventral direction. Additionally, the flow on the left side of the vessel was reduced greatly and even reversed slightly. At 42 -48 ms, there was little through-plane flow and the in-plane flow was reduced to weak clockwise vortices.   Figure 6 also demonstrates that using realistic but not subject-specific velocity conditions (mouse averaged, figure 6d) resulted in WSS patterns that were very similar to those obtained by using subject-specific conditions.
Oscillatory shear index (OSI) is a measure of how much the shear stress reverses direction over the cardiac cycle. It was defined by He & Ku [15] as where t is the WSS vector and T is the duration of the cardiac cycle. OSI values on the inner and outer wall were higher than elsewhere (including the lateral walls), with a marginally higher value over a larger area on the inner wall than the outer (figure 6e). Thus, both the inner and outer walls of the arch before the first branch were subjected to low levels of WSS and broadly similar levels of OSI. For reference, the locations of high (HP) and low (LP) probability of lesion formation, as defined by Iiyama et al. [4], are shown in figure 6f. The location of markers and the corresponding regions of high-and low-probability of lesion development are approximate, given inter-mouse variation in aortic arch geometry and the schematic nature of the diagram used to locate the HP and LP zones.  Shear stresses were further quantified within 'patches' on the luminal surface of the inner and outer walls (figure 7) that were chosen to coincide with the HP and LP areas identified by Iiyama et al. [4]. Both patches were made slightly larger than these areas in order to ensure overlap with them, because it was difficult to precisely map the areas onto our geometry.
The inner wall patch was approximately 50 per cent larger than the HP area of Iiyama et al. and the outer wall patch was approximately 25 per cent larger than the LP area.
A histogram of the WSS levels occurring within each patch (figure 7) confirmed that WSSs on the inner wall increased as inlet boundary conditions were changed  from full MRI-measured velocities to through-plane only MRI measurements to the flat inlet profile, with the opposite holding true on the outer wall. Figure 7 also confirms that the WSS on the outer patch was lower in magnitude than the WSS on the inner patch. This was quantified by calculating average values: shear stress was 0.882 (53 dynes cm 22 ; range 0.68 -1.46) for the inner wall patch and 0.664 (40 dynes cm 22 ; range 0.36 -1.00) for the outer wall patch, using the full MRI inlet conditions. Finally, figure 7 also confirms quantitatively that using the mouse-averaged inflow velocity data gives results that closely followed the full MRI data. Similar histograms of OSI ( figure 7) show that the dependence of OSI on inlet conditions is greater than the dependence of WSS on these conditions. We also observe that there is lower OSI in the outer curvature patch compared with the inner curvature patch, an observation that seems to support the idea that high OSI (i.e. low and oscillating WSS) promotes lesions. However, further examination from a slightly different orientation (figure 6e) shows that the relatively low OSI patch seen in figure 7d is immediately adjacent to a region of high OSI on the outer curvature, which is not a site of lesion formation. This emphasizes the fact that our data do not show a consistent association between low and/or oscillatory WSS and sites known to have a predilection for lesion formation.

DISCUSSION
The region of the mouse ascending aorta with the highest probability of developing lipid deposition is the inner wall (lesser curvature) of the distal segment, before the first branch, while the area with the lowest probability is the outer wall (greater curvature) of the same segment [4] (figure 6f ). This variation has been used to investigate the possible dependence of lesion prevalence on mechanical stress. Endothelial cells within the high probability (HP) region are less elongated and more randomly aligned than those in the low probability (LP) region, consistent with an expected difference in shear stresses [16]. Previous studies [5 -10] have systematically mapped the pattern of WSS in these regions. However, these studies have assumed idealized aortic root blood velocity profiles. This work investigated the effect of using more realistic inflow velocity conditions on computed shear stress patterns. Throughplane and in-plane velocities were obtained from phase contrast MR images and high-resolution anatomical information was obtained from micro-CT of    1 1.2 1.3 1.4 1  vascular corrosion casts. This approach is subject to certain limitations, which should not affect the overall conclusions [17]. WSS patterns were generally insensitive to inflow conditions at locations distal to the first branch with the exception of localized changes around the second branch ostium. However, in regions between the inlet and the first branch, the choice of inlet boundary conditions affected the level of WSS on both the inner and outer wall. When using a flat inlet velocity profile, the WSS was under-estimated on the outer wall and overestimated on the inner wall, compared with using the full MRI inlet velocity boundary conditions. The absence of in-plane velocity components contributed to this effect; excluding them again caused an under-estimation of WSS on the outer wall and over-estimation on the inner wall. In fact, the presence of swirling (in plane) velocities at the aortic root appears to reduce spatial and temporal variations in WSS, consistent with previous studies [18,19] where swirl induced by out-of-plane curvature reduced spatial variations in WSS.
These results demonstrate the importance of using realistic inflow data when computing shear stresses in the mouse aortic arch. Such data are not easy to obtain. It may therefore also be important that the non-subject specific (average mouse) inlet velocity boundary conditions contained coherent throughplane and in-plane flow structures and resulted in a very similar WSS pattern to those computed using the full (mouse-specific) boundary conditions. This may justify the use of our mouse-average inlet velocity field (see the electronic supplementary material) as a reasonable approximation of the actual velocity field in future studies. However, it should be borne in mind that velocities are likely to depend on factors such as age, weight, sex and strain as well as the anaesthetic used.
As in previous studies [5 -10], low WSS was found on the inner wall (normalized shear stresses less than one), coinciding with the HP region. The presence of this patch of WSS was, as expected, attributable to the development of counter-rotating vortices, much like the Dean vortices forming in curved tubes under steady flow conditions. Another feature of the WSS patterns was the presence of an area of low values on the outer wall of the arch, between the aortic root and the first branch. By analysing the average WSS in patches on the inner and outer walls, we unexpectedly found that the timeaveraged WSS in the LP region was actually lower than in the HP region ( figures 6 and 7).
The low WSS zone on the outer wall resulted from a recirculation zone on the upstream wall of the first branch. Recirculation occurred at this location (and at the equivalent sites in the other branches) at peak systole because the flow failed to change direction to enter the branch and became detached from the wall (see the electronic supplementary material, supplementary results). The recirculation zones were themselves characterized by relatively slow-moving fluid near the wall and hence low WSS. Additionally, the recirculation zone provided a barrier to flow entering the first branch along the outer wall of the arch. The fluid decelerated and was forced outwards in the dorsal and ventral directions away from the central plane of the arch. The deceleration and secondary flow also caused low WSS.
This feature of the WSS pattern was generally consistent across different mice when using subjectspecific geometry and inlet velocity data. However, individual variation in geometry could disrupt it: the patch of low shear was not seen in one mouse (070909PM) that had tighter aortic arch curvature than the other mice, an arch that started to curve closer to the aortic root and a more distal position of the first branch. These features caused earlier and stronger formation of the Dean-type secondary flow, skewing the highest velocity toward the outer wall and hence causing higher levels of WSS there.
Other authors have not explicitly remarked on the presence of a similar low WSS region on the outer curvature, although because this region has not previously been the focus of attention, many existing studies do not show clear images of the WSS distribution in this area. However, the study of Trachet et al. [20] does demonstrate a remarkably similar low WSS region (see their fig. 2 for young adult mice), and is entirely consistent with our findings. In reviewing studies that do show cycle-averaged WSS distributions on the outer curvature, it became clear that the existence of a low shear region on the outer curvature can be suppressed by sharp arch curvature and a more distally located brachiocephalic ostium, which is consistent with our findings among our cohort of mice and may explain why this region was not observed and/or commented on by others. Another factor affecting this low WSS zone may be the flow split at the brachiocephalic trunk. In a recent study, Trachet et al. [21] carefully measured volumetric flow rate from the major branches of the aorta, observing more flow entering the brachiocephalic than we specified (22.5% versus our value of 15.3%). Using this higher value would in fact make the separation effect that we observe more pronounced, suggesting that our assumed flow splits may have underestimated the extent of the low WSS zone in our mice.
The coincidence of low WSS and high lesion probability on the inner wall is consistent with the low WSS theory of atherogenesis [22]. On the other hand, the coincidence of low WSS and low lesion probability suggests that factors in addition to low and oscillatory WSS could be important in atherogenesis. The apparent importance of other localizing factors in atherogenesis could be explained by strain-dependent differences in aortic arch geometries. Different strains can have different aortic arch geometries [8], and arch geometry can affect WSS patterns on the outer wall as evidenced by mouse 070909PM. The LDLR 2/2 mice used by Iiyama et al. when mapping lesion prevalence were of mixed 129-C57BL/6 background while the strain for our MRI velocity measurements was C57BL/6; Zhu et al. [8] found differences in arch geometry between C57BL/6 and 129/SvEv mice. We must also keep in mind the possibility that the haemodynamic environment could differ between the wild-type (WT) mice used in this study and the knockout (KO) mice used in lesion mapping studies. Potential sources for such differences include arch geometry, the magnitude and waveform of cardiac output, blood viscosity and arterial stiffness. These issues have been studied by others. Zhu et al. [8] did not observe large differences in arch geometry between WT and KO mice. One study reported a higher cardiac output in older KO mice versus WT [23], while another [24] reported the opposite for younger mice, and flow rates in the KO mice varied by a factor of two; despite this wide range, there is a clear predilection for lesions to form on the inner, as opposed to the outer, curvature of the aortic arch in Apo-E 2/2 and LDLR 2/2 mice [25], indicating that disease localization in KO mice is not particularly sensitive to cardiac output. Finally, KO mice were reported to have an 11 per cent lower haematocrit and a 13 per cent higher pulse-wave velocity versus WT mice [23], but these are expected to have only modest effects on arch haemodynamics, where flow patterns are dominated by curvature, branch outflow rates and inlet velocity profile. We therefore believe that our haemodynamic results are valid for KO mice, but further confirmatory haemodynamic studies in KO mice are warranted.
Another possibility, recently suggested by Zhou et al. [26], is that retrograde flow rather than low shear stress is the localizing factor leading to lesion development. This suggestion was based on Doppler ultrasound measurements that showed late systolic/early diastolic retrograde flow on the inner but not the outer curvature of the ascending aorta, and by the pro-atherogenic effect of experimentally increasing retrograde flow. It is consistent with the earlier numerical results of Suo et al. [6], which indicated a greater change in flow direction during the cardiac cycle at the inner wall than the outer wall. Our own results (figure 6) also show greater changes in direction (characterized by the OSI) on the inner than those on the outer wall. However, both these regions had much higher OSI values than other parts of the ascending aorta. Hence, the lowest OSI regions were not those associated with the lowest lesion prevalence, and areas associated with the lowest lesion prevalence had moderately high values of OSI. The apparent discrepancy with the results of Suo et al. and Zhuo et al. may derive from the fact that reverse flow was characterized at only a few points in the ascending aorta in these earlier studies; our maps show significant variation of OSI over short distances ( figure 6). The discrepancy may also reflect the limited resolution of the Doppler measurements, the use of a flat inlet profile in the earlier numerical study, or other technical differences such as the compensation for effects of anaesthesia in the present study (see the electronic supplementary material, supplementary methods).
Alternatively, low WSS and oscillatory flow may be important, but the measures of WSS level and oscillations that we used, namely cycle-averaged WSS and OSI, respectively, may not be the indicators most closely associated with disease propensity, even though they are widely used in the literature. For example, atherogenesis might depend on the amount of time in the cardiac cycle that endothelial cells are exposed to WSS below a certain threshold rather than on the cycle-averaged WSS level, or on fluctuations in flow direction other than the complete changes in sign captured by the OSI.
Of course, a further possibility is that atherogenesis is not exclusively due to low and/or oscillating shear stress, however parametrized. For example, atherogenesis may depend on the mass transfer of blood-borne molecules (e.g. LDL and NO) between blood and endothelium, which is usually, but not inevitably, related to WSS [17,27]. In rabbit models, there is evidence for a dependence of lesions on wall strain [28] and, as here, for association with variables in addition to low shear [29]. It is relevant in this context that wall strain is likely to vary between the inner and outer walls of the aortic arch, and that it can influence endothelial cell elongation [30].
Scaling laws dictate that many aspects of mouse haemodynamics differ from those in people [31]. Indeed, the mechanisms driving atherogenesis may differ between mouse and human arteries. Nonetheless, we suggest that further study of the role in atherogenesis of mechanobiological factors in addition to low/oscillatory WSS is warranted.