Key influencers in an aneurysmal thrombosis model: A sensitivity analysis and validation study

Intracranial aneurysm (IA) is a type of cerebrovascular pathology, which is a localized dilation or ballooning of the cerebral blood vessel caused by the weakness of the wall of a cerebral artery or vein.¹ There are three main treatment options for patients with IAs: observation, surgical therapy, and endovascular therapy.² The goal of treating patients with unruptured IAs is to maximize their duration of high-quality life by optimally balancing the risks of aneurysm rupture with those of treatment-related adverse outcomes.³ In the literature, 24.2% (244/1009) IAs⁴ failed to obtain aneurysm occlusion after endovascular or surgical treatment; aneurysm reopening and retreatment after endovascular coiling occurred in 20.8% (1697/8161) and 10.3% (840/8161), respectively.⁵ Some aneurysms fail to develop a stable clot even with sufficient levels of flow reduction and may end up with post-treatment rupture, leading to high risks of mortality and morbidity.^6–8 From autopsy studies of aneurysms, researchers found that there are two different types of thrombi: organized white thrombus, rich in fibrin and platelet, and non-organized red thrombus, rich in fibrin and erythrocyte.^8–10 They can be found in stable clots and unstable clots, respectively. Red thrombi are the result of stagnation of blood flow, resulting in a clot containing all elements of normal blood, and they contain more enmeshed erythrocytes among sparse fibrin strands compared to precipitation or white thrombi. The red thrombi are expected to progress to organized white thrombi; otherwise, they may promote an inflammatory reaction, eventually leading to the disintegration of the aneurysm wall with subsequent rupture. Achieving organized white thrombi may reduce the probability of post-treatment rupture, and non-organized red thrombi have also been suggested as a potential predictor for unsatisfactory treatment results.^8–10 

Thrombosis is the process of the formation of a blood clot inside a blood vessel that obstructs the flow of blood through the circulatory system, which is different from the natural hemostasis process.¹¹ Under physiological conditions, the formation of a blood clot is a well-regulated process that includes three elements: (1) primary hemostasis, (2) secondary hemostasis/coagulation, and (3) fibrinolysis.¹² The key difference between natural hemostasis and intracranial aneurysm thrombosis relates to how to trigger primary hemostasis. In natural hemostasis, the onset is triggered by blood exposure in endothelial tissue caused by injury to the vascular wall, while in cerebral aneurysm thrombosis it has been linked to endothelial damage present in the aneurysm sac, wall inflammation, blood-borne tissue factor, and contact with artificial surfaces after treatment.^13–16 The formation of thrombi within the aneurysm is an important determinant of the outcome after endovascular therapy,¹⁷ and controlled thrombosis leading to stable clotting is also the main post-intervention goal of flow diverter treatments.^1,14,18 However, the mechanisms that govern the initiation and evolution of thrombosis are not fully understood¹³ as the regulation of thrombosis is extremely complicated.¹⁹ 

In unruptured aneurysms, thrombosis (spontaneous or device-induced) can stabilize the aneurysm or accelerate the path to rupture.¹³ Currently, in vivo or image-based analysis of thrombosis hemodynamics in realistic anatomies and physiologies is very difficult, if not impossible.²⁰ In recent decades, significant effort has been directed toward computational predictions of hemodynamics in aneurysms.¹³ However, computational prediction of thrombosis within aneurysms is relatively unexplored. The integrated thrombosis model originally developed by Sarrami-Foroushani et al.⁸ incorporates biochemical reactions, platelet activity, and hemodynamics. Briefly, the model combines platelet activation and transport with fibrin generation and defines a flow-induced platelet index (FiPi) as a quantitative measure of thrombus stability. According to aneurysm autopsy studies, two different types of thrombus have been identified: the unstable red thrombus (rich in fibrin and erythrocytes) and the stable white thrombus (rich in fibrin and platelets). FiPi quantifies the effect of blood flow on the transport of platelets to and from the site of thrombus formation, and thus on the final platelet content of the formed thrombus. FiPi is related to the initial concentration of the resting platelets, the initial concentration of the activated platelets, and the concentration of the bound platelets. During the thrombus formation process, resting platelets become activated by exposure to thrombin or other activated platelets. The activated platelets adhere to the fibrin network aggregate to form bound platelets. Both activated and bound platelets are derived from resting platelets, and thus FiPi is highly related to the resting platelet concentration. Sarrami-Foroushani et al.⁸ set FiPi >0.15 as a threshold for the formation of fibrin and platelet-rich white thrombi. This thrombosis model is not only capable of predicting both the hemodynamic changes and the thrombus formation process but is also able to predict long-term thrombus stability by investigating the thrombus composition. However, such a comprehensive computational thrombosis model that considers both hemodynamics and biochemical reactions is very complex with a large number of uncertain model parameters. Previously, we calibrated the hemodynamic thresholds, residence time (RT), and shear rate (SR), of thrombosis initiation against real population-specific data,²⁰ but how the rest of the model parameters affect the final formed thrombus has not yet been assessed. To have a robust computational thrombosis model for possible clinical use in the future, it is essential to assess the model reliability through comprehensive sensitivity analysis (SA) of the model parameters and validation studies based on clinical information from real patients. Uncertainty quantification of a computational model is crucial in in silico trials to ensure the accuracy and reliability of predictions (model credibility),²¹ thus improving confidence in regulatory submissions. It helps identify and manage potential risks, ensuring robust and credible simulation outcomes that can effectively replace or supplement traditional clinical trials.

This paper aims to first identify the most influential factors in our previously developed thrombosis model⁸ through a comprehensive global SA. We then validate the thrombosis model based on a real patient case (partial thrombosis before treatment and residual neck after immediate post-treatment) using patient-specific parameters for those identified as influential and two previously calibrated trigger thresholds²⁰ of thrombosis initiation. In addition, we improve our thrombosis modeling for untreated aneurysms by narrowing the thrombosis initiation in areas near the wall, as in the real situation, the thrombus is difficult to be suspended in an aneurysm lumen on its own without any anchors to the surrounding aneurysm wall. The novelty of this study is that we not only identify the most influential factors in the modeling of aneurysmal thrombosis but also demonstrate for the first time the ability of our thrombosis model in predicting the thrombus formation both before and after treatment based on a clinical case of the patient.

Using the reduced lumped parameter model, we conducted 14 000 simulations for −50% to +200% variation of 12 kinetic parameters (Table S1) using the thrombus formation duration as the output metric. As shown in Fig. S1, the SA results indicate that the duration of thrombus formation is sensitive to 4 kinetic parameters: $CFI,50, KFITH, KTHAP$⁠, and $KTHAT$⁠. Platelet recruitment and deposition were assumed to depend on the concentration of free platelets and the value of a second-order Hill function $ϕPBFI$ with $CFI,50$ = 60 nM. $KFITH$ is a kinetic constant related to thrombin-mediated fibrin generation. $KTHAP$ is the kinetic constant of the kinetic reaction of thrombin generation on the surface of activated platelets. Thrombin inhibition by antithrombin was modeled as a second-order reaction with kinetic constant, $KTHAT$⁠. The other 8 kinetic-related parameters have limited or negligible effects on the lumped parameter model.

We identified 4 key parameters from the SA results of the lumped parameter model. These 4 kinetic-related parameters and another 14 parameters (Table S1) were assessed based on a 3D aneurysm geometry (Fig. S2: a spontaneous thrombosis case from our previous study;²⁰ male; 51 years old; aneurysm size, 6.8 mm; aspect ratio, 1.5; location, PComA) using the 3D full model and EE method. For 18 parameters and 5 randomly generated paths, we conducted $r(k+1)=95$ model runs for the 3D full thrombosis model using the EE method. According to previous studies,^8,22,23 the clot was assumed to be formed in regions where the fibrin concentration is greater than 600 nM. As shown in Fig. S3, the space-averaged fibrin concentration in the aneurysm sac is converged at 100 cardiac cycles when using three different resting platelet concentration values. Therefore, we ran all 95 simulations for 100 cardiac cycles in this SA study using EE method. The output metric is the measure of thrombus composition (FiPi⁸ and the bound platelet concentration). The EE method allows us to classify the inputs into three groups: (1) small $μi*$⁠: inputs have negligible effects; (2) large $μi*$ and small σ_i: inputs having large linear effects without interactions; (3) large $μi*$ and large σ_i: inputs having large non-linear and/or interaction effects.

As shown in Figs. 1 and S4, the SA results show that the concentration of resting platelets has the greatest effect on the final formed thrombus composition (FiPi and bound platelet concentration), while the other 17 parameters have limited or negligible effects on the thrombosis model.

High RT and low SR are widely used in computational models to characterize flow stasis that triggers the thrombus formation process in aneurysms. In our previous in silico observational study,²⁰ we calibrated these trigger thresholds as RT 1.9 s and SR 11 s⁻¹. Using these calibrated trigger thresholds, we ran the flow simulation model (about 10 hours per run using 128 cores) to obtain hemodynamics inside the aneurysm sac to predict the possible thrombosis region. The flow simulation with all the reaction terms switched off is an initialization simulation for the following coupled flow and thrombosis simulation. From hemodynamics (Fig. 2), this patient-specific aneurysm is shown to be a spontaneous thrombosis case, and the main possible thrombosis region is located in the left side of the aneurysm sac, which coincides with the clinical ground truth before treatment.

We then performed the coupled flow and thrombosis model (about 10 days per run using 128 cores) to investigate how the selection of patient-specific and non-patient-specific parameters affects the final formed thrombi (Table I). The settings of non-patient-specific parameters were from Sarrami-Foroushani et al.,⁸ while the settings of patient-specific parameters were from the patient's clinical record. As shown in Fig. 3, the red part (located in the lower-left corner of the aneurysm sac) in the clinical ground truth figure is the manually labeled partial thrombosis region. There is no significant difference between the patient-specific simulation results and the non-patient-specific simulation results in terms of the location of the thrombi and the shape of the thrombi. However, a higher percentage of stable thrombi (44.6%) was obtained in the patient-specific model compared to the non-patient-specific model (40.6%). The patient-specific and non-patient-specific models predicted that the thrombosis region was located in the middle and left sides of the aneurysm sac and the thrombus grew mainly from the middle of the lumen to the aneurysm wall, which is unrealistic, as the thrombus is difficult to be suspended in the middle of the sac without any anchors. To address this model limitation, we constrained the thrombosis initiation to occur only near the wall or other thrombosed regions. As shown in the third row of Fig. 3, the thrombus grew mainly near the left wall, rather than in the center of the sac, which better matches the clinical ground truth prior to treatment. Although the wall-constrained initiation model successfully predicted the primary thrombosis initiation site before treatment, the extent of thrombosis was significantly underestimated compared to the clinical ground truth. As introducted in Sec. IV, the comparison with clinical observations is complicated by the unstable and unconverged state of the thrombus. Nevertheless, successfully predicting the main thrombosis initiation site is considered a significant achievement.

As shown in Fig. S5, we virtually deployed patient-specific coils and a flow diverter using GIMIAS (version 1.8.r1).²⁴ We then performed a post-treatment simulation with patient-specific virtual treatments and patient-specific concentrations of resting platelets and fibrinogen for 200 cardiac cycles of simulation time using CFX's adaptive time-stepping with minimum, maximum, and initial time steps of 0.0001, 0.05, and 0.01 s. It took about 2 months to obtain the results of the post-treatment simulation using 128 cores. The volume of the formed thrombi in the aneurysm sac over time is presented in Fig. S6. The thrombus grows very slowly, and it is noteworthy that the thrombus volume at 175 cardiac cycles reaches 95% of the final volume observed at 325 cycles.

According to the O'Kelly–Marotta^25,26 (OKM) grading scale, aneurysms are assigned grades based on the amount of contrast filling of the aneurysm sac (filling grades, A, B, C, and D) and how long contrast persists in the aneurysm sac with respect to angiographic phase (stasis grades, 1, 2, and 3). For this patient-specific aneurysm case treated with a flow diverter and ten coils, the aneurysm incompletely filled its lumen with contrast that persists within the lumen into the capillary phase of the angiogram [Fig. 4(a)]. It was assigned grade 2B according to the clinical record. As shown in Fig. 4, our simulation results show good agreement with the clinical immediate post-treatment angiographic result: minimal residual flow in the neck of the aneurysm after treatment.

Thrombosis is a biological response closely linked to intracranial aneurysms. The thrombus formation process is usually slow and complex as it is associated with blood flow and the net result of a series of biochemical reactions. The multi-scale and multi-physics nature of thrombosis has inspired a wide range of modeling approaches applied to various phenomena that aim to address how a thrombus forms.²⁷ Although different modeling methods can be coupled as informed by the scale and physics, the development of an all-encompassing computational model of thrombosis, combining all relevant underlying phenomena for patient-specific applications, remains impractical, and, instead, it is necessary to simplify models and to focus on specific questions.^27,28 The main interest of our thrombosis model is to investigate the thrombus composition/stability. Our novel model combines platelet activation and transport with fibrin generation, which is key to characterizing stable and unstable thrombus.⁸ Our model does not consider the fibrinolysis process, which involves the breakdown of a fibrin clot.²⁹ Therefore, the thrombus would keep growing according to the hemodynamics until a constant volume is reached and will not dissolve as the fibrinolysis process is not included in our model. The clot is to be formed in regions where the fibrin concentration is greater than 600 nM.^8,22,23 Based on this, we could assume that the fibrinolysis takes place where the concentration drops below this threshold¹³ if the fibrinolysis process needs to be included in our model in the future.

Given the complexity of thrombosis, with at least 80 coupled reactions that regulate thrombus growth,^30,31 a comprehensive computational thrombosis model considering both the hemodynamics and biochemical reactions is usually very complex and time-consuming with a large number of uncertain model parameters even after simplification. Our thrombosis model was originally developed by Sarrami-Foroushani et al.,⁸ where they built computer models of the in vitro phantom experiments and compared computational simulations of the flow diverter-induced thrombosis against in vitro observations reported in Ref. 32. A good agreement was achieved in that study. There are 31 model parameters in our thrombosis model⁸ with 8 biochemical reactions coupled to the transport of the blood flow. Previously, we calibrated the trigger thresholds (RT and SR thresholds) of thrombosis initiation as there is no consensus on the trigger thresholds, with different values used throughout the literature.^8,20,33 Building on this threshold calibration study,²⁰ in the present study, we performed a global SA to identify the most influential parameters and further validate our thrombosis model based on a real patient case. The unique most influential model parameter identified by the whole SA workflow is the resting platelet concentration, which means the concentration of resting platelets has the biggest effect on the final formed thrombus composition.

Our flow simulation model successfully predicted the spontaneous thrombosis status before treatment. The flow simulation is efficient (about 10 hours per run for the untreated aneurysm) compared with the time-consuming thrombosis model (about 10 days per run for the untreated aneurysm), but it can only provide hemodynamic information. To investigate the details of the formed thrombi, we ran the coupled hemodynamics and thrombosis model. As shown in Fig. 3, even using literature average values, our model is robust in predicting the main thrombi formation region. The concentration of patient-specific resting platelets primarily affects the composition/stability of the thrombi (Table I), which is consistent with our SA results. The resting platelet count distribution obtained from more than 400 000 cases from the UK Biobank can be found in Fig. S7, where 90% of the population was within $1.14×108$/ml and $4.38×108$/ml (Table S1). Given there is such a large range in the resting platelet concentrations for the general population, it is important to use patient-specific resting platelet concentration information when investigating the composition/stability of thrombi for individual cases.

We assumed the thrombus formation was triggered by blood flow stasis, which was characterized by high RT and low SR. Here, we further constrained the thrombosis initiation to only happen in regions near the wall or other thrombosed regions by adding a Hill function into the trigger mechanism for all internal points in the untreated aneurysm sac. This makes the thrombus formation process more realistic, as for untreated aneurysms, the thrombus is difficult to be suspended in the lumen on its own without any anchors to the surrounding aneurysm wall. As shown in Figs. 2–4 the simulation results of our thrombosis model show good agreement with the clinical ground truth both before treatment (spontaneous thrombosis) and immediately post-treatment (residual neck). It has been shown that aneurysmal thrombi form or at least deposit in regions of slow flow and low shear stress.^14,34–36 As illustrated in Fig. 4(c), there are flow stasis regions in the vicinity of the device within the main vessel that contribute to thrombus formation. In practice, this patient was treated with dual antiplatelet therapy consisting of aspirin and prasugrel (“75 mg of aspirin indefinitely and 10 mg of prasugrel for 6 months”). Dual antiplatelet therapy helps prevent stent-related thromboembolic events in cardiac patients and is commonly used during neurointerventional procedures.³⁷ This may explain why there were no thromboembolic events reported in the clinical records for this patient. The reported incidence of thrombus formation at the interface of the coil and the parent vessel is approximately 7%, based on retrospective analyses.^38,39

The post-treatment simulation is very time-consuming, taking months. The suggested future work from this study is to accelerate the thrombosis simulation to investigate the long-term post-treatment thrombus formation after patient-specific virtual endovascular treatments in an efficient way. Despite limitations due to the excessive run time, our results have demonstrated that our calibrated model can accurately predict the formed thrombus regions. The model could, therefore, be considered as a useful tool in clinical decision-making after further population-level and patient-specific validation studies, particularly when the run times are reduced and it becomes viable to use the model when planning treatments.

Limitations: (1) We assumed the thrombus formation in the aneurysm sac was triggered by blood flow stasis, and high RT and low SR were widely used in computational models to characterize the flow stasis. There may be other trigger mechanisms involved that were not included in our analysis. (2) The partially thrombosed regions were manually labeled with ITK-SNAP 3.8.0. There may be unavoidable errors due to subjective factors and the technical limitations of precise labeling of thrombosed areas. Although our model successfully predicted the main thrombosed region located in the left side of the aneurysm sac, the simulation results also showed that there is a small piece of thrombi formed in the right side of the aneurysm. In the raw images, we did not see apparent thrombi there. (3) As our thrombosis model is very time-consuming, we performed the global SA of our model only using a single geometry. Although our simulation results have demonstrated that our calibrated model can accurately predict the formed thrombus regions for the validation case, this study is limited to only one case. Further population-level validation studies should be investigated to increase the model's credibility.

In this comprehensive SA into all thrombosis model parameters and further clinical patient case validation study, we identified the unique most influential factor in aneurysmal thrombosis modeling, the resting platelet concentration. We also demonstrated that our thrombosis model is effective in predicting the thrombus formation both before and after treatment based on a clinical patient case, thereby further validating our model. Further large-scale validation studies across multiple patients are required to build additional trust in the model, but our results suggest there is significant value in using computational models to aid clinical decision-making.

The flow stasis-induced thrombosis model in our group was originally developed by Sarrami-Foroushani et al.⁸ Sarrami-Foroushani et al.⁸ assumed thrombosis to initiate and progress in regions where RT is greater than a threshold (e.g., 2.0 s) and SR is less than a threshold (e.g., 25 s⁻¹). These trigger thresholds were calibrated by Liu et al.²⁰ as 1.9 and 11 s⁻¹, respectively. As shown in Fig. 5, four main biochemically coupled events that result in a thrombus of fibrin mesh and aggregated platelets were considered, with five biochemical species: prothrombin (PT), thrombin (TH), antithrombin (AT), fibrinogen (FG), and fibrin (FI), and three categories of platelets: resting platelets (RP), activated platelets (AP), and fibrin bound aggregated platelets (BP). (1) Thrombin generation—conversion of prothrombin to thrombin on the surface of resting, activated, and bound platelets. Thrombin inhibition by antithrombin was also considered; (2) Fibrin generation—thrombin can convert fibrinogen into (insoluble) fibrin; (3) Platelet activation—resting platelets become activated by exposure to thrombin or other activated or bound platelets; (4) Platelet aggregation—activated platelets attach to the fibrin network and aggregate to form bound platelets. Details of the simulation specifications and equations for describing the blood flow transport and the biochemical reactions can be found in the supplementary material.

As shown in the supplementary material Table S1, there are 31 input parameters in our thrombosis model. The default values of these 31 parameters are obtained from literature.⁸ To perform comprehensive SA, we identify the upper and lower bounds of each parameter from the literature (e.g., upper and lower bounds of the fibrinogen concentration) or UK Biobank (e.g., upper and lower bounds of the resting platelet concentration). Where this information was not available, the upper and lower bounds of each parameter used in the following SA are set as a 200% (upper bound) or 50% (lower bound) variation of the literature default value.

There are five parameters that do not need to be included in the SA (Fig. 6 and Table S1). Three of these parameters are the initial concentrations of thrombin, fibrin, and bound platelets, as we assume that no thrombus formed before the thrombus formation process began. These three biochemical species are also the product of the associated biochemical reactions; for example, thrombin is generated by the conversion of prothrombin to thrombin on the surface of resting and activated platelets. Therefore, it is reasonable to set the initial concentrations of these three species to 0. The other two parameters that are not included in the SA are the RT and SR thresholds, as these were previously calibrated by Liu et al.²⁰ using the prevalence of clinical spontaneous thrombosis.

The thrombosis model is complex and time-consuming, requiring 5–20 days per case using 256 cores, as it combines hemodynamics with eight coupled biochemical reactions. To efficiently identify the most influential parameters, we first developed a lumped parameter model and used the least squares fitted linear model⁴⁰ to screen for key parameters from the 12 kinetic parameters. For those identified as key parameters using the lumped parameter model, along with the other 14 parameters, we then employed the 3D full thrombosis model and the elementary effect (EE) method⁴¹ to identify the most influential parameters. Finally, we validated the thrombosis model using a real patient case, which exhibited partial thrombosis before treatment and a residual neck after immediate post-treatment, by applying patient-specific values for the parameters identified as the most influential.

Previous studies have found that inlet flow boundary conditions of CFD models affect intracranial aneurysm hemodynamics, with inter-subject variability in cerebral blood flow found to be 10%–20%.^42–45 With in vivo measurement-derived boundary conditions⁴⁶ unavailable for our simulations, we used a previously developed multivariate Gaussian model (MGM)⁴⁷ to generate patient-specific boundary conditions as internal carotid flow waveforms. The MGM model was trained and calibrated using data from 17 healthy young adults. In this study, the patient-specific age (52 years old), gender (male), heart rate (68 bpm), systolic blood pressure (117 mm Hg), and diastolic blood pressure (73 mm Hg) information were used to generate the inflow boundary conditions. Poiseuille's law was used to scale the MGM-generated waveforms to achieve a time-averaged wall shear stress of 1.5 Pa at the inlet. At the outlets, the commonly used approaches are zero-pressure, Murray's law (divide outflows according to the cube of the diameter), and reduced-order models, with zero-pressure imposed by the majority of CFD teams.^48–50 In our previous posterior communicating artery (PComA) aneurysm study,⁵¹ we compared imposing 80:20, 60:40, 50:50, 40:60, and 20:80 mass flow splits at the middle cerebral artery (MCA) and the anterior cerebral artery (ACA) outlets with the default zero-pressure condition. We found that downstream conditions are not very important for aneurysm flow but that conditions applied to any branch vessels originating near the aneurysm are crucial (e.g., the PComA). For the validation case in this study, as the outlet branch vessels do not originate near the aneurysm, it is reasonable to impose zero-pressure boundary conditions.

In the simplified lumped parameter model, we assume that each species is fully diffused at its initial concentration. We then apply the least squares fitted linear model⁴⁰ to investigate the lumped parameter model using the duration of thrombus formation (the simulation time from 0% to the final 95% thrombosed) as the output metric.

The matrix X_nk has 1s in the first column and experimental values for the k parameters in the n simulations in the remaining columns. B_k contains the k + 1 unknown coefficients corresponding to the intercept b₀ and the k parameters. Y_n contains the n output values from the n simulations. Unless otherwise stated, we used n = 1000. If n is strictly greater than k + 1, it will not be possible to solve the equations exactly, unless the model is in fact linear, as the system of equations is overdetermined.⁴⁰ However, a least squares solution will generally be available. As solving an overdetermined system of equations for a least squares solution is computationally expensive and random samples tend to be poorly conditioned for large k because of clustering, we used Sobol sequences⁵² to generate the samples. Sobol sequences are a particularly common example of low-discrepancy sequences and exhibited faster convergence in comparison with random and Latin hypercube sampling (LHS) sampling, as has been demonstrated with correlated normal distributions in finance applications.⁵³ 

The computational cost to implement the EE method is $r(k+1)$⁠. Unless otherwise stated, we used p = 4 and r = 5 in this study.

After the comprehensive SA study, we identified the most influential parameters of our thrombosis model. To improve the model's performance and credibility, it is necessary to obtain patient-specific values for the parameters identified as influential ones. After anonymization, we used a patient case (partial thrombosis before treatment and residual neck after immediate flow diverter and coiling treatment; male; 52 years old; normotension; platelet count, $2.07×1011/ml$⁠; aneurysm size, 16.5 mm; aspect ratio, 2.2; location, MCA) from Leeds General Infirmary as a validation study case by comparing the thrombus regions predicted by our model with the clinical ground truth both before and after treatment. We collected the 3D rotational angiography images and detailed clinical records from Leeds General Infirmary [Fig. 7(a)], segmented the vasculature and aneurysm with a deep learning-based approach VASeg,⁵⁶ manually labeled the partially thrombosed regions with ITK-SNAP 3.8.0 [Fig. 7(b)], obtained the patient-specific vascular surface mesh, deployed the virtual stent and coils with GIMIAS (version 1.8.r1),²⁴ generated the volume mesh using ANSYS ICEM CFD v19.3 (Ansys Inc. Canonsburg, PA, USA), imposed the patient-specific inlet flow waveform generated from a multivariate Gaussian model,^42,47 ran the thrombosis model on an HPC cluster ARC4, and post-processed the simulation results with ANSYS CFD-POST and Paraview 5.10.0-RC1. The patient-specific clinical flow diverter model and coil information can be found in Table S2.

When modeling the stent, we are often concerned with alterations in the flow pattern and hemodynamics in the aneurysm rather than detailed flow fields near the walls of the parent vessel.⁵⁷ Given this, only the portion of the stent that crosses the neck of the aneurysm is modeled to reduce the computational expense of the simulations. The resolution of the mesh in the vicinity of the stent wires was set according to Sarrami-Foroushani et al.,⁸ where the independence of the mesh was obtained for the maximum edge size of 0.01 mm on the wires. The packing density, defined as the ratio of the volume of the physically inserted coil to the volume of the aneurysm, for this patient-specific case is 21.5%. The coils were discretized with a mesh resolution of 1.5 × the diameter of the primary coil.^58,59 The above settings resulted in volumetric meshes with 27 × 10⁶ total number of elements for this patient-specific case with a single flow diverter and 10 coils. The details of the mesh convergence analysis were previously described.^8,20

Rayz et al.³³ found that thrombus forms in layers, with the initial layer adhering to the arterial wall in regions of increased flow RT and then gradually expanding into the aneurysmal bulge. Moreover, imaging studies of untreated aneurysms often show thrombus formations that are closely associated with areas of altered flow patterns near the aneurysm wall, rather than suspended freely within the aneurysm sac.⁶⁰ For untreated aneurysms, the thrombus is difficult to be suspended in an aneurysm lumen on its own without any anchors to the aneurysm wall. Previously, we calibrated the trigger thresholds based on the prevalence of clinical spontaneous thrombosis (RT threshold 1.9 s and SR threshold 11 s⁻¹).²⁰ Here, we further constrain the thrombus to only initiate and progress near the wall or other thrombosed regions (where the fibrin concentration, C_FI, is greater than 600 nM.⁸).

In this study, we mainly consider saccular aneurysms, as the saccular type accounts for 90% of IAs.² Saccular aneurysms are spherical in shape, and the aneurysm sac can be approximated with a least squares fit ellipsoid or sphere.⁶¹ Therefore, we can separate all aneurysm points into near-wall points and internal points by a virtual ellipsoid or sphere inside the aneurysm sac. We set the center of the virtual ellipsoid or sphere to coincide with the center of the least squares fit ellipsoid or sphere. The radius of the virtual ellipsoid and sphere is set as half the value of the least squares fit ellipsoid or sphere radius. Then we used a Hill function to constrain the initiation of thrombus formation for all internal points. The Hill function is a sigmoidal activation function of the form $ϕpFI=CFIn/(CFIn+CFI,50n)$⁠, where the rate of occurrence of an event, p, requires an appropriate concentration of fibrin, $Ci,50$ is the fibrin concentration where the half maximum activation (half saturation) occurs, and the Hill coefficient (the exponent n) reflects the steepness of the response curve. In this study, we set $Ci,50=$ 600 nM and n = 4 in the trigger mechanism.

Spontaneous thrombosis of unruptured intracranial aneurysms is a common event that can be detected incidentally during advanced neuroradiological studies before treatment.^62–64 These spontaneously thrombosed aneurysms are considered unstable dynamic structures that may grow, recanalize, bleed, compress, or cause thromboembolic events.^63–65 The validation case may eventually become complete spontaneous thrombosis. However, it was treated when the partial thrombosis was fresh and unstable, as noted by the clinician. Complete spontaneous thrombosis can sometimes stabilize the growth of the lesion; however, 33% (7/21) of the completely thrombosed aneurysms presented recanalization at follow-up.⁶² The spontaneously formed thrombus was neither stable nor mature before treatment, making it inappropriate to compare a converged simulation result with an unstable and unconverged clinical ground truth. For the pre-treatment simulations, our approach is to select the simulation time that best matches the clinical ground truth. For instance, in the non-constrained initiation model, the thrombus formed at 30 s of simulation time closely resembles the clinical ground truth, so we used this time point. Consequently, for the constrained initiation model, we also ran the simulation for 30 s.

See the supplementary material for the governing equations, simulation specifications, and other additional images and data that support the findings of this study.

Professor Frangi acknowledges support from the Royal Academy of Engineering under the RAEng Chair in Emerging Technologies (INSILEX CiET1919/19) and the ERC Advanced Grant - UKRI Frontier Research Guarantee (INSILICO EP/Y030494/1). A.F.F. acknowledges seminal funding from the European Commission to @neurIST (FP6-2004-IST-4-027703) and the @neurIST Consortium. The NIHR—Manchester Biomedical Research Centre (BRC) (NIHR203308) and the British Heart Foundation (BHF) Manchester Centre of Research Excellence (RE/24/130017) support the work of A.F.F. The views expressed in this publication are those of the author(s) and not necessarily those of any of the funders. M.M. was funded by the EPSRC Centre for Doctoral Training in Fluid Dynamics (EP/L01615X/1) while at the University of Leeds, Leeds, UK.

The authors have no conflicts to disclose.

Ethics approval is not required.

Qiongyao Liu: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). Toni Lassila: Conceptualization (equal); Formal analysis (equal); Methodology (equal); Supervision (equal); Writing – review & editing (equal). Fengming Lin: Methodology (equal); Resources (equal); Software (equal); Visualization (equal). Michael MacRaild: Methodology (equal); Writing – review & editing (equal). Tufail Patankar: Conceptualization (equal); Data curation (equal); Methodology (equal); Resources (equal). Fathallah Islim: Data curation (equal); Resources (equal). Shuang Song: Methodology (equal); Software (equal). Huanming Xu: Methodology (equal); Software (equal). Xiang Chen: Writing – review & editing (equal). Zeike A. Taylor: Conceptualization (equal); Methodology (equal); Supervision (equal). Ali Sarrami-Foroushani: Conceptualization (equal); Methodology (equal); Supervision (equal). Alejandro F. Frangi: Conceptualization (equal); Funding acquisition (equal); Methodology (equal); Project administration (equal); Resources (equal); Supervision (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.