Sang-Wook Lee Biomedical Simulation Laboratory,University of Toronto,5King’s College Road Toronto,Toronto,ON M5S3G8Canada;School of Mechanical and AutomotiveEngineering,University of Ulsan,Ulsan680-749,South KoreaLuca AntigaDepartment of Bioengineering, Mario Negri Institute for PharmacologicalResearch,24020Ranica(BG),Italy David A.Steinman1Biomedical Simulation Laboratory,University of Toronto,5King’s College Road Toronto,Toronto,ON M5S3G8Canada e-mail:steinman@mie.utoronto.ca Correlations Among Indicators of Disturbed Flow at the Normal Carotid BifurcationA variety of hemodynamic wall parameters(HWP)has been proposed over the years to quantify hemodynamic disturbances as potential predictors or indicators of vascular wall dysfunction.The aim of this study was to determine whether some of these might,for practical purposes,be considered redundant.Image-based computationalfluid dynamics simulations were carried out for Nϭ50normal carotid bifurcations reconstructed from magnetic resonance imaging.Pairwise Spearman correlation analysis was performed for HWP quantifying wall shear stress magnitudes,spatial and temporal gradients,and harmonic contents.These were based on the spatial distributions of each HWP and,harmonic(DH)parameter were found to depend on how the wall shear stress magnitude was defined in the presence offlow reversals.Many of the proposed HWP were found to provide essentially the same information about disturbedflow at the normal carotid bifurcation.RRT is recommended as a robust single metric of low and oscillating shear. On the other hand,gradient-based HWP may be of limited utility in light of possible redundancies with other HWP,and practical challenges in their measurement.Further investigations are encouraged before thesefindings should be extrapolated to other vas-cular territories.͓DOI:10.1115/1.3127252͔Keywords:wall shear stress,atherosclerosis,hemodynamic wall parameter,carotid bifurcation1IntroductionThere is much evidence suggesting that initiation and progres-sion of atherosclerotic disease is influenced by“disturbedflow”͓1͔.Notwithstanding the imprecise nature of this term͓2͔,variousmetrics have been proposed over the years to quantifyflow dis-turbances.Originally focused on the magnitudes of wall shear stress͑WSS͓͒3,4͔these hemodynamic wall parameters͑HWP͒have since incorporated spatial and temporal gradients of WSS ͓5–8͔and,more recently,the harmonic content of time-varying WSS waveforms͓2,9͔.In a recent computationalfluid dynamics͑CFD͒study of the relationship between geometry and disturbedflow at the carotid bifurcations of young adults͓10͔,we noted that ourfindings were relatively insensitive to the choice of either time-averaged wall shear stress magnitude͑TAWSS͒or oscillatory shear index͑OSI͒as metrics of disturbedflow.This was found to be explained by a strong and significant inverse correlation between these two quan-tities.Such correlations among HWP are not unexpected,as rec-ognized early by Friedman and Deters͓11͔;however,they have been little-investigated in light of the growth in the number and complexity of candidate HWP.With this in mind,the objective of the present study was to use a representative sample of normal carotid bifurcation geometries to comprehensively test for correlations among established and recently-proposed HWP.Especially in the context of large-scalestudies of so-called geometric and hemodynamic risk factors inatherosclerosis,we aimed to determine whether a subset of HWP,or even a single HWP,might serve as a sufficiently robust markerof disturbedflow.2Materials and Methods2.1Computational Fluid Dynamics.N=50anatomically re-alistic carotid bifurcation geometries were digitally reconstructedfrom black blood magnetic resonance imaging͑MRI͒of25osten-sibly healthy young adults,as described previously͓12͔.CFDsimulations were carried out using a well-validated in-housefinite-element-based CFD solver͓13–15͔.Quadratic tetrahedral-element meshes were generated by a commercial mesh generator ͑ICEM-CFD;ANSYS,Berkeley,CA͒using a nominally uniform node spacing of0.2mm,previously shown to be sufficient forresolving wall shear stresses to within10%accuracy͓16͔.Rigidwalls and Newtonian rheology were assumed.Pulsatileflowboundary conditions were prescribed based on representativewaveform shapes and allometrically-scaled inlet and outletflowrates.Further details of the CFD simulations are provided else-where͓10͔.For each tetrahedral element the vector WSS,␶w,was calcu-lated as the projection of the stress tensor onto the element’s sur-face at each node,using the element’s quadratic shape functions. As nodes are connected to multiple elements,contributions to each nodal␶w were averaged together.From these time-varying nodal WSS vectors,a variety of HWP were computed,as summa-rized in Table1,and detailed below.1Corresponding author.Contributed by the Bioengineering Division of ASME for publication in the J OUR-NAL OF B IOMECHANICAL E NGINEERING.Manuscript received August12,2008;finalmanuscript received January1,2009;published online May11,2009.Review con-ducted by Fumihiko Kajiya.Paper presented at the2008Summer Bioengineering Conference͑SBC2008͒,Marco Island,FL,June25–29,2008.2.2Magnitude-Based HWP.Time-averaged wall shear stress was calculated by integrating each nodal WSS magnitude over the cardiac cycle.For each CFD model,the nodal TAWSS were normalized by the fully-developed ͑i.e.,Poiseuille ͒value,based on the model’s imposed cycle-averaged flow rate,the as-sumed viscosity,and the mean diameter ofthe model’s commoncarotid artery ͑CCA ͒at a location three radii upstream of the bifurcation ͑i.e.,CCA3,as defined in Ref.͓10͔͒.Oscillatory shear index,a dimensionless metric of changes in the WSS direction,originally introduced by Ku et al.͓6͔,was malized with respect to its fully-developed value at CCA3.2.3Gradient-Based HWP.Originally proposed by Lei et al.͓19͔,the wall shear stress spatial gradient ͑WSSG ͒may be con-sidered a marker of endothelial cell tension.As shown in Table 1,it is calculated from the WSS gradient tensor components parallel and perpendicular to the time-averaged WSS vector ͑m and n ,respectively ͒.Here the WSS gradients were calculated directly from the velocities of each element,taking advantage of their quadratic shape functions.As with WSS itself,elemental contri-butions to nodal WSSG values were averaged together.Since the WSSG of fully-developed flow in a straight uniform pipe is zero,here it is normalized by the fully-developed TAWSS divided by the mean diameter at CCA3.Subsequently,Longest and Kleinstreuer ͓20͔proposed the WSS angle gradient ͑WSSAG ͒to highlight regions exposed to large changes in WSS direction,irrespective of magnitude.As indicated by the formula in Table 1,this was done by calculating,for each element’s node ͑index j ͒,its direction relative to some reference WSS vector ͑index i ͒,here chosen to be that at the element’s centroid.The WSS temporal gradient ͑WSST ͒,originally suggested by Ojha ͓7͔as a factor in distal anastomotic intimal hyperplasia,is simply the maximum absolute rate of change in WSS magnitude over the cardiac cycle.Here it is normalized to WSST at location CCA3,determined from Womersley’s solution of fully-developed pulsatile flow based on the CCA3diameter,the heart rate,and the imposed flow rate Fourier coefficients.2.4Harmonic-Based HWP.Recently,Himburg and Fried-man ͓2͔suggested the harmonic content of the WSS waveform as a possible metric of disturbed flow,subsequently linking this to the frequency-dependent responses of endothelial cells ͓21͔.Fol-lowing those authors,the time-varying WSS magnitude at each node was Fourier-decomposed,with the dominant harmonic ͑DH ͒simply defined as the harmonic with the highest amplitude.Around the same time,Gelfand et al.͓9͔defined the harmonic index ͑HI ͒as the relative fraction of the harmonic amplitude spec-trum arising from the pulsatile flow components.Table 1Definitions of hemodynamic wall parameters …HWP …2.5Data Analysis.Correlations among the various HWP were evaluated in two ways:locally,whereby spatial distributionsof HWP were compared;and globally,whereby overall“burdens”of disturbedflow were compared.Local correlations were tested byfirst discretizing each model’ssurface into afinite number of contiguous patches,within each ofwhich the respective HWP was averaged͑for the ordinal DH,themedian was used͒.This amounted to mapping each distribution ofHWP onto an objective template plane,fixed with respect to abifurcation-specific coordinate system͓22͔.Remembering thatmodel dimensions have previously been normalized with respectto their respective CCA3radius,patches were nominally0.5unitsin length along the direction of the vessel axis,with12patchesdistributed circumferentially at each axial level.To facilitate pool-ing of patches from all cases,and in light of differences in eachCFD model’s relative length,only those patch locations presentfor all CFD models were included,resulting in288patches permodel.Although absolute patch sizes may have varied betweenand within models,this mapping procedure ensured that all sur-faces were spatially discretized in a consistent way͓22͔.Global comparisons were made by calculating,for each CFDmodel,the fraction of its surface area exposed to“abnormal”val-ues of a particular HWP.As detailed in Ref.͓10͔,“abnormal”wasdefined,for each HWP,as the upper quintile͑for TAWSS,thelower quintile͒of the HWP data pooled from all CFD models,thus providing an objective threshold between normal and abnor-mal.To ensure a consistent spatial extent across all cases,thesurfaces were clipped at planes three andfive maximally inscribedsphere radii along the common and internal carotid arteries ͑CCA3and ICA5,respectively͓10͔͒.Then,for each case,the area of the surface experiencing HWP above͑for TAWSS,below͒thethreshold was calculated.To factor out the influence of vesselsize—for the same shape,a larger vessel will experience a largerarea of disturbedflow—this absolute surface area was divided bythe total͑clipped͒surface area.In this way,for each HWP,eachCFD model was assigned a single value characterizing how muchof its surface was exposed to disturbedflow.For both local and global comparisons,a Spearman rank corre-lation coefficient͑r͒and significance͑p-value͒were computed for each of the28unique pairs of HWP using PRISM version4͑GraphPad Software,San Diego,CA͒.Correlations having p Ͻ0.05were deemed strong for͉r͉Ͼ0.8,weak for͉r͉Ͻ0.5,and moderate in between.Spearman correlation analysis was chosen in part because HWP distributions are unlikely to be normal.Also, in assessing these correlations based on rank,we sought to iden-tify,in the local correlations,whether the sites of extrema for one HWP would be reflected in the sites of extrema for another HWP. Global correlations sought to identify whether the ranking of cases from low to high burden of“disturbedflow,”based on the threshold of a given HWP,would be the same as that obtained based on another HWP.3ResultsAs depicted in Fig.1for a representative case,disturbedflow based on WSS magnitude quantities͑TAWSS,OSI,and RRT͒was concentrated around the outer walls of the bifurcation,consistent with many previous observations.For gradient-based HWP ͑WSSG,WSSAG,and WSST͒elevated values were concentrated around the bifurcation apex and,to a lesser but more variable extent,around the external and internal carotid artery͑ECA and ICA͒branches.Distributions of harmonic HWP͑DH and HI͒were more distinctive:higher DH was concentrated away from the outer walls of the bifurcation,whereas elevated HI reflected the general locations,if not the specific spatial extents,of the magnitude-based HWP.Overall,these observations hint at the correlations among the patched HWP distributions,detailed in Table2and depicted graphically for selected HWP pairs in Fig.2.Strong correlations were seen between TAWSS and both RRT͑r=−0.99͒and WSSG ͑r=0.86͒,albeit for different reasons:regions of elevated RRT correlated well with those experiencing low TAWSS͑r=−0.99͒, whereas elevated WSSG correlated with elevated TAWSS͑r =0.86͒.Moderate inverse correlations were found between TAWSS and both OSI and HI͑r=−0.66and r=−0.72,respec-tively͒,whereas TAWSS was positively correlated with WSST ͑r=0.63͒.Although many of the correlations were weak,all were statistically significant.Of all HWP,only DH was neither stronglynor moderately correlated with any other HWP.It is also worthnoting that the correlations identified by pooling the cases werefairly consistent across the50cases analyzed individually,as evi-denced by the relatively narrow confidence intervals shown in thesame table.According to the global correlations summarized in Table3,TAWSS,OSI,or RRT would rank vessels in similar order fromlow to high burdens of disturbedflow,as indicated by the strongpositive global correlation coefficients.Conversely,if disturbedflow was defined as elevated WSSG,vessels would be ranked inreverse order to this,as indicated by the moderate negative corre-lations with TAWSS,OSI,and RRT.As with the local correla-tions,HI was moderately correlated with the trio of magnitude-based HWP,while DH was only weakly correlated with otherHWP.Overall,corresponding local and global correlations wereof similar strength.A notable exception was OSI versus WSSAG, for which the local correlation was moderate͑r=0.73͒,while the global correlation was weak͑r=0.05͒and not significant.Reasons for this are given in Sec.4.It is worth noting that the above results were found to be rela-tively insensitive to the choice of data analysis method.For ex-ample,the Spearman correlation analysis of the continuous distri-butions of HWP͑i.e.,the CFD nodal values prior to patching͒revealed correlations similar to those obtained after patching.Glo-bal correlations based on a90th percentile threshold for disturbedflow were similar to those reported here using the80th percentilethreshold,afinding consistent with Ref.͓10͔.Finally,in light ofthe concentrations of HWP extrema around the bifurcation region,we repeated the local correlation analysis using patches extendingaxially only halfway along each branch͑i.e.144patches permodel͒to exclude the distal parts of the branches.Again thetrends were the same,although some of the moderate correlationsactually increased in strength,such as TAWSS versus OSI͑from r=−0.66to r=−0.77͒and TAWSS versus WSST͑from r=0.62to r=0.75͒.The detailed results are omitted in the interest of space, and because they do not affect the implications and conclusions discussed below.4Discussion4.1Summary and Implications of Findings.This compre-hensive evaluation of correlations among HWP at the normal ca-rotid bifurcation clearly demonstrates that some of these param-eters may,for practical purposes,be considered redundant.By virtue of their definitions in Table1,correlations among WSS, OSI,and RRT were expected,although not necessarily at the strengths observed.While WSS and OSI were moderately corre-lated,Fig.1suggests the OSI captures apparentflow disturbances at the ECA branch.Notwithstanding whether these are significant in the context of atherosclerosis—plaques do tend to occur at the ICA—our results would suggest that RRT can replace WSS and OSI as a single marker of“low and oscillatory”shear.In fact,as noted earlier,RRT is,by definition,the inverse of the magnitude of the time-averaged WSS vector.This explains its near-perfect correlation with TAWSS,which,remember,is the time-average of the WSS magnitude.In other words,RRT is simply another type of time-averaged WSS,but inverted and with a more tangible connection to the biological mechanisms underlying atherosclero-sis͓18͔.To appreciate the practical implications of replacing TAWSS and OSI with RRT,consider our recent study in which exposure todisturbed flow was found to be significantly predicted by a com-bination of bifurcation area ratio and tortuosity ͓10͔.There,the findings were shown to be independent of the choice of TAWSS or OSI as the metric of disturbed flow.Here,repeating the multiple regressions using RRT above the 80th percentile as the criterion for disturbed flow,we found near-identical—if anything,slightly stronger—coefficients:R adj 2=0.341͑p =0.0001͒,␤tortuosity =−0.498͑p =0.0001͒,and ␤AR1=0.459͑p =0.0007͒.In other words,exposure of the vessels to “disturbed flow”is the same,whether defined by extrema of TAWSS,OSI,or RRT.The strong positive local correlations ͑and consequent strong negative global correlations ͒between TAWSS and its spatial and temporal gradients likely reflect the fact that all of these quantities are highest around the apex of the bifurcation.As pointed out by Ojha ͓7͔,the use of WSS spatial gradients as risk indicators for intimal thickening is questionable in light of their concentration about the bifurcation apex,a region usually spared of plaques.As suggested recently by Goubergrits et al.͓23͔,regions elsewhere experiencing elevated WSSG may represent a consequence of ath-erosclerosis rather than a cause.Either way,for the normal carotid bifurcation at least,our findings would suggest that TAWSS could be used instead of WSSG,which is anyway more susceptible to measurement uncertainty ͓24,25͔,owing to its reliance on spatial gradients.A similar conclusion may be drawn from the moderate correla-tions between TAWSS and WSST,although it is worth remember-ing that in this study all CFD models were exposed to the same waveform shape .Intersubject variations in waveform shape are reported to be on the order of 10%͓26,27͔.Such variations in flow rate dynamics have been found to have a relatively minor influence on variations in the distributions of a variety of HWP,at least relative to the influence of uncertainty in the reconstructed geometry ͓25͔.Thus,it is reasonable to conclude that our findings are robust to our assumptions about waveform shape.It was also observed that the spatial distributions of OSI and WSSAG were moderately correlated ͑r =0.73͒,consistent with previous qualitative observations for carotid bifurcations ͓25͔and coronary arteries ͓23͔.As pointed out by Goubergrits et al.͓23͔,WSSAG may be thought of as an extension of OSI;however,being based on differential versus integrated quantities,WSSAG distributions tend to be noisier and more sensitive to uncertainty,something evident here and also in Ref.͓25͔.Nevertheless,for the representative case presented in Fig.1,this correlation is less obvious.Elevated values of WSSAG do coincide with the periph-ery of those regions exposed to elevated OSI;however,the core region of elevated OSI at the carotid bulb is characterized by low WSSAG and,as with the other gradient-based HWP,elevated WSSAG is concentrated at the bifurcation apex,a regiontypicallyFig.1HWP distributions for a representative case.Except for the ordinal DH,contour levels de-picted in each frame’s legend correspond to the 80th,85th,90th,and 95th percentile values based on the HWP distribution pooled over all cases.Note identification of CCA3and ICA5clip planes in the upper left …TAWSS …panel.spared of atherosclerosis.Moreover,the global correlation of these quantities was much weaker ͑r =0.05͒.This may be ex-plained in reference to the respective scatter plot in Fig.2,which clearly shows that the local ͑i.e.,patchwise ͒Spearman rank cor-relation was biased by the preponderance of patches having low OSI and WSSAG values.Focusing only on those patches having OSI Ͼ0.1,it can be seen that there is no obvious correlation with WSSAG.Because the global correlations focused only on those regions exposed to the upper quintile of HWP values,they better reflect the correlation,or lack thereof,of these extrema.Having said this,it is worth noting that,for most of the other pairwise comparisons,global and local correlation coefficients were in much closer agreement.Of all HWP,only DH was found to be essentially independent of the other HWP.This was somewhat surprising,since Himburg and Friedman ͓2͔,in introducing the use of WSS harmonics as metrics of disturbed flow,reported an inverse correlation between DH and TAWSS ͑Pearson r =−0.62͒.By way of explaining this,we note that their study was carried out on porcine iliac arteries,which are nominally straight vessels experiencing largely axial flows.On the other hand,flow at the carotid bifurcation is decid-edly nonaxial,and likely subject to more reverse flow.In such regions,rectification of the time-varying WSS vector—remember that DH was derived here from the WSS vector magnitude,per its original definition ͓28͔—could alter its harmonic content.To appreciate the impact of this,we recomputed DH and HI using instead the time-varying “axial”WSS,namely,the compo-nent of the instantaneous WSS vector projected onto a unit vector defined by the direction of its time-averaged value.In this way,flow reversals relative to the nominal axial direction are pre-Table 2Spearman rank correlation coefficients for pairwise local comparisons of the pooled HWP distributions …N =14,400patches total ….Shown in brackets are 95%confidence intervals,drawn from pairwise comparisons of the 50datasets individually …N =288patches each ….All correlations significant to p <0.001,except WSSAG versus DH,p =0.03.OSIRRT WSSG WSSAG WSST DH HI TAWSS Ϫ0.66͓Ϫ0.76,Ϫ0.53͔Ϫ0.99͓Ϫ0.99,Ϫ0.98͔0.86͓0.80,0.91͔Ϫ0.29͓Ϫ0.41,Ϫ0.16͔0.63͓0.46,0.79͔Ϫ0.33͓Ϫ0.45,Ϫ0.21͔Ϫ0.72͓Ϫ0.81,Ϫ0.60͔OSI 0.72͓0.59,0.82͔Ϫ0.38͓Ϫ0.62,Ϫ0.16͔0.73͓0.68,0.78͔Ϫ0.27͓Ϫ0.47,Ϫ0.04͔0.06͓Ϫ0.08,0.16͔0.53͓0.39,0.66͔RRT Ϫ0.81͓Ϫ0.88,Ϫ0.72͔0.38͓0.25,0.52͔Ϫ0.57͓Ϫ0.74,Ϫ0.41͔0.30͓0.19,0.43͔0.74͓0.65,0.82͔WSSG 0.09͓Ϫ0.11,0.26͔0.67͓0.52,0.82͔Ϫ0.31͓Ϫ0.47,Ϫ0.14͔Ϫ0.51͓Ϫ0.67,Ϫ0.26͔WSSAG 0.07͓Ϫ0.13,0.24͔Ϫ0.02͓Ϫ0.15,0.10͔0.45͓Ϫ0.32,Ϫ0.55͔WSST Ϫ0.14͓Ϫ0.26,Ϫ0.04͔Ϫ0.08͓Ϫ0.38,0.15͔DH0.27͓0.14,0.41͔Fig.2Scatter plots for selected pairwise comparisons of HWP .Note that the local …patched …data are plotted using a log-log scale to better depict the full dynamic range of data.served.As can be seen by comparing Fig.3to Fig.1,this had a marked effect on the spatial distributions of DH and HI.The rea-son for this is also given in Fig.3:rectification of the instanta-neous WSS served to break up the clear fifth harmonic oscillation of the time-varying axial WSS in favor of the lower frequencies.As summarized in Table 4,this served to strengthen the correla-tions between the harmonic and other HWP.Of note is the local correlation coefficient for TAWSS versus DH ͑r =−0.69͒,now close to that originally reported by Himburg and Friedman.Whether DH and HI should be defined based on magnitude or axial WSS cannot be answered by the present study.It is also not clear whether DH should be considered a monotonic HWP in light of possible limits to the temporal response of endothelial cells to shear ͓21͔.Nevertheless,our findings do bring to attention a heretofore-underappreciated issue in the harmonic analysis of WSS in the presence of strongly nonaxial flow.4.2Potential Limitations.This study has made the custom-ary assumptions of rigid walls,Newtonian rheology,and fully-developed inlet boundary conditions,previously shown to be of relatively minor influence on the distribution of WSS ͓16,29,30͔.The use of allometrically-scaled flow rates was recently shown to have little impact on the relative burdens of disturbed flow among the 50cases considered here ͓10͔;however,as noted earlier,theassumption of a constant waveform shape might have served to underestimate variations in WSST,as well as the harmonic HWP.An obvious limitation of this study is the relatively narrow scope of vascular configurations considered,namely,normal ca-rotid bifurcations.We do note,however,that recent work from Goubergrits et al.͓23͔similarly reported possible redundancies among gradient and magnitude-based HWP for the case of a nor-mal coronary artery.Huo et al.͓31͔also noted a significant power-law relationship between TAWSS and OSI on the outer walls of the common carotid and celiac arteries,based on a CFD model of flow along the length of the mouse aorta.While a power-law relationship is expected based on the definition of OSI,those au-thors did note a difference in the power-law coefficients derived from the carotid and celiac sites,an observation consistent with the branch-specific clustering of data points in the OSI versus TAWSS scatter plot in Fig.2.Those authors also noted a corre-spondence between elevated TAWSS and elevated WSSG,al-though no quantitative relationship was found.Nevertheless,we encourage further investigation before our findings should be ex-trapolated to other vascular territories.It is also important to appreciate that our study has made no attempt to prioritize any of these HWP in terms of their purportedTable 3Spearman rank correlation coefficients for each pair-wise global comparison of the relative surface area exposed to the HWP beyond its 80th percentile value …N =50cases ….OSIRRT WSSG WSSAG WSST DH HI TAWSS 0.80a0.97aϪ0.63a Ϫ0.16Ϫ0.250.040.68a OSI 0.88aϪ0.50a 0.05Ϫ0.260.040.66a RRT Ϫ0.63aϪ0.10Ϫ0.250.020.74a WSSG 0.58a0.54a Ϫ0.33b Ϫ0.27WSSAG 0.43cϪ0.240.20WSST Ϫ0.20Ϫ0.04DH0.01ap Ͻ0.05.bp Ͻ0.01.cp Ͻ0.001.Fig.3Distributions of DH and HI based on the axial WSS component rather than WSS magnitude,shown for same case depicted in Fig.1.The arrows indicate the site of the time-varying WSS waveforms and corresponding spectra shown to the right.Table 4Spearman rank correlation coefficients for local and global comparisons of DH and HI,as derived from the time-varying axial WSS rather than the WSS magnitude.Local GlobalDHHI DH HI TAWSS Ϫ0.69a Ϫ0.77a 0.76a 0.83a OSI 0.55a 0.65a 0.78a 0.89a RRT 0.70a 0.81a 0.78a 0.89a WSSG Ϫ0.58a Ϫ0.54a Ϫ0.54b Ϫ0.46a WSSAG 0.31a 0.51a Ϫ0.180.12WSST Ϫ0.36aϪ0.14a Ϫ0.34cϪ0.15DH0.63a0.64aap Ͻ0.05.bp Ͻ0.01.cp Ͻ0.001.links to the underlying biological mechanisms.Thus,some of the HWP we regard here as redundant might be shown to have closer mechanistic links in studies of individual biological responses to the local hemodynamic environment.Rather,we suggest that our findings are most applicable to large-scale studies of hemody-namic factors in atherosclerosis,which are more concerned with quantifying overall burdens or identifying patterns of localization of disturbedflow,whatever this vague term may prove to mean precisely.5ConclusionsFor the normal carotid bifurcation at least,many of the pur-ported indicators of disturbedflow are significantly correlated. Based on thesefindings we recommend the use of relative resi-dence time͑RRT͒as a robust single metric of low and oscillatory shear.In light of possible redundancies,any perceived benefits of gradient-based HWP are likely outweighed by practical challenges with their measurement.Dominant harmonic͑DH͒is a promising new HWP,but issues related to its definition in nonaxialflow need to be resolved.AcknowledgmentThe authors thank Dr.Mort Friedman for helpful discussions. The authors also thank the anonymous reviewers for their valu-able suggestions.DAS acknowledges the support of Grant No. MOP-62934from the Canadian Institutes of Health Research. SWL and DAS were supported by,respectively,a postdoctoral fellowship and a career investigator award from the Heart and Stroke Foundation.References͓1͔Chien,S.,2008,“Effects of Disturbed Flow on Endothelial Cells,”Ann.Biomed.Eng.,36͑4͒,pp.554–562.͓2͔Himburg,H.A.,and Friedman,M.H.,2006,“Correspondence of Low Mean Shear and High Harmonic Content in the Porcine Iliac Arteries,”ASME J.Biomech.Eng.,128͑6͒,pp.852–856.͓3͔Caro,C.G.,Fitz-Gerald,J.M.,and Schroter,R.C.,1971,“Atheroma and Arterial Wall Shear.Observation,Correlation and Proposal of a Shear Depen-dent Mass Transfer Mechanism for Atherogenesis,”Proc.R.Soc.London,Ser.B,177͑1046͒,pp.109–159.͓4͔Fry,D.L.,1968,“Acute Vascular Endothelial Changes Associated With In-creased Blood Velocity Gradients,”Circ.Res.,22͑2͒,pp.165–197.͓5͔Friedman,M.H.,and Ehrlich,L.W.,1975,“Effect of Spatial Variations in Shear on Diffusion at the Wall of an Arterial Branch,”Circ.Res.,37͑4͒,pp.446–454.͓6͔Ku,D.N.,Giddens,D.P.,Zarins,C.K.,and Glagov,S.,1985,“Pulsatile Flow and Atherosclerosis in the Human Carotid Bifurcation.Positive Correlation Between Plaque Location and Low Oscillating Shear Stress,”Arteriosclerosis ͑Dallas͒,5͑3͒,pp.293–302.͓7͔Ojha,M.,1994,“Wall Shear Stress Temporal Gradient and Anastomotic Inti-mal Hyperplasia,”Circ.Res.,74͑6͒,pp.1227–1231.͓8͔Hyun,S.,Kleinstreuer,C.,and Archie,J.P.,Jr.,2000,“Hemodynamics Analy-ses of Arterial Expansions With Implications to Thrombosis and Restenosis,”Med.Eng.Phys.,22͑1͒,pp.13–27.͓9͔Gelfand,B.D.,Epstein,F.H.,and Blackman,B.R.,2006,“Spatial and Spec-tral1heterogeneity of Time-Varying Shear Stress Profiles in the Carotid Bifur-cation by Phase-Contrast MRI,”J.Magn.Reson Imaging,24͑6͒,pp.1386–1392.͓10͔Lee,S.W.,Antiga,L.,Spence,J.D.,and Steinman,D.A.,2008,“Geometry of the Carotid Bifurcation Predicts Its Exposure to Disturbed Flow,”Stroke, 39͑8͒,pp.2341–2347.͓11͔Friedman,M.H.,and Deters,O.J.,1987,“Correlation Among Shear Rate Measures in Vascular Flows,”ASME J.Biomech.Eng.,109͑1͒,pp.25–26.͓12͔Thomas,J.B.,Antiga,L.,Che,S.L.,Milner,J.S.,Hangan-Steinman,D.A., Spence,J.D.,Rutt,B.K.,and Steinman,D.A.,2005,“Variation in the Carotid Bifurcation Geometry of Young Versus Older Adults:Implications for Geo-metric Risk of Atherosclerosis,”Stroke,36͑11͒,pp.2450–2456.͓13͔Ethier,C.R.,Steinman,D.A.,and Ojha,M.,1999,“Comparisons Between Computational Hemodynamics,Photochromic Dye Flow Visualization and Magnetic Resonance Velocimetry,”The Haemodynamics of Arterial Organs—Comparison of Computational Predictions With In Vivo and In Vitro Data,X.Y.Xu,and M.W.Collins,eds.,WIT,Southampton,UK,pp.131–183.͓14͔Minev,P.D.,and Ethier,C.R.,1999,“A Characteristic/Finite Element Algo-rithm for the3D Navier–Stokes Equations Using Unstructured Grids,”Com-put.Methods Appl.Mech.Eng.,178,pp.39–50.͓15͔Ethier,C.R.,Prakash,S.,Steinman,D.A.,Leask,R.L.,Couch,G.G.,and Ojha,M.,2000,“Steady Flow Separation Patterns in a45Degree Junction,”J.Fluid Mech.,411,pp.1–38.͓16͔Moyle,K.R.,Antiga,L.,and Steinman,D.A.,2006,“Inlet Conditions for Image-Based CFD Models of the Carotid Bifurcation:Is It Reasonable to Assume Fully Developed Flow?,”ASME J.Biomech.Eng.,128͑3͒,pp.371–379.͓17͔He,X.,and Ku,D.N.,1996,“Pulsatile Flow in the Human Left Coronary Artery Bifurcation:Average Conditions,”ASME J.Biomech.Eng.,118͑1͒, pp.74–82.͓18͔Himburg,H.A.,Grzybowski,D.M.,Hazel,A.L.,LaMack,J.A.,Li,X.M., and Friedman,M.H.,2004,“Spatial Comparison Between Wall Shear Stress Measures and Porcine Arterial Endothelial Permeability,”Am.J.Physiol.Heart Circ.Physiol.,286͑5͒,pp.H1916–1922.͓19͔Lei,M.,Kleinstreuer,C.,and Truskey,G.A.,1996,“A Focal Stress Gradient-Dependent Mass Transfer Mechanism for Atherogenesis in Branching Arter-ies,”Med.Eng.Phys.,18͑4͒,pp.326–332.͓20͔Longest,P.W.,and Kleinstreuer,C.,2000,“Computational Haemodynamics Analysis and Comparison Study of Arterio-Venous Grafts,”J.Med.Eng.Tech-nol.,24͑3͒,pp.102–110.͓21͔Himburg,H.A.,Dowd,S.E.,and Friedman,M.H.,2007,“Frequency-Dependent Response of the Vascular Endothelium to Pulsatile Shear Stress,”Am.J.Physiol.Heart Circ.Physiol.,293͑1͒,pp.H645–H653.͓22͔Antiga,L.,and Steinman,D.A.,2004,“Robust and Objective Decomposition and Mapping of Bifurcating Vessels,”IEEE Trans.Med.Imaging,23͑6͒,pp.704–713.͓23͔Goubergrits,L.,Kertzscher,U.,Schoneberg,B.,Wellnhofer,E.,Petz,C.,and Hege,H.C.,2008,“CFD Analysis in an Anatomically Realistic Coronary Artery Model Based on Non-Invasive3D Imaging:Comparison of Magnetic Resonance Imaging With Computed Tomography,”Int.J.Card.Imaging, 24͑4͒,pp.411–421.͓24͔Glor,F.P.,Long,Q.,Hughes,A.D.,Augst,A.D.,Ariff,B.,Thom,S.A., Verdonck,P.R.,and Xu,X.Y.,2003,“Reproducibility Study of Magnetic Resonance Image-Based Computational Fluid Dynamics Prediction of Carotid Bifurcation Flow,”Ann.Biomed.Eng.,31͑2͒,pp.142–151.͓25͔Thomas,J.B.,Milner,J.S.,Rutt,B.K.,and Steinman,D.A.,2003,“Repro-ducibility of Image-Based Computational Fluid Dynamics Models of the Hu-man Carotid Bifurcation,”Ann.Biomed.Eng.,31͑2͒,pp.132–141.͓26͔Ford,M.D.,Alperin,N.,Lee,S.H.,Holdsworth,D.W.,and Steinman,D.A., 2005,“Characterization of V olumetric Flow Rate Waveforms in the Normal Internal Carotid and Vertebral Arteries,”Physiol.Meas,26͑4͒,pp.477–488.͓27͔Holdsworth,D.W.,Norley,C.J.,Frayne,R.,Steinman,D.A.,and Rutt,B.K., 1999,“Characterization of Common Carotid Artery Blood-Flow Waveforms in Normal Human Subjects,”Physiol.Meas,20͑3͒,pp.219–240.͓28͔Friedman,M.H.,2008,personal communication.͓29͔Lee,S.W.,and Steinman,D.A.,2007,“On the Relative Importance of Rhe-ology for Image-Based CFD Models of the Carotid Bifurcation,”ASME J.Biomech.Eng.,129͑2͒,pp.273–278.͓30͔Younis,H.F.,Kaazempur-Mofrad,M.R.,Chan,R.C.,Isasi,A.G.,Hinton,D.P.,Chau,A.H.,Kim,L.A.,and Kamm,R.D.,2004,“Hemodynamics and Wall Mechanics in Human Carotid Bifurcation and Its Consequences for Atherogenesis:Investigation of Inter-Individual Variation,”Biomech.Model.Mechanobiol.,3͑1͒,pp.17–32.͓31͔Huo,Y.,Guo,X.,and Kassab,G.S.,2008,“The Flow Field Along the Entire Length of Mouse Aorta and Primary Branches,”Ann.Biomed.Eng.,36͑5͒, pp.685–699.。