Analysis of Hemodynamics and Aneurysm Occlusion after Flow-Diverting Treatment in Rabbit Models

BACKGROUND AND PURPOSE: Predicting the outcome of flow diversion treatment of cerebral aneurysms remains challenging. Our aim was to investigate the relationship between hemodynamic conditions created immediately after flow diversion and subsequent occlusion of experimental aneurysms in rabbits. MATERIALS AND METHODS: The hemodynamic environment before and after flow-diversion treatment of elastase-induced aneurysms in 20 rabbits was modeled by using image-based computational fluid dynamics. Local aneurysm occlusion was quantified by using a voxelization technique on 3D images acquired 8 weeks after treatment. Global and local voxel-by-voxel hemodynamic variables were used to statistically compare aneurysm regions that later thrombosed to regions that remained patent. RESULTS: Six aneurysms remained patent at 8 weeks, while 14 were completely or nearly completely occluded. Patent aneurysms had statistically larger neck sizes (P = .0015) and smaller mean transit times (P = .02). The velocity, vorticity, and shear rate were approximately 2.8 times (P < .0001) larger in patent regions—that is, they had larger “flow activity” than regions that progressed to occlusion. Statistical models based on local hemodynamic variables were capable of predicting local occlusion with good precision (84% accuracy), especially away from the neck (92%–94%). Predictions near the neck were poorer (73% accuracy). CONCLUSIONS: These results suggests that the dominant healing mechanism of occlusion within the aneurysm dome is related to slow-flow-induced thrombosis, while near the neck, other processes could be at play simultaneously.

D espite the increased use of flow-diverting (FD) devices for treating intracranial aneurysms, the exact effects and processes responsible for the evolution of the aneurysm, and the ultimate outcome of these procedures are poorly understood. [1][2][3][4] Intrasaccular thrombosis and endothelial cell growth across the aneurysm neck have been observed after flow diversion and have been proposed as the main mechanisms driving the healing process. [5][6][7] However, the interaction, dominance, relative importance, and time precedence of these processes are not well-under-stood. Knowledge of these mechanisms is important to facilitate the development of future devices and therapies.
The purpose of this study was to investigate the relationship between the hemodynamic conditions created immediately after deployment of FD devices and the subsequent aneurysm occlusion in a rabbit model of saccular aneurysms.

Animal Models and Imaging
Twenty-three elastase-induced aneurysms were created in New Zealand white rabbits, following the approach described in Altes et al. 8 Four weeks after their creation, the aneurysms were treated with a FD device (Pipeline Embolization Device; Covidien, Irvine, California). Two days before treatment, the subjects were premedicated with aspirin (10 mg/kg by mouth) and clopidogrel (10 mg/kg by mouth) and were continued on these for 1 month after treatment. Immediately before treatment, 3D rotational angiography (3DRA) images were acquired, and velocities in the surrounding vessels were measured with Doppler sonography. Six animals were sacrificed before 1 week after treatment; all others, after 8 weeks. Immediately before sacrifice, 3DRA imaging was repeated. Some of the rabbits used in this study were part of another investigation in which we analyzed the mechanism of the endothelialization after flow diverter implantation. This article is entirely unrelated to the previous study.

Hemodynamics Modeling
Subject-specific computational fluid dynamics (CFD) models were constructed from pretreatment 3DRA images. 9 Unstructured grids were generated with a resolution of 0.2 mm. Models of the FD devices used to treat the aneurysms were created and virtually deployed within the reconstructed vascular models. 10 Blood flows were modeled by solving the unsteady 3D incompressible Navier-Stokes equations. 11 Physiologic flow conditions were derived from the Doppler sonography velocity measurements and imposed as boundary conditions in the computational models. Wall compliance was not included in the model. Blood attenuation was set to ϭ 1.0 g/cm 3 , and blood viscosity, to ϭ 0.04 poise. The governing equations were numerically solved by using an efficient finite-element solver with a time-step of 0.01 seconds. 12 After deployment of the FD device, the mesh was adaptively refined to resolve the stent wires and a new CFD simulation was performed by using immersed grids. 13,14

Occlusion Modeling
Regions of the aneurysm that remained open or patent and regions that were occluded at the time of sacrifice were identified by constructing a second vascular model from the 3DRA image acquired before sacrifice. The pretreatment model and this premortem or follow-up model were manually aligned by using rigid registration. A new grid was generated filling the volume of the follow-up model. The aneurysm neck was interactively delineated on the pretreatment model by connecting selected points along the path of minimum geodesic distance. 15 The aneurysm orifice defined by the delineated neck was triangulated, and grid points on each side of this surface were labeled as belonging to the aneurysm or the vessel. 15 The bounding box of the aneurysm was found and voxelized with isotropic voxels of the same resolution as the 3DRA images. Mesh points in the aneurysm region were then labeled as "open" or "occluded," depending on whether they were inside or outside the aligned follow-up model. The methodology is illustrated in Fig 1.

Global (per Aneurysm) Data Analysis
The degree of aneurysm occlusion at the time of sacrifice was quantified by the percentage of occluded voxels identified as described above. Aneurysms were then classified into 2 groups: 1) patent group: Ͻ80% of the voxels occluded at follow-up; 2) occluded group: Ͼ80% of the voxels occluded at follow-up. This 80% threshold was chosen because aneurysms with Ͼ80% occlusion of their volume looked almost completely closed on DSA images with perhaps a small "bump" in the parent artery-that is, they did not exhibit remnant necks of persistent filling of the sac.
Using the 3D meshes, we calculated geometric variables, including volume, maximum size, neck area, and maximum neck size for each aneurysm. Similarly, the following hemodynamic quantities were calculated over the aneurysm region and averaged over the cardiac cycle: mean aneurysm inflow rate ϽQϾ, mean aneurysm velocity ϽvϾ, mean aneurysm wall shear stress ϽϾ, and mean aneurysm shear rate Ͻ␥Ͼ. These quantities were computed before and after treatment. Additionally, a new variable denoted "mean aneurysm transit time" (MATT) was defined as

MATT ϭ ͑V/A͒ ϽvϾ
where V is the aneurysm volume, A is the neck area, and ϽvϾ is the mean aneurysm velocity. The rationale is that this variable, which has units of time, should do the following: 1) Increase with the aneurysm volume (ie, it takes longer to traverse a large aneurysm) 2) Decrease with the neck area (ie, a larger neck implies more chances of getting out of the aneurysm, while a smaller neck implies more chances of staying within the aneurysm and recirculating) 3) Decrease with the blood velocity (ie, slower velocities imply longer transit times).
The values of geometric and hemodynamic variables were then averaged over the patent and occluded groups. Cases belonging to the incompletely occluded group sacrificed before 1 week after treatment (between brackets in the Tables) were excluded from the statistical analysis. The Wilcoxon nonparametric test was used to test whether the means were statistically different between the occluded and patent groups.

Local (Voxel-by-Voxel) Data Analysis
Velocity (v), vorticity (), and shear rate (␥) were calculated on the CFD mesh and interpolated and averaged over each aneurysm voxel. The mean values of these quantities over open or occluded voxels of each aneurysm were calculated and compared.
Predictive statistical models of regional occlusion were then constructed. We defined 2 classes or outcomes: 1) occluded voxels (ϩ1), and 2) open voxels (Ϫ1). Seven attributes or features were considered: 1) ϽvϾ pre , 2) ϽϾ pre , 3) Ͻ␥Ͼ pre , 4) ϽvϾ post , 5) ϽϾ post , 6) Ͻ␥Ͼ post , and 7) d ϭ path distance to the neck. These values were calculated for each voxel from the pretreatment and posttreatment CFD simulations. A logistic regression model 16 was trained by using a "leave-one-out" methodology. Specifically, it was trained with all open and occluded voxels from all aneurysms except 1. The aneurysm left out was then used to test the accuracy of the model. The probabilities of each voxel of this aneurysm belonging to each class were calculated, and the voxel was "predicted" to belong to the class of the highest probability. The accuracy of the prediction was calculated by counting the number of correctly and incorrectly classified voxels. The process was repeated by leaving each aneurysm out one at a time, and the total predictive accuracy of the model was computed.
Next, the voxels of each aneurysm were divided into 3 groups according to their distance to the aneurysm neck. Denoting by d the distance from a voxel to the neck and d max the maximum distance from any voxel to the neck, we assigned voxels to the following 3 groups: The predictive statistical analysis described above by using the logistic regression model was repeated for voxels in each of these 3 regions. The corresponding predictive accuracies were calculated and compared.
Finally, we computed and compared the predictive power of the following statistical models typically used in machine-learning studies 17 : 1) logistic regression, 2) neural network, and 3) support vector machine. Each of these models was trained with data from all aneurysms except 1 and was tested on the one left out. The process was repeated by leaving each aneurysm out in turn. The total predictive accuracy of each model was computed and compared.

Hemodynamics
Flow visualizations revealed that velocities within the aneurysm are substantially reduced after treatment, that flow structures may change and, in particular, become smoother and simpler (ie, less swirling), and that the location of the inflow stream may shift from the distal part of the neck to the proximal part. To illustrate these observations, we present 2 examples in Fig 2, corresponding to incomplete and complete occlusions. The placement of FD devices thus causes both quantitative and qualitative changes in the hemodynamic environment within the aneurysm.

Global (per Aneurysm) Characteristics
Geometric characteristics of the aneurysms, the parent arteries, and the FD devices are presented in Tables 1 and 2. Aneurysms have been ordered by percentage occlusion. The aneurysm vol-ume was not statistically different between the open and occluded aneurysm groups. The neck area (P ϭ .0006), maximum neck size (P ϭ .0015), and maximum aneurysm size (P ϭ .02) were significantly larger in the patent than in the occluded groups.
Results from the global hemodynamics analysis are presented in the On-line Table. Only the mean inflow rate and mean aneurysm transit time were significantly different between the patent and occluded groups. The inflow rates before (P ϭ .0015) and after treatment (P ϭ .0305) were approximately 2.8 times larger in the patent group. The mean aneurysm transit time after treatment (P ϭ .02) was approximately 3 times larger in occluded aneurysms than in patent aneurysms. The mean posttreatment shear rate achieved marginal significance (P ϭ .0757) and was approximately 1.6 times larger in the patent group. All other quantities were not statistically different between the open and occluded groups.
To assess the possibility of predicting whether the aneurysm will be completely occluded at follow-up, we performed a receiving operating curve analysis. The area under the curve was calculated for each geometric and hemodynamic variable (listed in Tables 1 and 2 and the On-line Table). The best predictor was neck area and size, with an accuracy of approximately 70%-73%. However, our sample included a small number of aneurysms. and it is likely that this accuracy could be improved with larger samples.

Regional (Voxel-by-Voxel) Characteristics
Results from the local hemodynamics analysis are presented in Fig  3. This figure shows the mean value of posttreatment hemodynamics variables (velocity, vorticity, and shear rate) averaged over the open (red) and occluded (green) voxels of each aneurysm. Aneurysms in the x-axis are ordered by percentage occlusion. Open voxels have larger posttreatment velocity, vorticity, and shear rates than occluded voxels. These differences are less pronounced in aneurysms that had large patent regions (toward the left of Fig 3). The differences seem to increase with the percentage of occluded voxels (to the right of Fig 3). Ratios of posttreatment hemodynamic variables in open voxels over occluded voxels are presented in Table 3 for aneurysms that were not completely occluded. These results indicate that velocity, vorticity, and shear rate are on average larger by a factor of approximately 2.6 in open regions than in occluded regions. The ratios, however, are larger in aneurysms that were nearly completely occluded than in incompletely occluded ones. The P values indicate that these differences are statistically significant. Table 4 presents the accuracy and number of correctly predicted voxels over the total number of voxels for logistic regression models created for the entire aneurysm (sac), the dome, the body, and the neck regions. In addition, a range of accuracy values is provided for each region. These ranges were computed by calculating the accuracy of all possible combinations of the 7 features used to build the models (255 models for each region). Overall, the logistic regression model is capable of predicting which voxels will be occluded with very good accuracy and that voxels near the dome are more easily predicted than voxels near the neck.

Predictive Models
A comparison of different statistical classifiers typically used in machine learning is presented in Table 5. Overall, all statistical models had similar accuracies. The best performance was obtained with the neural network. The cases that were most difficult to predict were the ones with patent regions that were approximately 60%-70% of the aneurysm volume. The predictive accuracy for aneurysms that were approximately 30%-46% occluded ranged from approximately 40% to 50% (data not shown).

DISCUSSION
Flow diverters cause a significant disruption of the aneurysm inflow stream, which induces both qualitative and quantitative changes in the hemodynamic environment within the aneurysm. Aneurysms completely or almost completely occluded after treatment had statistically smaller necks than aneurysms that re-mained patent. They also had larger mean aneurysm transit times (related to the aneurysm residence time). Immediately after flow diverter implantation, the local velocity, vorticity (speed of rotation or swirling of fluid elements), and shear rate (rate of deformation of fluid elements) were all statistically smaller in regions later observed to occlude than in regions that remained patent. Statistical models based on local hemodynamic conditions were able to predict local occlusion with a very good accuracy. For a given aneurysm, regions close to the neck were more difficult to predict than regions that were close to the aneurysm dome. The mechanism of intrasaccular thrombosis and subsequent aneurysm occlusion is intimately related to the hemodynamic environment created after treatment. This environment depends mainly on the geometry of the aneurysm and parent artery and on the hemodynamic effects of the FD device. Aneurysm morphol-    ogy, in particular the location of the orifice on the parent artery and the neck size, determines the inflow stream into the aneurysm. The characteristics of the FD device, in turn, determine the disruption of this inflow stream and subsequent recirculation within the aneurysm sac. Thus, associations between both geometric and hemodynamic variables with outcome are expected. However, in this series, all aneurysms were created at the same location and were treated with similar FD devices. Thus, it is not surprising that the best predictor of complete occlusion was the neck size. Further studies with larger samples are needed to better understand the relationships among geometry, hemodynamics, and outcomes. Previous studies have largely focused on the changes of hemodynamic quantities from pre-to posttreatment. Pereira et al 18 used dynamic DSA images to estimate the mean aneurysm velocity (mean aneurysm flow amplitude) before and after FD treatment in a pilot clinical study (n ϭ 21), and found that at a later follow-up examination, aneurysms that were completely occluded (n ϭ 18) had larger reductions of the mean aneurysm flow amplitude than aneurysms still patent (n ϭ 3). Zhang et al 19 used CFD porous media to model the flow in 2 stented aneurysms, one remaining patent and the other completely occluded. They reported a larger reduction of aneurysm inflow velocity in the occluded than in the patent aneurysm. Huang et al 20 analyzed the flow dynamics in 14 rabbit aneurysm models treated with flow diverters. They used micro-CT images of harvested arteries with the flow diverter to create a model of the stent and place it in CFD models and found larger reductions in aneurysm inflow and larger increases of relative residence time in occluded aneurysms (n ϭ 10) than in patent aneurysms (n ϭ 4), though these differences did not reach statistical significance. They also observed changes in the location of the inflow stream, from the distal part of the neck to a more central/proximal position. In our study, we used CFD models with subject-specific flow conditions and FD models that reproduced the geometric characteristics of the implanted devices and compared global and local hemodynamic quantities in both occluded and patent aneurysms as well as between occluded and patent regions of each aneurysm. Our study focused on hemodynamic conditions created immediately after treatment rather than on changes from pre-to posttreatment.
The results of our study suggest that global hemodynamic conditions may not be sufficient to understand or predict aneurysm occlusion following flow-diversion treatment. They also suggest that thrombosis induced by the local hemodynamic conditions seems to be the main driver of aneurysm occlusion at the dome and body of the aneurysm, while other mechanisms such as endothelial tissue growth may be important near the neck. It is also possible that the thrombosis starts at the dome and progresses toward the neck, altering the intra-aneurysmal hemodynamic environment and inducing thrombus formation near the neck later followed by endothelialization. The regional analysis described in this work provides local information that could potentially be used to establish thresholds for thrombosis initiation and develop quantitative predictive models of local aneurysm occlusion.
The current study has a number of limitations that should be considered when interpreting the results. Unlike human cerebral arteries, flow reversals were observed in the parent arteries of our animal models. This could affect the values of the time-averaged hemodynamic quantities investigated. Alignment of follow-up and pretreatment vascular models was performed manually. This could result in inaccuracies of the quantification of the patent and occluded regions of the aneurysm. Our analysis did not distinguish between voxels that were close to the wall or near the middle of the aneurysm, which may have different hemodynamic environments. Likewise, the changes in the hemodynamic environment as the aneurysm progressively thromboses were not taken into account, only the hemodynamic conditions immediately posttreatment were studied. More refined models incorporating thrombus formation and proximity to the wall should be considered in future studies. The virtual deployment, the estimation of flow conditions from Doppler sonography images, as well as the many assumptions and approximations made during the CFD modeling process may affect the hemodynamic results. However, previous studies have indicated that these models could reproduce in vivo velocity measurements and observed flow characteristics. 21 Finally, the number of animal models included in the study is limited. Studies with a larger number of subjects are needed to confirm the trends observed in this work and to achieve stronger statistical significance.
Despite these limitations, this study confirms that local hemodynamic conditions created immediately after flow diverter implantation play an important role in determining whether different regions of the aneurysm will thrombose and occlude or remain patent. The next set of studies should focus on the connection between hemodynamic conditions and biologic processes responsible for not only thrombosis and aneurysm occlusion but also endothelial tissue growth and parent artery reconstruction.

CONCLUSIONS
Aneurysms completely occluded at follow-up had statistically smaller necks and longer mean transit times than aneurysms that remained patent. Posttreatment velocity, vorticity, and shear rate were significantly larger in regions that remained patent (ie, had larger "flow activity") than in regions that occluded. Machinelearning models based on local hemodynamic variables were capable of predicting local occlusion with good precision, especially away from the neck. This suggests that the dominant healing mechanism at the dome and body of the aneurysm is related to flow-induced thrombosis, while near the neck, other processes such as flow-modulated endothelialization may also play a fundamental role.