Computation in the rabbit aorta of a new metric - the transverse wall shear stress - to quantify the multidirectional character of disturbed blood flow.

Spatial variation of the haemodynamic stresses acting on the arterial wall is commonly assumed to explain the focal development of atherosclerosis. Disturbed flow in particular is thought to play a key role. However, widely-used metrics developed to quantify its extent are unable to distinguish between uniaxial and multidirectional flows. We analysed pulsatile flow fields obtained in idealised and anatomically-realistic arterial geometries using computational fluid dynamics techniques, and in particular investigated the multidirectionality of the flow fields, capturing this aspect of near-wall blood flow with a new metric - the transverse wall shear stress (transWSS) - calculated as the time-average of wall shear stress components perpendicular to the mean flow direction. In the idealised branching geometry, multidirectional flow was observed downstream of the branch ostium, a region of flow stagnation, and to the sides of the ostium. The distribution of the transWSS was different from the pattern of traditional haemodynamic metrics and more dependent on the velocity waveform imposed at the branch outlet. In rabbit aortas, transWSS patterns were again different from patterns of traditional metrics. The near-branch pattern varied between intercostal ostia, as is the case for lesion distribution; for some branches there were striking resemblances to the age-dependent patterns of disease seen in rabbit and human aortas. The new metric may lead to improved understanding of atherogenesis.

Introduction

The patchy distribution of atherosclerotic lesions within the arterial system appears to be a consequence of local variation in blood flow characteristics (Giddens et al., 1993; Steinman, 2004). Disturbed flow may play a key role; it has most commonly been identified using a low mean wall shear stress (WSS) criterion (Caro et al., 1971; Cecchi et al., 2011) or from the degree of flow reversal, quantified by the oscillatory shear index (OSI) (Ku et al., 1985). A variety of other parameters have been proposed, including the relative residence time (RRT) (Himburg et al., 2004), the WSS spatial gradient (Lei et al., 1995), the WSS angle gradient (Longest and Kleinstreuer, 2000), the WSS angle deviation (Hyun et al., 2000), the peak WSS temporal gradient (Ohja, 1994) and the dominant harmonic (Himburg and Friedman, 2006). A study by Lee et al. (2009) in human carotid bifurcations indicated that strong correlations exist between many of these metrics, implying that they capture essentially the same flow features. Similar conclusions have been reached for other vascular segments (Goubergrits et al., 2002; Huo et al., 2008).

Although a large number of studies comparing distributions of these indices – particularly low WSS, high OSI and high RRT – to the localisation of atherosclerotic disease report a good spatial correlation (e.g., Chatzizisis et al. (2011), Samady et al. (2011)), several studies (e.g., Steinman et al. (2002), Gijsen et al. (2007), Peiffer et al. (2012)), including all those making systematic point-by-point comparisons of WSS and disease, did not find the expected relation (Peiffer et al., 2013a). This discrepancy would be plausibly explained if the traditional metrics fail to identify some characteristic of disturbed flow that is of pathobiological importance. One feature haemodynamic metrics generally do not capture is the multidirectionality of the flow field. The idealised examples in Fig. 1 demonstrate that the frequently-used time-averaged WSS (TAWSS), OSI and RRT indices are unable to distinguish between uniaxial and multidirectional flows.

Fig. 1

Three different flow environments to which an endothelial cell could be exposed. The black arrows represent WSS vectors at various times in the cardiac cycle. The grey arrows indicate their evolution with time. The table lists TAWSS (Pa), OSI, RRT (relative to purely forward flow) and transWSS (Pa) for the three environments.

A multidirectional shear stress model of atherogenesis was put forward more than 25 years ago to reconcile the prevailing low and high shear stress theories (McMillan, 1985). The hypothesis was supported by biological and fluid mechanical arguments but was never substantiated. Interest in the multidirectional character of disturbed arterial flow has recently re-emerged. A potential aneurysm formation indicator (AFI), defined as the cosine of the angle between the instantaneous WSS vector and the temporal mean WSS vector, has been used to quantify changes in WSS direction over time (Mantha et al., 2006), whilst a directional oscillatory shear index (DOSI) has been introduced to investigate for in vitro flows the relative oscillatory character of shear stresses in two directions perpendicular to each other (Chakraborty et al., 2012).

The aim of the present study is to improve understanding of multidirectional disturbed flow and its relevance to atherosclerotic disease. We propose a time-averaged scalar haemodynamic metric to capture flow multidirectionality, and analyse distributions of this metric in idealised and anatomically-realistic geometries. We compare the distributions with those of the more traditional TAWSS, OSI and RRT, and comment briefly on the spatial correlation between high levels of the new metric and the occurrence of arterial disease.

Methods

Four haemodynamic metrics, defined below, were computed for the pulsatile flow fields in idealised and realistic arterial geometries previously presented by Kazakidi et al. (2011) and Peiffer et al. (2012), respectively. Geometries and boundary conditions for these computational fluid dynamic analyses are first briefly reviewed.

Idealised branching geometry

The computational domain (Kazakidi et al., 2009, 2011) was a simplified model of a single intercostal artery branching off the aorta: a cylindrical tube, representing the intercostal artery, emerged perpendicularly from a flat surface, representing the dorsal wall of the descending thoracic aorta. The assumption of a flat aortic wall around the junction is justified by the large difference between the diameters of the intercostal artery and the thoracic aorta. Three combinations of boundary conditions at the aortic inflow and side-branch outflow were considered.

Case I: Steady aortic inflow with pulsatile but non-reversing side-branch flow (side-branch peak-to-mean flow ratio of 1.5).

Case II: Steady aortic inflow with pulsatile reversing side-branch flow (side-branch peak-to-mean flow ratio of 2.75, making the flow reverse for approximately one-third of the cycle).

Case III: Pulsatile aortic inflow that reversed near the wall but not at the centreline, with pulsatile non-reversing side-branch flow (peak-to-mean flow ratio of 1.5 in the aorta and the side-branch). Aortic and side-branch flows were in phase.

Reversing intercostal artery flow is seen in almost all young adults but disappears during the 4th decade of life (Sloop et al., 1998). Case I and Case II correspond in this regard to older and younger people, respectively; in these Cases, aortic flow was kept steady to simplify the evaluation of the effects of flow reversal in the branch. The effect of reversing near-wall aortic flow was added in Case III, which is otherwise identical to Case I.

A parabolic velocity profile was imposed at the aortic inlet and along the lateral planes of the computational domain. For pulsatile aortic inflow conditions this profile was complemented by a time-dependent component, consisting of a single sinusoidal harmonic. Further details of the imposed velocity profiles are given by Kazakidi et al. (2011). The mean inflow Reynolds number Rein (based on the aortic diameter) was 500 for all cases. A flow split of 0.79%, based on the average inflow rate, was applied at the side branch outlet, where a Womersley profile was imposed; flow did not reverse near the wall. At the aortic outlet a zero pressure gradient was specified, and a no-slip condition was applied on the aortic and side-branch walls.

Realistic aortic geometry

Peiffer et al. (2012) simulated pulsatile flow in geometries obtained by vascular corrosion casting of two immature and two mature rabbit aortas. A pulsating parabolic velocity profile was imposed at the aortic root – the waveform is shown in the Supplemental Data of Peiffer et al. (2012). Rein (based on the diameter at the aortic root) was set to 300 for all cases, and the vessel walls were non-compliant. A sensitivity analysis for steady flow conditions in this earlier study showed that patterns of WSS are insensitive to whether the inflow velocity profile is blunt, parabolic or skewed, and to secondary flows at the inlet.

Flow splits to large branches were based on experimental measurements or were estimated, as explained in Peiffer et al. (2012), with 14.7% of the flow in the ascending aorta going to the first branches of the aortic arch (the brachiocephalic trunk, left common carotid artery and their branches), 7.1% to the left subclavian artery and its branches, and 2% of the descending aortic flow to the 9 or 10 intercostal arteries. Further flow divisions to the individual branches were calculated using Murray's law (Murray, 1926). Parabolic rather than Womersley velocity profiles were imposed at the outlets of the branches; they pulsated in phase with the aortic flow. Parabolic and Womersely profiles have the same mean component; the higher harmonic content will not have affected WSS patterns in the aorta because of the flow extensions added to the branch arteries. At the outlet of the descending aorta, below the fifth pair of intercostal branches, a zero velocity gradient boundary condition was applied. A no-slip condition was again prescribed at the arterial wall.

Computational methods

Pulsatile flow fields were calculated using the computational fluid dynamics techniques described by Kazakidi et al. (2011) for the idealised branching geometry and by Peiffer et al. (2012) for the anatomically realistic aortic geometries. Volume meshes for use with the spectral/hp element solver Nektar were generated using prismatic boundary elements and tetrahedral internal elements. The mesh for the idealised branching geometry, generated with Felisa non-proprietary software, contained 452 prismatic elements and 2116 tetrahedral elements, corresponding to 266,336 local degrees of freedom per variable when used in combination with a 6th order polynomial expansion within each element. Meshes for the aortic geometries, obtained using TGrid 4.0.24 (ANSYS, Inc.), contained a total of ~105 spectral elements and were used in combination with a 5th order polynomial expansion; equivalent meshes in a finite-element environment would contain ~107 linear tetrahedral elements. Incompressible Newtonian flow was simulated for all cases. At least two cardiac cycles were simulated to dampen initial transients; cycle-to-cycle speed differences below 1% were considered acceptable. Data were extracted at equal time intervals of 20/cycle for the idealised branching geometry and 111/cycle for the anatomically realistic geometries.

Traditional haemodynamic metrics

The TAWSS, OSI (He and Ku, 1996) and RRT (Ku et al., 1985) were calculated using the following definitions:where represents the instantaneous WSS vector, t the time and T the cardiac cycle. Since RRT is a relative concept, the metric was normalised using its value at the aortic inlet for the idealised geometries and using the aortic average for the anatomically representative geometries.

Definition of the transverse WSS

A characteristic of multidirectional disturbed flow is that the WSS vector does not remain parallel to a single axis throughout the cardiac cycle. Only such flows have non-zero WSS components in the direction normal to the temporal mean vector. The calculation of the time-average of the magnitudes of these components results in the transverse WSS (transWSS), according to the following definition:where represents the normal to the arterial surface.

The second form of the definition highlights that transWSS is constructed as the time average of the wall shear stress vector projected in the direction . The quantity physically corresponds to the unit vector within the surface plane that is perpendicular to . TransWSS can take any value between 0 and the TAWSS, so long as there is a preferred flow direction and hence is not the null vector. (Note that the original definition of the OSI considered wall shear stress components that were at either 90 or 180 degrees to the TAWSS vector; whether transverse or reverse components were used depended on the circumferential location of the region of interest within the artery (Ku et al., 1985). However, the OSI has not subsequently been used in this way.)

A low transWSS indicates that the flow, which can be steady, pulsatile or oscillatory, remains approximately parallel to a single axis throughout the cardiac cycle. A high transWSS indicates large changes in flow direction, smaller changes in direction but over a larger portion of the cardiac cycle or small changes in the direction of high-speed near-wall flow. As the examples in Fig. 1 suggest, transWSS does not differentiate between purely forward and reversing unidirectional flow. It therefore does not replace the TAWSS, OSI or RRT, but rather complements them.

Results

The distributions of TAWSS, OSI, RRT and transWSS for the three Cases are shown in Fig. 2. As detailed by Kazakidi et al. (2011), for steady aortic flow in combination with pulsatile but non-reversing side-branch flow (Case I), TAWSS was low in a patch slightly downstream of each side of the orifice and was elevated around the ostium (especially upstream) and in a separate patch downstream of it. When side-branch flow reversed (Case II) or when near-wall aortic flow reversed (Case III), the location and/or extent of the downstream patch of high TAWSS changed somewhat but the overall pattern of the metric remained essentially the same.

Fig. 2

Distributions of four haemodynamic metrics on the aortic wall of the idealised branching geometry presented by Kazakidi et al. (2009, 2011). The different rows represent different sets of boundary conditions: steady aortic flow with non-reversing (Case I) or reversing (Case II) side-branch flow, or pulsatile aortic flow with non-reversing side-branch flow (Case III). Net aortic flow is from top to bottom in each plot and the white disc represents the branch mouth, viewed en face. TAWSS⁎ and transWSS⁎ were normalised using the TAWSS at the aortic inlet; RRT⁎ was normalised using the RRT at the aortic inlet.

The map of RRT looked like the inverse of the map of TAWSS; the overall pattern was again generally insensitive to the changes in boundary conditions, with the same caveat about the location and relative extent of the patch downstream of the ostium. The maps of OSI for Case I and Case II were similar in character to each other, although the reversing side-branch flow in Case II increased the area where this metric was elevated. The highest values occurred in three patches in a transverse line downstream of the ostium; a four-lobed area of intermediate values surrounded this line. In Case III, the reversal of near-wall aortic flow led to non-zero values of OSI almost everywhere in the map. Even in this Case, however, a superimposition of the three peaks in a transverse line, seen in the other two maps, is apparent. (Note that the OSI pattern for Case II shown in Fig. 6(b)(ii) of Kazakidi et al. (2011) was based on a limited number of steps in the cardiac cycle and gave a different pattern from the more rigorous one shown here.)

For all three Cases, the map of transWSS showed an overall “butterfly” pattern, albeit of varying size and detailed features. However, within this shape the location of the highest values was substantially different in Case II from Cases I and III. In the latter Cases, transWSS was elevated in a circumscribed region downstream of the ostium; the distance from the branch was similar to that seen for the maxima in OSI. In Case II, on the other hand, the reversal in side-branch flow displaced the transWSS peaks to the lateral margins of the branch. Thus not only did the pattern of transWSS differ from that of more traditional metrics, as expected from its definition, but it also seemed more sensitive to the alteration in velocity waveform.

Rabbit aortic flow

Fig. 3 shows the path of the tip of the WSS vector during the cardiac cycle for four points in a rabbit aorta. The directionality of the WSS vector varied greatly with location. The flow was reversing but confined to a single line in the distal thoracic aorta (point d). A similar pattern was seen for the proximal region of the same aortic segment, but there were slightly greater deviations from the line of mean flow (point c). In the aortic arch region the flow was truly multidirectional (points a and b).

Fig. 3

Traces showing changes in the direction and magnitude of the WSS vector over the cardiac cycle at four points (a–d) in a rabbit aorta. (The locations of the four points are indicated on the TAWSS map in Fig. 4.) Each line shows the path traced during the cardiac cycle by the tip of the vector, which has its origin at the centre of the circle. For example, the vector at point “d” in the descending thoracic aorta oscillates between forward and reverse directions along the aortic axis (0–180°), whereas the vector at point “b” in the aortic arch is more multidirectional in character and reaches greater magnitudes. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

In Fig. 4 the distributions of the TAWSS, OSI, RRT and transWSS are shown for the same rabbit aorta. As described by Peiffer et al. (2012), TAWSS was low, and OSI and RRT were high, along the lesser and greater curvature of the ascending aorta, on the inner wall of the aortic arch and along two stripes in the descending aorta proximal to the first intercostal branch pair. (An inverted colour scale has been used for TAWSS to make this clearer.) The latter regions were narrower in the OSI map than in the TAWSS and RRT maps. TransWSS was high in areas of the ascending aorta characterised by high OSI and RRT and low TAWSS, but after that its pattern diverged from those of the other metrics. Thus in the aortic arch, transWSS was highest on the lateral walls, and the two stripes proximal to the first intercostal branch pair seen for the other metrics were replaced by a single dorsal stripe in this animal. In the descending thoracic aorta, multidirectional flow was found in two stripes that included the left and right intercostal branch ostia, and not in the single ventral stripe seen for high OSI, high RRT and low TAWSS; it was highest at the level of the most proximal intercostal ostium. A close-up of this region is shown in Fig. 5. To a reasonable approximation, the patterns of OSI and RRT are the inverse of the pattern for TAWSS. TransWSS showed a strikingly different distribution both around and away from the branch.

Fig. 4

Distributions of four haemodynamic metrics in a rabbit aorta. The luminal surface of each segment, opened ventrally whilst preserving tapering, is shown en face. The colour scale of the TAWSS map was inverted for ease of comparison. The letters on this map refer to the selection of points used in Fig. 3. Net aortic flow is from top to bottom. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 5

Close-up of the proximal left intercostal branch ostium in Fig. 4, viewed from outside the vessel, showing distributions of four haemodynamic metrics. The colour scale of the TAWSS map was inverted for ease of comparison. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 6 shows the transWSS distributions for two immature and two mature animals. (The map on the left corresponds to the transWSS map in Fig. 4.) In all four animals, two streaks of elevated transWSS are clearly visible in the descending section. They appear to start further upstream in the mature than the immature rabbits. The variable spatial relation between the streaks and the intercostal ostia gave rise to a range of different transWSS patterns around the ostia. Of particular interest, when an ostium was situated well within a longitudinal stripe of high transWSS, its lateral and downstream margins were often surrounded by an arrowhead-shaped region of even higher transWSS. When an ostium was situated between the two stripes of high transWSS, lobes of high transWSS connected the ostium to the stripes.

Fig. 6

Distribution of transWSS in aortas of two immature and two mature rabbits. The luminal surface of each segment, opened ventrally whilst preserving tapering, is shown en face. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Discussion

The concept that disturbed blood flow plays a role in the initiation and development of atherosclerosis is well established; however, the imprecise nature of the term is acknowledged (e.g., by Himburg and Friedman (2006)). As noted in the Introduction, two recent studies have developed new indices to capture the multidirectional nature of disturbed flow. The DOSI was developed to predict the patchy elongation of endothelial cells exposed to complex flow patterns in vitro (Chakraborty et al., 2012). It focuses on the relative oscillatory character of the flow in two directions perpendicular to each other, unlike the transWSS metric introduced here which considers the absolute magnitude of oscillations perpendicular to the predominant flow direction. The AFI (Mantha et al., 2006), developed to predict aneurysm formation, is defined as the cosine of the angle between the instantaneous and temporal mean WSS vectors. (The cosine of this angle has also been used in a modified OSI definition by Papaharilaou et al. (2002).) Note that the AFI considers only the direction of the shear, whilst its magnitude is ignored; transWSS accounts for both. Moreover, because cosines are monotonic between 0° and 180°, the AFI gives more weight to reverse flow than to flow which is perpendicular to the mean vector. Hence purely pulsatile flow could have an extreme AFI value, although its transWSS would always be zero. Finally, the value of AFI varies during the cardiac cycle and it is unclear which value, or what combination of values, should be compared to aneurysm prevalence or other disease indices.

A theoretical limitation to the transWSS metric is that it cannot be defined unless there is a dominant flow direction. Purely oscillatory flow and purely circular flow are two examples where it could not be defined. If such cases were analysed computationally, a dominant direction would be assigned randomly due to numerical errors in the calculations, and hence an artefactual transWSS would be generated. The definition could be modified to eliminate this mathematical incompleteness. For real arterial flows, however, the probability of this problem arising is not sufficiently high to justify the additional complication. A more plausible problem can be illustrated by the hypothetical example of a large forward flow along the aortic axis in systole, a slightly smaller but still large reverse axial flow in diastole, and a small circumferential flow occurring in a period when the axial flows are halted. In this case, the two large axially-aligned flows would nearly cancel when averaged, leading to the small circumferential component being identified as the mean flow direction; the axial flows would then give rise to a large transWSS. Examination of our data suggests that errors of this type are also unlikely to occur. Fig. 4 shows that the large majority of OSI values in the realistic geometry fall below 0.25 (a purely oscillatory flow would give 0.5), whilst comparison of the OSI and transWSS columns of Fig. 2 does not show any artefactual peaks of transWSS in regions of peak OSI, despite the fact that OSI values in this idealised model reach higher values.

Evaluation of transWSS in an idealised branching geometry showed that the pattern of transWSS differs from the patterns of traditional haemodynamic metrics. Moreover, the pattern of transWSS but not the patterns of the other metrics was influenced by the boundary condition imposed at the side branch: for pulsatile but non-reversing side-branch flow the transWSS reached its maximum a short distance downstream of the branch opening, whereas the maximum was on the lateral margins of the ostium for reversing side-branch flow. (An inference from this observation is that patterns of transWSS obtained for the realistic geometry may be more sensitive than patterns of the other metrics to the necessary approximations made when setting inflow and outflow boundary conditions). By combining the information available in the maps of the TAWSS, OSI, RRT and transWSS, the governing flow environment in different regions of the branching geometry can be characterised. For example, a circumscribed region of multidirectional flow was seen downstream of the ostium (high transWSS, high OSI and low TAWSS); this is where flow stagnation occurs during at least part of the cycle (Kazakidi et al., 2011). To the sides of the ostium, the near-wall flow deviates from its mean direction without reversing (high transWSS but low OSI).

TransWSS also showed distinctive patterns in anatomically-realistic models of the thoracic aorta. Maps for all four aortas were dominated by two stripes of elevated transWSS running along or adjacent to the left and right intercostal branch ostia. This pattern was not seen in maps of TAWSS, OSI or RRT. Another difference is that, whereas distributions of the traditional haemodynamic metrics did not vary notably from one intercostal branch to another, a range of different near-branch transWSS patterns was found. This may reflect variation in the location of the branch with respect to the two stripes, or a more critical dependence of transWSS on branch geometry. Lesion distributions also vary from ostium to ostium within the same subject (Cremers et al., 2011).

Lesion distributions around branches additionally change with age: diet-induced lipid deposition occurs most frequently in an arrowhead-shaped region centred on the downstream half of intercostal branch mouths in immature rabbit aortas but to the sides of branches at later ages. Typical lesions at each age are shown in Fig. 7, as are two near-branch patterns of transWSS obtained in the present study. The striking similarity between patterns of lesions and patterns of transWSS suggests that the latter may play a causal role in atherogenesis. Since comparable patterns of lipid deposition have been reported in children and young adults (Sinzinger et al., 1980; Sloop et al., 1998), the findings may have implications for human as well as rabbit disease.

Fig. 7

Examples of diet-induced lesions around intercostal branch ostia in immature (top left) and mature (bottom left) rabbit aortas, lipid deposition being stained red and normal tissue counterstained blue (Cremers et al., 2011), and equivalent maps of transWSS (immature, top right; mature, bottom right), enlarged from Fig. 6 of the present study. Net aortic flow is from top to bottom in all figures. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Although it can be seen in Fig. 6 that neither the downstream nor the lateral pattern of transWSS is unique to a single age group, there is a trend for the former to occur in immature vessels and the latter in mature ones. The dependence of the two patterns on boundary conditions and the possibility that critical boundary conditions change with age require further investigation. Results obtained in the idealised geometry suggest that the degree of flow reversal may be important; this feature is known to change with age (Sloop et al., 1998; Padilla et al., 2011).

Conclusion

A number of metrics have been developed to quantify disturbed blood flow. The transWSS metric described here does so in a novel way by time averaging the magnitude of WSS components perpendicular to the mean shear vector at each location on the vessel wall. Its distribution in an idealised model and in the rabbit aorta differs noticeably from that of the TAWSS, OSI and RRT, which are themselves highly correlated. Given the postulated importance of disturbed flow in the predilection of atherosclerosis for branch sites, the striking correspondence between selected lesion and transWSS patterns around intercostal ostia is noteworthy. A recent statistical comparison (Peiffer et al., 2013b) did not show a reliable correlation between lesion prevalence and TAWSS patterns, consistent with a review of the literature that also questioned this widely-assumed association (Peiffer et al., 2013a). Future studies using large samples will establish whether there is a stronger correlation between lesion prevalence and transWSS patterns, and will determine whether the change with age seen in lesion patterns can be explained by influences of altered boundary conditions (Sloop et al., 1998) or aortic geometry (Peiffer et al., 2012) on transWSS. If an association is found, it will be important to develop in vitro methods, such as those based on orbital flows (Potter et al., 2011; Chakraborty et al., 2012), for investigating influences of transWSS on endothelial cells.

Conflict of interest statement

None.