Intraluminal thrombus (ILT) is present in over 75% of all abdominal aortic aneurysms (AAAs) and probably contributes to the complex biomechanics and pathobiology of these lesions. A reliable predictor of thrombus formation in enlarging lesions could thereby aid clinicians in treatment planning. The primary goal of this work was to identify a new phenomenological metric having clinical utility that is motivated by the hypothesis that two basic haemodynamic features must coincide spatially and temporally to promote the formation of a thrombus on an intact endothelium—platelets must be activated within a shear flow and then be presented to a susceptible endothelium. Towards this end, we propose a new thrombus formation potential (TFP) that combines information on the flow-induced shear history experienced by blood-borne particles that come in close proximity to the endothelium with information on both the time-averaged wall shear stress (WSS) and the oscillatory shear index (OSI) that locally affect the endothelial mechanobiology. To illustrate the possible utility of this new metric, we show computational results for 10 carotid arteries from five patients where regions of low WSS and high OSI tend not to be presented with activated platelets (i.e. they have a low TFP), consistent with the thrombo-resistance of the healthy carotid despite its complex haemodynamics. Conversely, we show results for three patients that high TFP co-localizes with regions of observed thin thrombus in AAAs, which contrasts with findings of low TFP for the abdominal aorta of three healthy subjects. We submit that these promising results suggest the need for further consideration of the TFP, or a similar combined metric, as a potentially useful clinical predictor of the possible formation of ILT in AAAs.
Abdominal aortic aneurysms (AAAs) are responsible for significant morbidity and mortality, particularly in older males. Over 75% of these lesions contain an intraluminal thrombus (ILT), which, in many cases, fills the lesion and yields a lumen of nearly normal size despite the extreme dilatation of the aortic wall. There remains a significant controversy over the precise roles an ILT may play in the enlargement and potential rupture of these aneurysms, but it appears that thrombus is important both biomechanically and biochemically [1–3]. There is, therefore, a pressing need to understand better the conditions that lead to the formation of an ILT within these potentially life-threatening lesions.
The natural history of an ILT can be thought to consist of three basic phases: an initial change in platelet activity, the formation of an associated insoluble fibrin clot and possible fibrinolysis. The first phase actually involves three steps as well: the activation, adhesion and aggregation of platelets. Many prior studies have focused on specific aspects of the mechanisms underlying thrombus formation within a general haemodynamic field as revealed by multiple papers [4–6]. In particular, there has been an appropriate move towards multiscale models wherein one can capture contributions ranging from molecular-level chemical reaction kinetics to macroscopic haemodynamics. Notwithstanding these many advances, Cito et al. [6, p. 122] correctly conclude that ‘there is not any validated tool capable of informing clinicians in terms of predictive diagnostics and surgical planning of diseases wherein blood clotting plays a relevant role, such as aneurysm’.
In this paper, we adopted a different, complementary, approach to understand conditions that can lead to thrombus development in AAAs. Rather than including complexities of the molecular mechanisms, we sought a phenomenological metric that could reveal a haemodynamic threshold that stratifies the risk of forming a thrombus within a complex unsteady flow field. Moreover, different from most prior studies that focus on stenoses, ruptured atherosclerotic plaques or implanted devices wherein thrombus tends to form, we focused first on a region of the normal vasculature that is not susceptible to thrombus formation despite its complex geometry and haemodynamics. That is, we suggest that one can unveil lower bounds of useful haemodynamic metrics by identifying values for which thrombus does not occur and then contrast these values with those associated with clinical cases wherein thrombus has been observed. Towards this end, we focused first on the non-thrombogenic human carotid artery, including the bifurcation and the sinus, which is characterized by a local dilatation of the wall, and then examined three counter-examples based on AAAs harbouring a thin ILT. To complete the study, we also examined three models of healthy abdominal aorta free of thrombus. From these three cases, we found a new ‘thrombus formation potential’ (TFP) that appears to predict the location of thrombus formation within AAAs while explaining why normal regions of the arterial tree typically remain thrombus free.
2. Material and methods
Retrospective de-identified computed tomography (CT) or magnetic resonance (MR) images were collected for 11 patients, five to study the paired carotid arteries, three to study non-thrombotic infrarenal aortas and three to study AAAs harbouring a hint of ILT.
(a) Model construction: carotid arteries
Luminal geometries of 10 carotid arteries (two per subject) were extracted from four CT and one MR angiography scans (table 1) of patients under observation for conditions not directly related to carotid disease (average common carotid radius approx. 3.09±0.20 mm). Modest wall calcifications were found in the right carotid artery of Patient C2, but all arteries were free of stenosis. All patients presented a low degree of tortuosity in the carotid bifurcation except Patient C3, but thrombus was absent in all five patients.
The five pairs of images were first processed using a three-dimensional level-set segmentation algorithm available in the open source Vascular Modeling Toolkit (VMTK) . In order to minimize the impact of boundary effects in the haemodynamic simulations, we consistently included extended segments of the common and internal carotid arteries; the external carotid artery was segmented up to its first branch. VMTKwas also used to measure select geometric features such as the radius of the common carotid artery, defined as the average of multiple measurements taken along the vessel centreline, and the bifurcation angle between the internal and external carotids . Subsequently, the processed level-set images were imported into a customized version of the open source code SimVascular , where two-dimensional segmentations were extracted at several cross sections via thresholding. Non-uniform rational B-spline (NURB) surfaces were then lofted for each vessel through the segmented contours and finally merged into a complete geometric model (figure 1). A blending procedure was then used to smooth the junctions between the internal and external carotid arteries and thereby to minimize segmentation artefacts in the simulations.
The geometric model was then meshed using tetrahedral elements of nominal size 0.3 mm using the MeshSim library (Simmetrix Inc., Clifton Park, NY, USA) included in SimVascular. To capture the steep gradients in velocity close to the arterial wall, which was modelled as rigid, we meshed such regions at an increased resolution using four boundary layers of gradually decreasing thickness. Curvature-based refinement was also used in regions of high curvature such as the junctions between the internal and external carotids. On average, the meshes included about 1.7 million elements (1 724 027±307 495 elements), with variations attributable to differences in vessel sizes and curvatures.
Given the similar sizes of the 10 carotid models and the lack of patient-specific measurements of blood flow, we applied the same inflow boundary conditions to all models to simulate both rest and moderate exercise conditions. Specifically, the inlet flow waveform reported by Marshall et al.  was scaled according to mean flow rates collected by Sato et al.  on 10 healthy patients under rest and moderate exercise (40% of peak oxygen uptake) conditions. The Womersley solution for pulsatile flow in rigid wall cylinders  was then used to prescribe the inlet profiles [13,14]. To mitigate uncertainty on the assumed heart rates, again due to the lack of patient-specific data, we constructed a sequence of four cardiac cycles whose durations varied by 0.1 s about the mean values reported by Sato et al. . The sequences of four variable consecutive cycles (0.9–1.0–0.9–0.8 s and 0.6–0.7–0.6–0.5 s for rest and stress conditions, respectively) was iterated twice, thus a total of eight cardiac cycles were simulated for each geometry at rest and moderate exercise. As a comparison, simulations were also performed over eight cycles at a constant duration of 0.9 and 0.6 s, respectively.
At the outlets, we applied a three-element Windkessel model to replicate resistance (R) and capacitance (C) effects of the distal vasculature that were not included directly in the three-dimensional domain . The total resistance was first determined to match the mean arterial pressures reported by Sato et al.  for rest and exercise (87 and 103 mm Hg respectively), and then distributed among internal and external carotid arteries to reproduce the flow splits reported in the same source (QECA/QICA=0.54 and 0.58 at rest and exercise, respectively). Similarly, the value for total capacitance was taken from the literature  and then distributed between the ECA and ICA using the inverse of the coefficient used to split the resistance (table 2).
(b) Model construction: healthy infrarenal aortas and abdominal aortic aneurysms
The luminal geometries of the abdominal aorta and its main branches were extracted from three image datasets of healthy patients and three abdominal CT angiography scans of patients with AAAs partially covered by a thin layer (less than 3 mm) of ILT. The procedure was similar to the segmentation of the carotid bifurcations (figure 1). Briefly, level-set images obtained through three-dimensional segmentation in VMTK were imported into SimVascular to build a NURBS computational model. The final models included part of the suprarenal aorta, healthy or aneurysmal infrarenal aorta, superior mesenteric artery, celiac trunk, and renal and iliac arteries. The computational domains were discretized into tetrahedral elements of 1.0 mm maximum edge length. We used three boundary layers to resolve the potentially sharp gradients in velocity in regions next to the wall, which was again modelled as rigid. Curvature-adaptive meshing was used to reduce the size of the elements in smaller vessels and at bifurcations. The average number of elements for these models was 5 949 585±346 503.
The boundary condition at the supraceliac inlet was prescribed by imposing the averaged volumetric flow rate and waveform obtained by Les et al.  from PC-MRI measurements on 36 patients affected by AAA. Stress conditions were simulated using the waveform and flow rate reported by Suh et al.  for mild exercise (table 2). Owing to the lack of data on flow splits in the abdominal region for the degree of exercise considered here, the Windkessel parameters at rest were kept unchanged for exercise conditions. Similar to the analyses for the carotids, simulations were performed for eight cardiac cycles of varying duration. The average cycle at rest and exercise was 0.9 and 0.6 s corresponding to heart rates of approximately 67 and 100 beats per minute, respectively. A Womersley velocity profile was prescribed at the inlet [17–19] and Windkessel RCR boundary conditions were imposed at all outlets using values reported in the full-body three-dimensional model of Xiao et al.  (table 2). To avoid instabilities in the numerical solution in situations of back flow, we weakly constrained the shape of the velocity profiles to be parabolic at the outlets of the iliac and superior mesenteric arteries using the augmented-Lagrangian technique developed by Kim et al. .
(c) Finite-element models
The time-dependent velocity and pressure fields were obtained by solving the incompressible Navier–Stokes equations with the stabilized finite-element method implemented in SimVascular [21,22]. Blood was modelled as a Newtonian fluid with constant viscosity (0.004 Pa s) and density (1060 kg m−3). All arterial walls were assumed to be rigid, a simplifying assumption dictated by the scarcity of information on regionally varying material properties of the carotid (common, internal and external carotids) and different abdominal arteries; if our proposed metric, TFP, were to be adopted clinically, such an assumption would invariably be made as well due to the general lack of patient-specific information on regional wall properties. As we discuss below, the rigid wall assumption actually leads to a conservative estimate of the TFP, which is also desirable clinically. Each cardiac cycle was subdivided in time steps of 0.0001 s. Snapshots of the computational results were saved every 50 time steps (i.e. every 0.005 s). To eliminate the effect of initial transients, only results collected during the last four cycles were used for post-processing and particle tracking. All the simulations were run on the Omega supercomputer of the Yale High Performance Computing facility, where we employed on average 512 cores.
(d) Haemodynamic wall parameters and endothelial cell activation
Considerable attention has been directed over the past few decades towards characterizing effects that different haemodynamic stimuli may have on the endothelial cell gene expression that leads to either pro-atherosclerotic or pro-thrombogenic phenotypes. For example, several indices based on wall shear stress (WSS) and its temporal and spatial variations have been proposed to capture these mechanobiological effects with different degrees of success . We evaluated multiple indices (e.g. WSS, relative residence times, particle residence times near the wall) and combinations thereof, but ultimately combined two of the most commonly employed ones to better localize endothelial regions exposed to both low and oscillatory shear flows, conditions that promote a pro-thrombogenic phenotype [24,25]. The time-averaged magnitude of the WSS (TAWSS) is defined as 2.1where T is the total duration of the sequence of the last four simulated cardiac cycles (e.g. 3.6 and 2.4 s, respectively, for rest and stress conditions), and |τw| is the Euclidean norm of the WSS (actually in-plane components of the traction vector) on the luminal surface. To facilitate the comparison between the different simulations, and to reduce the effects of imposing the same boundary conditions for all models, nodal values of TAWSS (having units of Pa) were normalized by the average value of the most proximal segment of the domain (i.e. the common carotid artery in carotid simulations and the supraceliac portion of the aorta in infrarenal abdominal aorta (IAA) and AAA simulations), then denoted as .
The oscillatory shear index (OSI) has been proved successful in identifying atheroprone regions of the vasculature [26,27]. It is defined as follows: 2.2This dimensionless scalar index assumes its maximum value (0.5) within regions where the WSS field changes directions due to complex flow patterns. For example, a purely oscillatory flow with equal forward and backward contributions will produce an OSI of 0.5; conversely, in unidirectional flows in regular geometries, the OSI will be identically zero.
In an attempt to localize regions of the wall exposed to both high OSI and low TAWSS, we propose using the ratio of these two indices to characterize the degree of ‘thrombogenic susceptibility’ of the vessel wall. We refer to this metric as the endothelial cell activation potential (ECAP), namely 2.3Higher values of the ECAP index will thereby correspond to situations of large OSI and small TAWSS, that is, of endothelial susceptibility.
(e) Particle tracking and platelet activation
Computational fluid dynamics (CFD) studies typically produce values of velocity and pressure at every node of a fixed Eulerian computational domain. Such techniques do not provide, however, any information on the loading history that a fluid particle experiences along its path. Because the mechanisms that trigger platelet activation can involve prolonged mechanical stimulation [28–30], a Lagrangian interpretation of the flow field can help discern which endothelial regions are more likely to receive flow rich in activated platelets . Hence, we implemented a particle tracking procedure to enrich our CFD results with information on the shear history of all particles throughout the flow domain. Following a strategy similar to the investigations of Lagrangian coherent structures in the vasculature [19,32], we first re-meshed the computational domain to enhance the resolution of the particle tracking analysis (nominal edge size of 0.15 and 0.5 mm for the carotid and the abdominal aorta simulations, respectively). The nodal coordinates of the newly refined meshes were used as initial particle injection sites in an in-house parallel particle tracking code built on top of VTK and MPI open source libraries. Specifically, the code integrated the Lagrangian trajectories of fluid elements from the CFD velocity field using a fourth-order Runge–Kutta scheme. VTK classes (e.g. vtkCellLocator), endowed with methods for linear basis interpolation on finite elements, were then used to extract the velocity vectors at arbitrary locations within the domain. Owing to the discrete size of the integration time step, a small percentage of particles could erroneously be lost through the wall boundaries. Following the procedure reported in Duvernois et al. , we added a small inward velocity component of magnitude 1 cm s−1 to the otherwise zero velocity vector at the wall nodes to limit the leakage of particles without significantly affecting the tracking.
For each fluid element and at each integration time step, we collected the values of the velocity gradient and thus flow-induced shear encountered along a particle's path. Such information was used to compute a PLatelet Activation Potential (PLAP) recently proposed by Shadden & Hendabadi  to investigate the flow through an arterial stenosis. Briefly, the PLAP is a non-dimensional scalar index that represents the magnitude of shear rates that a fluid particle accumulates while travelling throughout the fluid domain. In our analyses, we tracked particles by looping twice over the four cardiac cycle sequence prescribed at the inlet. PLAP is defined as 2.4where |D(x(τ),τ)| is the Frobenius norm of the symmetric part of the spatial gradient of velocity tensor, t is the time of injection of the particle and 2T indicates how long the particle has been tracked (i.e. eight cardiac cycles in our case). Collecting particle information for multiple cardiac cycles allows one to capture flow stagnation events that might be of importance in thrombogenesis. If a particle were to exit the domain before completing the eight cycles, integration of its trajectory was interrupted and its PLAP congealed at the value computed at the last time step before exiting the domain. Although the PLAP was calculated for all particles in the fluid domain, our primary interest was in the trajectories of those particles that at some time came close to the endothelial layer. Hence, we tracked all particles backwards in time and PLAP was thereby computed on fluid elements whose positions at the end of the CFD simulations were known (and thus used as initial locations for particle tracking). After averaging results obtained for different injection times, we assigned final PLAP values to their corresponding nodes for further post-processing and plotting. Given that we were primarily interested in investigating the presentation of activated platelets to the endothelium, we probed the PLAP field at each boundary node at 20 equispaced points lying on a segment starting at the wall node and pointing towards the inward normal direction. Each vessel wall node was then assigned the averaged segmental value of the PLAP. To take into account potential effects of the probing distance from the wall, we computed PLAP for segments of length equal to 5% and 10% of the local vessel radius. Lastly, similar to what done for TAWSS, PLAP values were also normalized individually for each CFD model by a nominal value obtained by averaging PLAP on the most proximal segment of the computational domain (common carotid artery in carotid bifurcation analyses, and supraceliac aorta in aortic simulations).
(f) Thrombus formation potential
Of the numerous combinations of metrics considered for WSS and flow-induced platelet activation, a simple multiplicative index, the TFP, was finally defined to combine on the one hand the WSS-based ECAP obtained via standard CFD computations and on the other hand the fluid shear history-based PLAP obtained by particle tracking. The main purpose of this new index was to identify, if present, local regions of the wall that at the same time were exposed to pro-thrombotic WSS stimuli and a flow rich in activated platelets. We defined TFP as 2.5
(a) Non-thrombogenic carotid arteries
The ECAP highlights those regions where WSS is characterized by both a low intensity and pronounced changes in direction during a cardiac cycle. As expected, we observed high values of ECAP near the carotid bifurcation, where the flow was disturbed further by the sinus (i.e. a local dilatation). Conversely, the proximal segment of the common carotid and the distal segments of the internal and external carotids exhibited values of ECAP close to zero, confirming that such regions are not exposed to disturbed flows. Figure 2 shows distributions of ECAP in all 10 carotid arteries, both at rest and under mild exercise conditions. In all vessels save one, regions of high ECAP were on the outer wall of the bifurcation; the highest values of ECAP were on the inner wall in C2R (patient C2, right carotid), which corresponded to a site of calcification in that particular artery.
In all cases, prescribing exercise conditions did not cause evident changes in terms of the magnitude or distribution of ECAP. Despite an increase in WSS under exercise, index normalization with respect to the proximal value of TAWSS was sufficient to maintain peak values of ECAP similar. Yet, exercise conditions led to a slight broadening of the region of high ECAP towards the internal carotid and along the direction of flow in the bifurcation of C5L. Interestingly, this vessel was characterized by the largest bifurcation angle (81°), consistent with the importance that geometric features have in determining the peculiar haemodynamics in carotids . In particular, large bifurcation angles shield some regions of the proximal internal carotid from high shear flows, which could give rise to complex haemodynamics in exercise conditions as well. Table 3 shows the 99th percentile values of the ECAP index for all 10 carotid simulations. We believe that analysing upper percentiles might be more convenient and robust than using peak nodal values in discerning those cases of major interest where ECAP (or any other index) maintains high values at least in a small non-negligible region rather than only at isolated computational nodes. Among the carotids, the largest 99th percentiles were observed in C1R (upper ECAP percentile=2.70 and 2.47 under rest and exercise conditions, respectively) and C2R (upper ECAP percentile=1.79 and 2.24 under rest and exercise conditions, respectively). Carotid C1R also presented the second largest nodal value (max ECAP=9.15 at rest), while the largest peak was observed in carotid C4L (max ECAP=11.4 at rest).
Tracking particle trajectories backwards, rather than forwards, in time provided several advantages in terms of analysis and visualization. By plotting results for PLAP at the initial (injection) position of the particles, one can easily locate regions of the luminal surface that were exposed to either high- or low-shear fluid elements. Moreover, the procedure to obtain a near wall PLAP was greatly simplified because we could use the same particle tracking simulation to probe and average the PLAP field at different distances from the wall boundaries. Figure 3 shows distributions of PLAP in all 10 carotid models averaged over a near wall distance equal to 5% of the local radius. Of note, general patterns could be observed in all cases: PLAP was high throughout the straight part of the common carotid artery (particularly on the external side), but lower in the vicinity of the bifurcation. Gradually higher values were also found when moving distally along the internal carotid artery, whereas the portion of the external carotid included in the computational domain was probably too short to observe a similar behaviour. These results suggest that large arterial segments (especially regions with relatively simple haemodynamics like the common carotid artery) might naturally be in contact with fluid elements that have encountered high shear flows along their trajectory, and thus may potentially be rich in activated platelets. Wall regions subjected to more complex haemodynamics, like the proximal bifurcation and the sinus, instead seem to be relatively spared and could present a lower PLAP. Particles injected in the vicinity or downstream of these locations were, in fact, advected back through the complex flow that characterized the bifurcation; hence, their distribution of PLAP values was likely dictated by other haemodynamic features (e.g. presence of vortices, recirculation events and stagnation).
Simulations of mild exercise did not significantly affect the magnitudes or distributions of PLAP. The higher velocity gradients (and thus higher shear-induced activation potentials) that could be expected in exercise due to increased flow rates were counterbalanced by the lower times that particles resided within the domain because of the higher heart rates and overall flow velocities. In general, exercise conditions flattened the PLAP distributions, with a slight reduction in the peak values and localized increases in PLAP in some of the low-index regions, as, for example, in C2R (on the interior side) and C4L (on the outer wall) proximal to the bifurcation. Index scaling was effective in removing a weak correlation that we initially observed between peak values of PLAP and common carotid radii. Prior to normalization, the maximum value of PLAP was observed in C5L (radius=3.35 mm, second largest of the 10 common carotids), while the minimum value was observed in C4L (radius=2.76 mm, the smallest of the 10). As shown in table 3, scaled 99th percentile values of PLAP were instead very similar for all simulations (average 99th percentile of PLAP=2.13±0.15 and 2.02±0.14 for a 5% probing distance under rest and exercise conditions, respectively).
Table 3 also shows effects of the probing distance on the upper percentile values of PLAP. While differences due to probing over a larger distance could be considered to be negligible on average (3% reduction and 1% increase under rest and exercise conditions, respectively), the effects were more apparent on a case to case basis. The most pronounced differences were observed in the smallest carotids, C1R (22% decrease from a 5% to 10% probing distance both at rest and exercise) and C2R (16% and 30% increase under rest and exercise conditions, respectively). The averaging distance did not change the regional distribution of PLAP, however.
Figure 4 shows computed distributions of TFP (recall equation (2.5)) for all 10 carotid bifurcations at rest and exercise conditions. Consistent with the above, the 99th percentile values of TFP (cf. table 3) were in all cases smaller than at least one of the corresponding values of the ECAP and PLAP indices, and in most cases (e.g. for all measurements taken with a 5% probing distance of the local radius) smaller than both despite the upper percentiles of PLAP always being larger than 1.75. The same observation remains valid when considering peak values at computational nodes. On average, the peak nodal values of TFP were 4.73±1.80 and 5.05±2.07 for a probing distance equal to 5% of the local radius under rest and exercise conditions, respectively, whereas the average peak nodal values of ECAP were 6.64±2.41 and 5.96±1.58 under rest and exercise conditions, respectively. Comparison of figures 2 and 3 highlights how regions of high ECAP consistently coincided with regions of low PLAP and vice versa, which thereby attenuated the values of TFP (figure 4) in the healthy carotid arteries. This finding might provide some new insight into the mechanisms that inherently protect regions of the normal vasculature that are subjected to complex haemodynamic stimuli, such as the carotid bifurcation, from de novo thrombus formation despite their otherwise high vulnerability for atherosclerosis.
Indeed, it was this observation that motivated our definition of the TFP index. ECAP and PLAP are distinct in simple geometries: if a high shear flow rich in activated platelets comes close enough to the wall (high PLAP), it increases the magnitude of the wall shear (i.e. traction vector) and hinders its directional oscillation (low ECAP). High ECAP and PLAP regions thus tend to not overlap in the normal vasculature. Figure 5 illustrates this concept further for the two carotid arteries that revealed, respectively, the minimum and maximum values of peak TFP (C3R and C5L). Locations of high ECAP and PLAP, where ‘high’ is defined as the upper 10th percentile on the local surface areas, were marked on the external surfaces but with little overlap between the two high index regions. The area of intersection between the two high index regions was less than 1% of the total area in both cases (0.7% and 0.2% for C1R and C3R, respectively). This figure also shows a detailed view of the carotid sinus, where values of TFP were the highest. Notably, regions of high TFP tended to coincide with regions where ECAP was high, mainly because ECAP was close to zero in most regions of the wall subjected to relatively simple haemodynamics (following the well-known behaviour of OSI, on which the definition of ECAP is based). Thus, the PLAP did not contribute to changing the boundaries of non-zero TFP regions (which are dictated primarily by ECAP); rather, it confined the locations of the areas of peak TFP. Note, also, that in both cases the peak value of TFP was outside the overlap area of high ECAP and PLAP regions.
In summary, despite complex geometries and blood flows due primarily to the bifurcation and the sinus, the human carotid artery tends not to harbour intraluminal thrombosis in the absence of an eroded or ruptured atherosclerotic plaque [25,35]. Such was the case for the 10 non-atherosclerotic carotid arteries studied herein. The present simulations suggest further that, despite the complex haemodynamics, regions within the carotid artery that may be presented with platelets that have recently experienced a high shear history are also regions that experience higher values of non-oscillatory WSS. Increased WSS tends to upregulate endothelial production of nitric oxide and prostacyclin, both of which are anti-thrombotic, hence suggesting that regions of high (i.e. normal, not pathologically high) WSS should be protected. The postulated TFP was found to be similar under rest and exercise conditions, with the maximum 99th percentile value of 1.58 (for C1R at rest, max TFP=8.06) for an averaging distance of 5% of the local radius and 2.15 (for C2R under exercise conditions, max TFP=7.86) for an averaging distance of 10% of the local radius. Although higher values may be found for other normal human carotids, the present set of 10 results suggest that upper percentile threshold values for TFP should be above approximately 2.0 or 2.5, depending on the computational scheme. With this in mind, let us now consider three models of healthy abdominal aortas free of thrombus and three aortic aneurysms that harboured only a localized, thin ILT.
(b) Healthy abdominal aortas
The fluid dynamic analyses conducted on three healthy infrarenal aortas suggested that ECAP values observed in this region are, in general, lower than in the carotid sinus (figure 6). Despite observing high OSI values in large portions of the infrarenal aorta, as expected due to the backflow induced by the renal flow during diastole, TAWSS was also comparatively higher in the same region. Average 99th percentile values of ECAP were then approximately 56% and approximately 66% less than their corresponding values in carotids under rest and exercise conditions, respectively (table 3). Such differences were partially reflected in TFP percentile values, which were on average approximately 18% (under rest) and approximately 44% (under exercise) larger in carotids than in the infrarenal aortas, regardless of the probing distance employed. Lastly, exercise conditions (figure 6) seemed to have a moderate effect in reducing TFP in infrarenal aortas (approx. 20% reduction in average upper percentile TFP value under exercise conditions). These findings substantiated our original idea to study the carotid bifurcation as a region of highly disturbed flow that is yet normally thrombo-resistant. Such regions provide important insight into possible threshold values of candidate metrics.
(c) Abdominal aortic aneurysms with thin thrombi
Figure 7 shows the geometric model used to simulate the haemodynamics associated with an approximate 4.4 cm diameter abdominal aortic aneurysm (AAA1), including major arteries of the abdominal vasculature. The lesion dilatation was more pronounced on the anterior side, with the proximal neck located approximately 3.5 cm below the renal arteries and shifted forward with respect to the path of the healthy aorta. Shown, too, are full-field distributions of ECAP, PLAP and TFP together with a detailed lateral view of the lesion, with the actual thin ILT visualized via a transparent overlay on the central surface. Not surprisingly, our computational analyses on AAA1 predicted haemodynamic conditions that are typically expected in AAAs. After detaching from the wall immediately below the neck, the flow promptly decelerated giving rise to recirculation and mixing. Accordingly, the ECAP index (a measure of low and oscillatory WSS) was high throughout the upper and central portion of the anterior portion of the lesion, with an upper percentile value higher than that found in carotid bifurcations C1R and C2R (upper ECAP percentile=4.13 in AAA1 versus upper ECAP percentile=2.70 and 1.79 under rest conditions in C1R and C2R, respectively). Conversely, the peak nodal values were significantly smaller (max ECAP=5.54 in AAA1 versus max ECAP=9.15 and 11.4 under rest conditions in C1R and C2R, respectively).
The PLAP index was also higher in the anterior and central portions of the lesion, in contrast with that which was observed for the carotid arteries wherein regions of high PLAP were systematically characterized by low values of ECAP. This finding for the aneurysm is likely due to the more complex haemodynamics, where secondary vortices and prolonged recirculation might favour the mixing of (and thus transfer of particles between) low and high shear flows. Although peak nodal values of both ECAP and PLAP were higher in some of the carotid simulations, the maximum TFP was higher in AAA1 (max TFP=8.65 in AAA1 versus max TFP=8.06 in C1R at rest and for a probing distance equal to 5% of the local radius). Differences were more pronounced in terms of 99th percentiles, with the aneurysmal value being more than three times larger than the largest carotid value (upper TFP percentile=6.27 in AAA1 versus upper TFP percentile=1.58 in C1R under rest conditions and for a probing distance equal to 5% of the local radius). The region covered by the actual ILT encompassed areas of high TFP, including the region wherein PLAP was found to reach its peak value.
Figure 8 shows spatial distributions of TFP in two other aneurysms (AAA2 and AAA3) of approximately 5.1 cm and approximately 5.2 cm diameter, respectively. Despite their larger size, close to the clinically accepted criterion for intervention (5–5.5 cm), these lesions harboured very thin (and in the case of AAA3 even fragmented) ILT. Interestingly, both AAAs presented calcified proximal necks with peculiar morphologies that proved to affect significantly their haemodynamics. In AAA2, a slight stenosis at the inlet of the aneurysm caused a prompt acceleration of the entering blood flow. The high shear jet was directed towards the anterior wall of the lesion, locally increasing the WSS and thus explaining the very low TFP values observed throughout the anterior region. It was only after prolonged recirculation that the newly injected flow reached the posterior wall, gradually losing intensity from the bottom towards the top. Two high TFP islands were visible in the upper posterior wall, one of which coincided with the region covered by the actual ILT. The peak value of TFP was found within the other region of high TFP; it was approximately 70% greater than the maximum value observed in the carotids (max TFP=13.54 in AAA2 versus max TFP=8.06 in C1R for a probing distance equal to 5% of the radius and rest conditions) and similarly the 99th percentile was about twofold higher (99th upper TFP percentile=3.22 in AAA1 versus 99th upper TFP percentile=1.58 in C1R for a probing distance equal to 5% of the radius and rest conditions). Nevertheless, differences in both percentile and peak values of ECAP or PLAP were not evident between this aneurysm and the carotid arteries, thus suggesting that it was critical to combine the information provided by the ECAP and PLAP to discriminate between regions protected from or prone to the formation of ILT.
AAA3 also presented a stenotic segment preceding its inlet, within which flow gained enough momentum to maintain shear stresses high along the wall of the lesion despite a complex haemodynamics characterized by multiple recirculation vortices. Thin ILT was present in 4 separate regions located in the upper and posterior wall. An analysis of the distribution of the TFP over the entire lesion showed high values only in one of these locations (max TFP=8.33 at rest for a probing distance equal to 5% of the radius, slightly above the maximum value found in the carotid bifurcations), while the index was uniformly low throughout the rest of the wall. As a result, the upper TFP percentile of AAA3 was the lowest among the aneurysms and comparable with the maximum value observed in the carotid analyses (upper TFP percentile=1.46 for AAA3 versus upper TFP percentile=1.58 for C1R at rest and for a probing distance equal to 5% of the local radius). The fact that three out of four thrombus pockets seemed to be in regions not haemodynamically prone to thrombus deposition (mainly due to high WSS) would suggest that these ILTs might not be the foci of a new formation event, but rather the remnants of a previously deposited (and more uniform) thrombus that had been excavated by the high WSS flow characterizing the lesion. Indeed, this was found to be the case by examining a prior image. Hence, this peculiar aneurysm could be considered as a cautioning case, with TFP distributions that do not fully represent a novel thrombus deposition, but rather provide valuable lower bounds for the desired threshold value of the TFP.
Figure 9 summarizes results from table 3 graphically and reveals that the TFP alone delineated well between all 13 normal thrombo-resistant arteries and the three pathological thrombo-prone aneurysms. These results well support our original hypothesis that one must consider both the condition of the endothelium and the flow-induced shear history experienced by platelets presented to the endothelium.
The primary goal of this work was to motivate a new phenomenological approach for identifying regions of possible formation of ILT on an intact but susceptible endothelium within AAAs. That is, in contrast to prior haemodynamic studies (see reviews [4–6]) of thrombus formation in regions of arterial injury (e.g. eroded stenoses or ruptured atherosclerotic plaques) or within medical implants (e.g. a bioprosthetic valve or metallic device), we focused on regions where blood would not be exposed to collagen or other highly thrombogenic material. In particular, this study was motivated by the hypothesis that two basic haemodynamic features must coincide to promote the formation of a thrombus on an intact endothelium—platelets must experience a sufficiently high shear history and then be presented to a susceptible region. We submit that the proposed TFP—which combines three prior useful metrics of low WSS, high oscillatory WSS and high flow-induced particle shear history—captures these two key features (cf. figure 9).
Simulations for 10 non-diseased carotid arteries from five patients revealed that regions that are usually susceptible to low and oscillatory shear stress (and hence atheroprone and potentially pro-thrombotic) did not overlap with regions that are likely to be presented with platelets that have recently experienced a high shear history. Although explored in detail for only one location of complex flow within the human vasculature, this finding seems reasonable in hindsight—there must be a mechanistic reason why the healthy carotid artery, with its geometric complexity, is not susceptible to thrombus formation under physiologic conditions despite being susceptible to atherogenesis. There is a need nonetheless to evaluate this new index and its lower bound value (e.g. 99th percentile of TFP∼2.0 for a 5% averaging distance) for other locations, perhaps including the aortic arch as well as other bifurcations. The carotid artery was selected because of the uniqueness of its geometry, which includes both a bifurcation and a distal sinus that effectively represents a local dilatation of the wall, and the availability of prior studies on the carotid that provided reliable information on inlet and outlet boundary conditions [10,11]. As confirmed by our analyses, this choice proved to be particularly insightful because the endothelium within the carotid sinus is constantly subjected to low and oscillatory WSSs and at the same time contains flows potentially rich in activated platelets. The carotid bifurcation was thus a convenient test case for our hypothesis of combined mechanisms. Conversely, three analyses conducted on thrombus-free normal infrarenal aortas (figure 6) revealed that albeit oscillatory, aortic blood flow was normally not stagnant and maintained a relatively high shear stress on the vessel wall. Such factors led to very low values of TFP index, in agreement with the observation that healthy aortas similarly do not normally develop thrombus, and that the carotid artery was a better comparator for identifying a safe (lower) bound value of the TFP.
By contrast, we also sought patient-specific cases of AAAs that harboured but a nascent ILT. Lesions with large thrombus would obviously not be useful in searching for conditions under which the ILT initiated because the overall geometry of the lumen and wall likely change following the accumulation of thrombus. Conversely, lesions without thrombus could provide additional information on the lower bound value of the TFP, but would not help in identifying how far above this bound the ILT might initiate. We were thus fortunate to find three rare cases wherein there was a small, thin ILT within a sizable aneurysm (figures 7 and 8). That the average 99th percentile of TFP in these aneurysms (approx. 3.60) was more than five times the corresponding value found in carotids (approx. 0.70) and healthy aortas (approx. 0.58), as shown in table 3 and figure 9, and that in all cases peak values of TFP co-localized with the regions of the actual ILTs (which were removed from the lesion numerically to perform unbiased simulations) provided additional confidence that the TFP merits further study. There is a need, therefore, to identify and analyse many more AAAs wherein there is but a nascent thrombus.
Albeit not shown (see the electronic supplementary material, tables and figures), we also reperformed all simulations for constant rather than variable heart rates. Although values of TFP changed slightly in both the arteries and the aneurysms, the co-localization of regions of high TFP with the actual thin ILT in the aneurysms was even stronger. Hence, the results presented visually in figures 7 and 8 may be conservative, which is desirable of a potential clinical metric. We also performed additional simulations on aneurysm AAA2 after numerically removing the slight stenosis at the proximal neck (i.e. inlet). Doing so resulted in higher values of TFP in the anterior portion of the lesion, where ILT is often found in similarly shaped lesions. In this case, one could speculate that the presence of the stenosis and associated inlet ‘jet-like’ flow may have been somewhat protective from the perspective of ILT formation. If so, categorization of such effects on the haemodynamics and associated platelet activation/adhesion could become a useful clinical aid.
There have been many prior studies of possible haemodynamic initiators of thrombus within AAAs. For example, Bluestein et al.  presented results many years ago for steady flow simulations within idealized axisymmetric lesions. They concluded that ‘the recirculation zone formed inside the aneurysm cavity creates conditions that promote thrombus formation’ [36, p. 280]. In particular, they suggested that ‘Once trapped in the recirculation zone, platelets with elevated shear histories and higher incidence of activation were readily deposited to the wall in areas of low wall shear stress’ [36, p. 284]. Although these observations were extremely insightful, they did not suggest a simple metric to capture these multiple effects. The TFP proposed herein represents such a metric that is consistent with their seminal findings.
More recently, attention has been directed more towards the strain history experienced by the platelets, particularly in vortical structures. For example, Biasetti et al.  focused on high fluid shear stresses generated by vortical structures within fusiform aneurysms, which move within the body of the lesion and potentially release activated platelets upon the break-up of the vortices. Shadden & Hendabadi  similarly focused on elevated fluid shear stresses and exposure times experienced by platelets, albeit for thrombus formation in stenoses. The PLAP used herein was taken directly from the latter work. Neither of these studies attempted to tie together the activation of the platelets due to elevated shear histories and their more likely adhesion and aggregation on intact but susceptible endothelium in regions of low WSS or high OSI. Again, the TFP is consistent with these important studies.
We emphasize again that our goal was to identify a potentially useful clinical metric, not to explore underlying biochemomechanical mechanisms of thrombus formation. The latter is yet important and necessitates multiscale, multiphysics approaches as being pursued by others (cf. [4–6]). Nevertheless, a phenomenological metric could also be very important, particularly in treatment planning. Such a clinical metric necessarily should be based on information usually available in the monitoring of AAAs, namely, medical images that reveal the patient-specific geometry. That is, even though large-scale fluid–solid interaction models are now practical , information on regional variations in patient-specific material properties remain wanting. We thus assumed rigid walls in all of our simulations, as would be expected of most computations based on current clinically available data. Albeit not shown, we performed a series of computations to determine TFP in idealized abdominal aortas with or without axisymmetric aneurysms that were generated using our growth and remodelling codes . In this way, we could compare directly the results based on a rigid versus a compliant wall (i.e. one endowed with appropriate regional variations in biaxial material stiffness and wall thickness). Our findings revealed, consistent with the definition of ECAP, that the computed value of TFP was lower in the deformable wall simulations for the healthy infrarenal aorta because the WSS maintained higher values during late systole and diastole (despite an overall lower peak value), thus resulting in generally higher TAWSS values when taking into account wall elasticity. In other words, a TFP computed for a healthy artery based on a rigid wall assumption (cf. figure 9 for the carotids and infrarenal aorta) assumes a higher (i.e. conservative) value, which is desirable clinically. Interestingly, the computed value of TFP for idealized, deformable aneurysms was similar to that obtained via the rigid wall assumption, in part because of the extreme stiffening that occurs following the loss of elastin in aneurysms . Hence, the rigid wall assumption adopted here was both consistent with expected future simulations of clinical cases and computationally conservative as desired when comparing normal versus pathologic cases.
In conclusion, for complex ‘pathologic’ geometries and haemodynamics with pronounced stagnation and recirculation (e.g. as in aneurysms), particles that have previously accumulated a high shear history could potentially be advected to the wall even by a lower shear flow that would result in low WSSs on the endothelium. In this case, the computed TFP may provide much more information than simply a WSS-based parameter like ECAP, a measure such as (near wall) residence time, or a single measure of platelet shear history [31,37,40,41]. That is, the present findings motivate further consideration of the TFP, or a similarly combined index, that accounts for the potential of activated platelets being presented to a quiescent (non-thrombogenic) versus a susceptible (thrombogenic) endothelium. Being able to predict when a lesion might develop an ILT could enable another level of computational sophistication (cf. ) and thereby contribute to treatment planning for AAAs.
This study was approved by the Institutional Review Boards of Yale University and the Veterans Affairs Connecticut Healthcare System.
This work was supported, in part, by the staff and facilities of the Yale University Faculty of Arts and Sciences High Performance Computing Center as well as grants from the NIH (R01 HL086418, U01 HL116323).
- Received February 27, 2014.
- Accepted August 13, 2014.
- © 2014 The Author(s) Published by the Royal Society. All rights reserved.