Uncertainty quantification of computational fluid dynamics-based predictions for fractional flow reserve and wall shear stress of the idealized stenotic coronary

Introduction: In clinical practice, hemodynamics-based functional indices, such as fractional flow reserve (FFR) and wall shear stress (WSS), are useful in assessing the severity and risk of rupture of atherosclerotic lesions. Computational fluid dynamics (CFD) is widely used to predict these indices noninvasively.

Method: In this study, uncertainty quantification and sensitivity analysis are performed for the computational prediction of WSS and FFR directly from 3D–0D coupled CFD simulations on idealized stenotic coronary models. Five geometric parameters (proximal, mid, and distal lengths of stenosis; reference lumen diameter; and stenosis severity) and two physiological parameters (mean aortic pressure and microcirculation resistance) are considered as uncertain input variables.

Results: When employing the true values of stenosis severity and mean aortic pressure, a discernible reduction of 25% and 9.5% in the uncertainty of the computed proximal WSS, respectively. In addition, degree of stenosis, reference lumen diameter, and coronary resistance contributed the uncertainty of computed FFR, accounting for 41.2%, 31.9%, and 24.6%, respectively.

Conclusion: This study demonstrated that accurate measurement of the degree of stenosis and mean aortic pressure is crucial for improving the computational prediction of WSS. In contrast, the reference lumen diameter, degree of stenosis, and coronary resistance are the most impactful parameters for FFR.

1. Introduction

Coronary artery disease (CAD) is the most common cause of mortality worldwide (1). Significant narrowing of the coronary lumen owing to endothelial plaque accumulation can interfere with the flow of oxygen-rich blood to the myocardium, thereby causing ischemic symptoms and myocardial infarction. Traditionally, a morphology-based index, the percentage reduction of arterial lumen diameter in coronary angiography images, has been used to define the severity of stenosis and as a decision-making tool for revascularization. However, multiple clinical trials have reported limited diagnostic accuracy and inferior correlation with prognostic outcomes. This finding emphasizes the need for a physiology-based functional index (2–4).

Pressure-wire-based fractional flow reserve (FFR) is defined as the maximal coronary flow to the myocardium in the presence of stenosis divided by the theoretical maximal blood flow in a normal coronary artery, and it is currently the reference standard for identifying hemodynamically significant stenotic lesions in coronary circulation (5–9). This was approximated as the ratio of the mean coronary pressure measured distal to the stenosis to the mean aortic pressure under hyperemic conditions. To measure FFR, an invasive pressure wire and pharmacological agents that induce maximal vasodilation of the microcirculation are required (10). In maximal vasodilation (hyperemia), the distal coronary pressure is directly proportional to the maximum vasodilated perfusion and coronary flow (7). FFR exhibits superior diagnostic performance compared with the morphology-based index (11, 12); however, the practical application rate of FFR in catheterization laboratories is only approximately 10%, possibly owing to potential risks during measurements (13). Noninvasive computational FFR methods have been developed as alternatives and showed promising potential with excluding the use of additional pressure wire and hyperemic agents (14). In this approach, FFR is computationally evaluated by solving the Navier-Stokes equations for coronary velocity and pressure fields using computational fluid dynamics (CFD) techniques based on 3D reconstructed vascular geometry derived from coronary CT angiography or x-ray angiography.

In addition, the local hemodynamic force, i.e., the wall shear stress (WSS) is widely accepted to play an important role in the development and progression of atherosclerosis (15, 16). In particular, high WSS is correlated with plaque rupture in severe stenosis (>70%); at minimal CAD, low WSS is associated with rapid plaque progression (17). Kumar et al. (18) reported that high-risk plaque rupture is more prevalent at the proximal segment of a stenotic lesion, which is presumed to be related to the increased time-varying structural strain caused by a high WSS. Furthermore, proximal WSS combined with FFR is shown to have an incremental prognostic value in the prediction of myocardial infarction compared with FFR alone.

With remarkable advancements in computer resources, computational modeling, and image processing techniques, CFD has been applied to predict patient-specific intravascular hemodynamics noninvasively. This approach requires input data, including patient-specific clinical data and the 3D lumen geometry of the coronary artery, which is typically reconstructed from computed tomography and x-ray angiography. However, owing to inevitable uncertainties associated with the inherent variability of physiological data, and the inaccuracy of in vivo measurements and lumen segmentation process of stenosis, the effect of these input uncertainties on the prediction accuracy of hemodynamics-based diagnostic indices must be quantified.

Although Monte Carlo simulation is a conventional and simple method for the analysis of uncertainty quantification (UQ), it is not feasible for complex problems because its convergence rate follows the principle of large numbers, i.e., its rate is approximately 1/√N, where N is the number of samples. Another technique developed in the recent decades to overcome the issue of long computing time is the polynomial chaos expansion (PCE) method (19, 20). The PCE approximates a random variable as a linear combination of polynomial functions of other random variables and is faster than the Monte Carlo method when the number of uncertain parameters is less than 20 (21). The PCE has been used to successfully solve various problems in fluid mechanics, including laminar boundary layer flow on a flat plate, supersonic flow over a convex corner, and inviscid flow around a three-dimensional wing, as illustrated in the study conducted by Hosder et al. (22).

Eck et al. (23) previously introduced the utilization of UQ and SA in the context of cardiovascular applications. They investigated on the variability of FFR using a 1D model with ten uncertain parameters: proximal length, stenosis length, distal length, proximal radius, stenosis radius, distal radius, hyperemic blood flow, arterial pressure, blood viscosity, and blood density. Concurrently, they applied both the Monte Carlo (MC) method and the Polynomial Chaos (PC) method to assess the merits and demerits of each technique. The findings revealed that among the uncertain inputs, stenotic radius, hyperemic flow, and arterial pressure played pivotal roles in contributing to the variance of FFR. Notably, the PC method demonstrated a substantial advantage over the Monte Carlo (MC) method. The PC method required a significantly lower number of samples—2,002 samples compared to the MC method’s 60,000 samples. Additionally, the distribution of FFR obtained through the PC method closely resembled that obtained via the MC method, but at a significantly lower computational cost.

Sankaran et al. (24) undertook a comprehensive investigation into the influence of uncertainties of minimum lumen diameter, lesion length, boundary conditions, and blood viscosity on blood flow and pressure in both an idealized stenosis model and a patient-specific model. They employed an adaptive stochastic method integrated with a data-driven approach to analyze the effects of these uncertainties. Their findings revealed that the minimum lumen diameter emerged as the most influential factor affecting hemodynamic simulations.

Gashi et al. (25) studied the effect of the model-order reduction approach (2D, 3D, steady, and unsteady) on the computationally predicted FFR and reported that stenosis severity is the dominant geometric parameter for all cases.

Although these studies provide valuable information regarding the influence of the uncertainty of geometrical and physiological features on computed FFR, they were based only on 1D model simulations or were limited to a systematic investigation of the interactive contribution of uncertain input parameters to variations in FFR output. Furthermore, information regarding the uncertainty of stenosis WSS predictions is lacking.

In this study, an investigation was conducted to assess the relative significance of input variables and the repercussions of their uncertainties on computational predictions of FFR and WSS. Utilizing uncertainty quantification (UQ) and sensitivity analysis (SA), the nonintrusive PCE method was employed. The analytical approach involved an unsteady 3D-0D coupled CFD analysis, employing idealized coronary artery models rather than a simplified 1D approach. Lumped parameter networks were used to model the blood flow within the coronary artery, with the assumption that the coronary resistance during hyperemia is a quarter of the coronary resistance in the normal state. To account for geometric uncertainties, an eccentric stenosis case was considered due to its prevalence in CAD. Seven prospective uncertain inputs were considered: proximal length, middle length, distal length, reference lumen diameter of the coronary artery (representing the healthy diameter), degree of stenosis, mean aortic pressure, and coronary resistance. Additionally, the convergence of the nonintrusive PCE method was addressed by comparing various statistical indices at different polynomial orders.

2. Materials and methods

2.1. Geometry modeling and uncertain input variables

To investigate the effects of uncertainty in the input variables on FFR and WSS, coronary models were generated based on a synthetically designed single conduit with various geometric parameters of stenosis, as shown in Figure 1. Five uncertain geometric features were considered: proximal, mid, and distal lengths of stenosis; reference lumen diameter of the coronary artery; and diameter-based stenosis severity. The eccentricity of the stenosis was fixed at 35% in all the models. In addition, the mean aortic pressure and normalized microcirculation resistance were included as uncertain physiological variables because the flow rate and cardiac pressure generally play important roles in coronary flow dynamics. The details of the mean value and uncertainty levels for each input variable considered in this study are presented in Table 1. The degree of uncertainty for the input variables was assumed based on the previous study, in which in vivo data typically observed in patients were considered (23). Mean values of vessel geometric parameters were determined by analyzing approximately 400 patients’ left anterior descending (LAD) coronary profiles from the literature (26–29). The mean aortic pressure (P_a) of 90 mm Hg was applied. The microcirculation resistance (R) was set to match the mean flowrate of 1 ml/s under baseline conditions and was reduced by four times (R/4) for the hyperemic condition (30).

FIGURE 1

Figure 1. Schematic of an idealized stenotic coronary model.

TABLE 1

Table 1. Descriptions of uncertain input variables of the computational model.

An idealized stenotic coronary model with a 3.6 mm-diameter was constructed based on the following equations for the geometric contours.

$d(z)=\{\begin{array}{c}D\left[1-s_0\sin {\displaystyle \frac{\pi }{2L_i}(z-z_0)}\right],z_0\le z\le z_1\\ D[1-s_0],z_1\le z\le z_2\\ D\left[1-s_0\cos {\displaystyle \frac{\pi }{2L_0}(z-z_2)}\right],z_2\le z\le z_3\end{array},(1)$

$e(z)=\{\begin{array}{c}{e}_0\sin {\displaystyle \frac{\pi }{2L_i}(z-z_0),z_0\le z\le z_1}\\ e_0\left(={ϵ}_0{s}_0{\displaystyle \frac{D}{2}}\right)z_1\le z\le z_2\\ e_0\cos {\displaystyle \frac{\pi }{2L_o}(z-z_2),z_2\le z\le z_3}\end{array},(2)$

$\left[\begin{array}{c}x(z)\\ y(z)\end{array}\right]=\left[\begin{array}{c}0\\ e(z)\end{array}\right]+{\displaystyle \frac{d(z)}{2}\left[\begin{array}{c}\cos \theta \\ \sin \theta \end{array}\right]0\le \theta \le 2\pi ,}(3)$

Here, d(z) and e(z) represent the vessel diameter at the stenotic segment and the central axis of stenosis along its length, respectively. e₀ represents the maximum off-center distance of the stenotic segment; D is the reference lumen diameter of the coronary artery; s₀ is the degree of stenosis.

The eccentricity of stenosis ε₀ is defined as the percentage measure of asymmetry or off-center narrowing within the arterial lumen. Previous studies have demonstrated a higher prevalence of eccentric stenosis compared to concentric stenosis. Mintz et al. (31) reported an occurrence of 795 out of 1,446 (55%) eccentric lesions based on angiography data. Similarly, Yamagishi et al. (32) investigated the morphological characteristics of 114 coronary plaques using intravascular ultrasound and observed that 92 (80%) cases of stenosis exhibited an eccentric pattern. However, Seo et al. (33) reported that lesion eccentricity has no statistically significant effect on FFR. Based on these characteristics, the coronary model in this study considered a fixed eccentricity of 35%, which significantly reduces the CFD computation time.

Eventually, the uncertainties in the variables are propagated to the CFD model via the Latin hypercube sampling method (33) with a uniform distribution.

2.2. Computational fluid dynamics

3D–0D coupled CFD simulations for stenotic coronary flow were performed (34, 35). The flow was assumed to be transient, and the blood was an incompressible Newtonian fluid with a density of 1,050 kg/m³ and kinematic viscosity of 3.5 × 10⁻⁶ m²/s. The unsteady Navier–Stokes equations describing the conservation of mass and momentum of the blood flow are as follows.

$\rho {\dot{u}}_i+\rho u_j{u}_{i,j}=-p_{,i}+({\tau }_{i,j}{)}_{,j^{\mathrm{\prime }}}(4)$

$u_{i,i}=0.(5)$

For the numerical calculations, these governing equations were split into four steps for time advancement based on a fully implicit fractional step method (36) as follows.

$\mathrm{Step}1:\rho {\displaystyle \frac{{\hat{u}}_i-u_i^n}{\mathrm{\Delta }t}+\rho {\displaystyle \frac{1}{2}({\hat{u}}_j{\hat{u}}_{i,j}+u_j^n{u}_{i,j}^n)=-p_{,i}^n+{\displaystyle \frac{1}{2}({\hat{\tau }}_{ij}+{\tau }_{ij}^n{)}_{,j}}}}(6)$

$\mathrm{Step}2:\rho {\displaystyle \frac{u_i^{\ast }-{\hat{u}}_i}{\mathrm{\Delta }t}=p_{,i}^n,}(7)$

$\mathrm{Step}3:p_{,ii}^{n+1}={\displaystyle \frac{\rho }{\mathrm{\Delta }t}{u}_{i,i}^{\ast },}(8)$

$\mathrm{Step}4:\rho {\displaystyle \frac{u_i^{n+1}-u_i^{\ast }}{\mathrm{\Delta }t}=-p_{,i}^{n+1},}(9)$

where Δt is the time increment; the superscript n indicates the time level; and ${\hat{u}}_i$ and $u_i^{\ast }$ are intermediate velocities. The second-order implicit Crank–Nicolson scheme was employed for the diffusion and convection terms.

To accommodate the specific time-varying pressure-flow characteristics of coronary flow, the 0D lumped parameter network (LPN) modeling technique for coronary microcirculation was integrated into 3D CFD at the outlet, as shown in Figure 2. The LPN model represents the coronary circulatory system by combining resistance and capacitance elements based on the similarity between the hydraulic and electric systems. The total resistance and capacitance of the distal coronary bed at rest (baseline) and under the hyperemic condition, which is the maximum flow condition achieved by minimizing microvascular resistance, were determined based on previous studies conducted by Sankaran (24) and Sharma et al. (37). The MPI parallel algorithm was applied to reduce the computation time to a simulation time within 1 h (about 20 min with 64 parallel cores) (35).

FIGURE 2

Figure 2. Schematic of a 3D–0D coupled CFD approach for an idealized stenotic coronary flow simulation.

A P2P1 finite element scheme (35) was employed with quadratic tetrahedral elements generated using a commercial mesh generator (ICEM-CFD, ANSYS) to achieve higher accuracy. A mesh convergence study was performed by comparing the stenotic pressure drop $\mathrm{\Delta }{p}_{stenosis}$ in the case of different mesh densities, shown in Figure 3. The difference of $\mathrm{\Delta }{p}_{stenosis}$ between the medium mesh (the number of nodes ∼1 million) and the densest mesh (the number of nodes ∼1.8 million) is within 0.5%; Meanwhile, the computation time significantly increases from 1,400 s to 3,000 s with the same computational resources (56 parallel cores, Intel Xeon E5-2630 v4). In this way, the mesh with the number of nodes of 1 million was chosen for this study.

FIGURE 3

Figure 3. Mesh independent test for the idealized stenotic coronary.

For the inlet boundary conditions, a typical aortic pressure waveform (39), scaled using a mean aortic pressure of 90 mmHg, was applied. The duration of each cardiac cycle is 1.0 s and is divided into 120 timesteps. The convergence criteria are $1\times {10}^{-6},5\times {10}^{-6}$ for velocity and mass, respectively. Figure 4 illustrates the blood flow rate under resting (baseline) and hyperemic conditions and the pressure waveforms for a cardiac cycle.

FIGURE 4

Figure 4. Aortic pressure and flow rate waveforms under resting and hyperemic conditions of a cardiac cycle. (Aortic pressure was applied at the inlet as a boundary condition and resting and hyperemic flow rates are computed results under corresponding conditions).

2.3. Point-collocation nonintrusive polynomial chaos method

In the nonintrusive polynomial chaos (NIPC) method, a surrogate model f_PC based on polynomials for output Y is generated, as follows.

$Y\approx f_{\mathrm{PC}}(\mathrm{X})=\sum _{i=1}^{N_P}{\alpha }_i{\mathrm{\Psi }}_i(X),(10)$

where Ψ_i is the expansion or polynomial (typically orthogonal polynomial) function in relation to the probability distribution of inputs $X;{\alpha }_i$ is a coefficient to be computed; and N_p is the total number of discretization terms.

In practice, N_p is expressed by the following equation in terms of the number of uncertain input parameters n and the order of the polynomial p:

$N_p={\displaystyle \frac{(n+p)!}{n!p!}}(11)$

The inputs $\mathit{X}((X_i),i=1,2,..,n)$ are treated as uncertain variables and hypothesized independent; therefore, the random space of uncertain inputs X is drawn by the joint probability distribution ρ_x, which is consisted of other random variables X_i:

${\rho }_X(x)=\prod _{i=1}^n{\rho }_{X_i}(x_i)(12)$

The point collocation NIPC method determines the components a_i by choosing N_s vectors $\{{\mathit{X}}_j={(X_1,X_2,\dots ,X_n)}_j,j=1,\mathrm{2..},N_s\}$ from the random space and subsequently evaluating the outputs of interest $\{Y_j,j=1,\mathrm{2..},N_s\}$ by deterministic CFD code. Corresponding to each chosen vector there will be an output, a linear system of equations can be established:

$\left(\begin{array}{cccc}{\mathrm{\Psi }}_1(X_1)& {\mathrm{\Psi }}_2(X_1)& \dots & {\mathrm{\Psi }}_{N_P}(X_1)\\ {\mathrm{\Psi }}_1(X_2)& {\mathrm{\Psi }}_2(X_2)& \dots & {\mathrm{\Psi }}_{N_P}(X_2)\\ ⋮& ⋮& \ddots & ⋮\\ {\mathrm{\Psi }}_1(X_{N_s})& {\mathrm{\Psi }}_2(X_{N_s})& \dots & {\mathrm{\Psi }}_{N_P}(X_{N_s})\end{array}\right)\left(\begin{array}{c}{\alpha }_1\\ {\alpha }_2\\ ⋮\\ {\alpha }_{N_P}\end{array}\right)=\left(\begin{array}{c}{Y}_1\\ Y_2\\ ⋮\\ Y_{N_s}\end{array}\right)(13)$

$\alpha =\mathrm{argmi}{\mathrm{n}}_{\alpha }\Vert \mathrm{\Psi }\alpha -Y{\Vert }_2(14)$

The number of chosen samples N_s should be equal or greater than number of discretization terms N_p; in the case of $N_s>N_p$, the deterministic problem becomes to over-determined system of equations (Eq. 13), then the least squares method can be used to solve (Eq. 14). Hosder et al. (22) suggested that twice as many samples as the minimum (practical) are required to better approximate the statistics at each polynomial degree $(N_s=2N_p)$. Here, $n_{ps}$ is denoted as the proportion of sample:

$n_{ps}={\displaystyle \frac{N_s}{N_p}}(15)$

To apply the NIPC method for FFR and WSS, the uncertain input parameters are predefined with a probability distribution, then these parameters are sampled and applied to a deterministic solver to obtain outputs of interest (FFR and WSS); finally, the calculation of UQ and SA was performed. In this study, uncertain input parameters were identified as geometric $(L_i,L_m,L_o,D,s_0)$ and physiological $(P_a,R)$ features. The uncertain parameters were hypothesized independently distributed and sampled with joint distribution by Latin hypercube to select the collocation points; well-validated in-house CFD code was used to calculate FFR and WSS (35).

The outputs of interest (here, FFR and WSS) are considered as random variables. Practically, the probability density functions of FFR and WSS $({\rho }_{FFR},\mathrm{and}{\rho }_{WSS})$ are dependent on the uncertain inputs $(L_i,L_m,L_o,D,s_0,P_a,R)$, which are governed by the CFD model. The characteristics of ${\rho }_{FFR}$ and ${\rho }_{WSS}$ were assessed using uncertainty quantification (UQ) technique, and the contribution of an uncertain input $X_i$ to the outputs of interest was evaluated using the sensitivity analysis (SA).

The convergence rate generally depends on both the number of samples N_s and the order of the polynomial p. Four different maximal polynomial orders were considered $(p=1,2,3\mathrm{and}4)$, with the proportion of sample $n_{ps}=1\mathrm{and}2$ were used to assess the convergence of UQ measures.

The Python package Chaospy (38) was used to generate stochastic samples and calculate the surrogate model for the simulation output of interest.

2.4. Uncertainty quantification and sensitivity analysis

Characteristics of unknown distribution of Y using UQ are described as statistical moments (e.g., mean and variance). In the NIPC method, a surrogate model f_PC consisting of orthogonal polynomials used to replace Y; therefore, the mean of unknown distribution of Y can be directly obtained using the coefficients and the orthogonality property. For example, the mean and variance of Y are given by:

$\mu [Y]=\mathbb{E}[Y]=\underset{\mathrm{\Omega }}{\int }y\rho (y)dy=\underset{\mathrm{\Omega }}{\int }\sum _{i=1}^{N_P}{\alpha }_i{\mathrm{\Psi }}_i(x){\rho }_X(x)dx={\alpha }_1(16)$

$\begin{array}{c}Var(Y)=\mathbb{E}[Y^2]-\mathbb{E}{[Y]}^2={\displaystyle \underset{\mathrm{\Omega }}{\int }{\left({\displaystyle \sum _{i=1}^{N_P}{\alpha }_i{\mathrm{\Psi }}_i(x)}\right)}^2{\rho }_X(x)dx-{\alpha }_1^2}\\ ={\displaystyle \sum _{i=1}^{N_P}{\alpha }_i^2\underset{\mathrm{\Omega }}{\int }{\mathrm{\Psi }}_i{(x)}^2{\rho }_X(x)dx-{\alpha }_1^2=\sum _{i=2}^{N_P}{\alpha }_i^2\underset{\mathrm{\Omega }}{\int }{\mathrm{\Psi }}_i{(x)}^2{\rho }_X(x)dx}\end{array}(17)$

The contribution of the uncertain input parameters to the variation in the output was investigated by sensitivity analysis. The final goal is determining how the uncertainties in the input variables contribute to the variance of the output, either individually or through interactions with other parameters. This is useful for modeling personalization to determine the parameters that need to be optimized to their true values (input prioritization) and those that can be fixed within their uncertainty domains (input fixing). The main and total sensitivity indices ($S_i$ and $S_i^T$, respectively) (23) were computed to represent input prioritization and input fixing, respectively.

The main sensitivity index, also called the first-order Sobol sensitivity index, is the variance of the conditional expectation of output Y, given the value of input $X_i$, normalized by the total variance:

$S_i={\displaystyle \frac{\mathrm{Var}[E[Y|X_i]]}{\mathrm{Var}[Y]}.}(18)$

Here, the index i varies from 1 to the number of uncertain input variables $n,1\le i\le n$.

The main sensitivity index $S_i$ represents the expected reduction in the total variance $\mathrm{Var}[Y]$ when the input variable $X_i$ is corrected to its true value.

The total sensitivity index $S_i^T$ includes the sensitivity of both the first-order effects and the interactions (covariance) between a given parameter $X_i$ and all other parameters:

$S_i^T=1-{\displaystyle \frac{\mathrm{Var}[E[Y|X_{-i}]]}{\mathrm{Var}[Y]},}(19)$

where $X_{-i}$ is the set of all uncertain input variables except $X_i$.

The convergence of the UQ measures was assessed by evaluating the mean value obtained with different polynomial orders p and the proportion of sample $n_{ps}$.

3. Results

The blood flow and pressure fields for the mean values of the input variables in the baseline and hyperemic conditions from the CFD simulations are shown in Figure 5. An increase of approximately 3.3-fold in the mean flow rate was observed when the microvascular resistance was decreased by a factor of four relative to the baseline condition to model the hyperemic condition. In addition, a significant difference of 2.2 mmHg and 12.7 mmHg was observed in pressure drop between baseline and hyperemic conditions. This difference is primarily due to increased flow separation at higher flowrates and possibly turbulence.

FIGURE 5

Figure 5. Representative blood flow and pressure fields obtained from CFD simulations for resting (baseline) and hyperemic conditions. (A) Baseline condition (R_base) (B) Hyperemic condition $(R_{hyper}=0.25\times R_{base})$.

3.1. UQ and SA for WSS

Slager et al. (39) reported a tendency to form high-risk plaques localized in stenotic proximal segments with a high WSS. Herein, UQ and SA were conducted to investigate the impact of uncertain input variables on the uncertainty of the average proximal WSS $(AWS{S}_{prox})$. A surrogate model was derived from 3D CFD data using PCE.

Figure 6 shows a comparison of the total sensitivity index $S_i^T$ of $AWS{S}_{prox}$ for various combinations of polynomial order p and sample size $n_{ps}$. The behavior of $S_i^T$ exhibited a similar trend for $(p,n_{ps})$ combinations higher than $(p=3,n_{ps}=2)$, which suggests suitable convergence.

FIGURE 6

Figure 6. Trend of total sensitivity indices of AWSS_prox for each combination of polynomial order and sample size (p, n_ps).

The degree of stenosis S_o and the mean aortic pressure $P_a$ were always in the highest-sensitivity group, with a main sensitivity greater than 50%, whereas the other variables ranged from 10% to 20%. In addition, the probability density distributions for each combination of $(p,n_{ps})$ were compared, as shown in Figure 7. The probability density distribution by $(p=3,n_{ps}=2)$ combination exhibits nearly the same asymmetric shape as the one by $(p=4,n_{ps}=2)$ combination with only minor difference of the mean and standard deviation, which indicates the results with the number of polynomials $p=3$ and $p=4$ approach the convergence, particularly, when oversampled as $n_{ps}=2$. However, it is clearly exhibited that the proportion of sample, $n_{ps}=1$ is not sufficient for the approximation of the statistics, even though high polynomial order $p=4$ is applied.

FIGURE 7

Figure 7. Probability density of AWSS_prox for each combination of polynomial order and sample size (p, n_ps).

The mean of $AWS{S}_{prox}$ obtained with different combinations of $(p,n_{ps})$ are given in Figure 8, showing that the obtained mean value tends to converge (relative error < 7.8%) as the polynomial order increases.

FIGURE 8

Figure 8. Mean of AWSS_prox for each combination of polynomial order and sample size (p, n_ps).

Table 2 lists the required number of samples (3D CFD realizations) and the corresponding computing times associated with each combination $(p,n_{ps})$. Based on this observation, the combination $(p=3,n_{ps}=2)$ is optimal in terms of the computational cost and prediction accuracy.

TABLE 2

Table 2. Number of samples (computing time; hours) for each combination of (p, n_ps).

As shown in Figure 9, the stenosis severity $s_0$ exhibited the highest main sensitivity index, followed by the mean aortic pressure $P_a$, whereas the remaining variables exhibited only a minor effect $(\le 5\text{\%})$. More precisely, these sensitivity indices indicate that if s_o and $P_a$ are changed to their true values, the uncertainty in $AWS{S}_{prox}$ prediction is reduced by 25% and 9.5%, respectively. The total sensitivity index exhibited a trend similar to that of the main sensitivity index, that is, the greatest influence was from stenosis severity $s_0$ and mean aortic pressure $P_a$. However, the difference between the main and total sensitivity indices was considerably large. This indicates that the interaction of the input variables in the computational model for WSS prediction is significant; thus, multivariate uncertainty analysis is crucial.

FIGURE 9

Figure 9. Sensitivity indices for computed AWSS_prox.

3.2. UQ and SA for FFR

For FFR computations, the same conditions for uncertain input variables were applied, except for microvascular resistance, which was decreased by a factor of 4 relative to the baseline condition, to ensure a hyperemic condition. The FFR was evaluated 20 mm downstream of the stenosis in the pressure fields from the CFD simulations.

Figure 10 shows the statistical characteristics of the probability density of the FFR due to the uncertainty of the input parameters, including the expected values, standard deviation, and 95% prediction interval. Evidently, the distribution is asymmetric and left-skewed toward lower FFR values, unlike $AWS{S}_{prox}$ which is right-skewed toward higher FFR values.

FIGURE 10

Figure 10. Probability density of computed FFR.

Figure 11 shows that only the coronary reference diameter D, stenosis severity $s_0$, and microvascular resistance R significantly affected the uncertainty in the FFR prediction. However, the influence of the mean aortic pressure was <0.05. The difference between the main and total sensitivity indices was nearly zero, thus indicating a negligible interaction of the input variables that are independent in the computational model for FFR prediction.

FIGURE 11

Figure 11. Sensitivity indices for computed FFR.

4. Discussion

Herein, UQ and SA were performed for the computational prediction of the WSS and FFR directly from 3D–0D coupled CFD simulations of idealized stenotic coronary models. The mean degree of stenosis was 50% based on diameter reduction. Five geometric parameters for the coronary stenosis model (proximal, mid, and distal lengths of stenosis; reference lumen diameter; and diameter-based stenosis severity) and two physiological parameters (mean aortic pressure and normalized microcirculation resistance) were considered as uncertain input variables for the CFD simulations. Although viscosity and density may be additional blood characteristics that influence the WSS and FFR values, they were not included in this study because the uncertainty of the parameters is not significant. Furthermore, Sankanran et al. (24) reported that viscosity and density have negligible effects on the uncertainty of the computed FFR.

Eck et al. (23) demonstrated that the stenotic radius, hyperemic flow rate, and arterial pressure contributed the most to the uncertainty in the FFR using 1D model-based simulation data. This finding is consistent with the results of the present study. Herein, the reference vessel diameter, stenosis severity, and microvascular resistance were found to exhibit the greatest impact on FFR prediction. Microvascular resistance under hyperemic conditions primarily determines the coronary hyperemic flow rate. In clinical practice, the FFR value is measured under the assumption of zero microvascular resistance after inducing a hyperemic condition with a drug; in reality, the value is decreased by three–five times owing to the resistance, which varies among patients. In the present study, arterial pressure influenced the uncertainty of FFR prediction; however, the relative sensitivity was low.

Sankaran et al. (24) showed that the minimum lumen diameter and outlet boundary resistance are the most important factors contributing to the variance of the computed FFR, and that their influence is independent and not interactive. Evidently, similar results were observed herein because the combined reference diameter and degree of stenosis were directly associated with the minimum lumen diameter. The interactive effect of input uncertainty on the FFR was also negligible.

For the average proximal WSS $(AWS{S}_{prox})$, the mean aortic pressure had the greatest impact on the prediction uncertainty. This is because the WSS was calculated in the baseline flow condition, and we expect that aortic pressure has a greater effect on the baseline flow rate with higher outlet resistance than in a hyperemic condition. In addition, because WSS is inversely proportional to the third power of the lumen radius, stenosis severity has a large impact on the uncertainty of $AWS{S}_{prox}$.

Although the uncertainty of segmental lengths of a stenosis, such as $L_i$, $L_m$, and $L_o$, appears to have a negligible influence on the uncertainty of the computed FFR (indicating that precise measurements of the stenosis length in the imaging and reconstruction processes are not crucial), they exhibited relatively higher sensitivity to the variation of the computed $AWS{S}_{prox}$. In particular, the total sensitivity index of the proximal length of stenosis $L_i$ was approximately the same as that of microvascular resistance in the computed $AWS{S}_{prox}$. This implies strong interactions between the length of the stenosis and other input parameters, which should not be excluded from the distinctive input parameters in the computational prediction of WSS.

Moreover, the wide inter-lesion variations in patients may produce a discernible difference in FFR. As depicted in Figure 12 as a representative case, doubling the stenosis length of a stenotic vessel while preserving the other input parameters produced a notably large pressure drop and in turn a lower FFR. This suggests that the stenosis length, with appreciation of the presence of wide interlesion variations, plays a crucial role in the determination of FFR and can be considered an independent geometric factor.

FIGURE 12

Figure 12. Comparison of computed FFR with different stenosis lengths (9 mm and 18 mm) (A) $s_0=50\text{\%}$ (B) $s_0=60\text{\%}$.

5. Conclusion

In this study, UQ and SA were performed to quantify the effect of uncertainty of the geometric and physiological input parameters including proximal, mid, and distal lengths of stenosis, reference lumen diameter, stenosis severity, mean aortic pressure and microcirculation resistance, on the computational prediction of the WSS and FFR with an idealized stenotic coronary model. The degree of stenosis and the mean aortic pressure exhibited highest sensitivity to variations in the computed WSS. When employing the true values of stenosis severity and mean aortic pressure, a discernible reduction of 25% and 9.5% in the uncertainty of the computed proximal WSS, respectively. In addition, degree of stenosis, reference lumen diameter, and coronary resistance are the primary contributors to the uncertainty of computed FFR, accounting for 41.2%, 31.9%, and 24.6%, respectively. In particular, the interactive effect of the input variables on the uncertainty of the computed WSS is significantly higher than that on the uncertainty of the computed FFR.

Data availability statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author contributions

NH, KL, and S-WL contributed to the conception and design of the study. NH performed simulations and statistical analyses. NH and S-WL wrote the first draft of the manuscript. KL and S-WL critically revised the manuscript for important intellectual content. All authors contributed to the article and approved the submitted version.

Funding

This research was supported by grants from the Basic Science Research Program through the National Research Foundation of Korea (NRF), funded by the Ministry of Education (NRF-2020R1I1A30666174 & NRF-2020R1C1C1010316).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

1. Global health estimates: Life expectancy and leading causes of death and disability. Geneva, Switzerland: World Health Organization (2019). Available at: https://www.who.int/data/gho/data/themes/theme-details/GHO/mortality-and-global-health-estimates

2. Adjedj J, Xaplanteris P, Toth G, Ferrara A, Pellicano M, Ciccarelli G, et al. Visual and quantitative assessment of coronary stenoses at angiography versus fractional flow reserve: the impact of risk factors. Circ Cardiovasc Imaging. (2017) 10(7):e006243. doi: 10.1161/circimaging.117.006243

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Van Belle E, Rioufol G, Pouillot C, Cuisset T, Bougrini K, Teiger E, et al. Outcome impact of coronary revascularization strategy reclassification with fractional flow reserve at time of diagnostic angiography: insights from a large French multicenter fractional flow reserve registry. Circulation. (2014) 129(2):173–85. doi: 10.1161/circulationaha.113.006646

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Toth G, Hamilos M, Pyxaras S, Mangiacapra F, Nelis O, De Vroey F, et al. Evolving concepts of angiogram: fractional flow reserve discordances in 4000 coronary stenoses. Eur Heart J. (2014) 35(40):2831–8. doi: 10.1093/eurheartj/ehu094

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Balanescu S. Fractional flow reserve assessment of coronary artery stenosis. Eur Cardiol. (2016) 11(2):77–82. Disclosure: The author has no conflicts of interest to declare. doi: 10.15420/ecr/2016:24:2

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Berry C, Corcoran D, Hennigan B, Watkins S, Layland J, Oldroyd KG. Fractional flow reserve-guided management in stable coronary disease and acute myocardial infarction: recent developments. Eur Heart J. (2015) 36(45):3155–64. doi: 10.1093/eurheartj/ehv206

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Hwang D, Lee JM, Koo BK. Physiologic assessment of coronary artery disease: focus on fractional flow reserve. Korean J Radiol. (2016) 17(3):307–20. doi: 10.3348/kjr.2016.17.3.307

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Hamilos M, Muller O, Cuisset T, Ntalianis A, Chlouverakis G, Sarno G, et al. Long-term clinical outcome after fractional flow reserve-guided treatment in patients with angiographically equivocal left main coronary artery stenosis. Circulation. (2009) 120(15):1505–12. doi: 10.1161/circulationaha.109.850073

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Cook CM, Petraco R, Shun-Shin MJ, Ahmad Y, Nijjer S, Al-Lamee R, et al. Diagnostic accuracy of computed tomography-derived fractional flow reserve: a systematic review. JAMA Cardiol. (2017) 2(7):803–10. doi: 10.1001/jamacardio.2017.1314

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Pijls NH, De Bruyne B, Peels K, Van Der Voort PH, Bonnier HJ, Bartunek JKJJ, et al. Measurement of fractional flow reserve to assess the functional severity of coronary-artery stenoses. N Engl J Med. (1996) 334(26):1703–8. doi: 10.1056/nejm199606273342604

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Tonino PA, De Bruyne B, Pijls NH, Siebert U, Ikeno F, van’ t Veer M, et al. Fractional flow reserve versus angiography for guiding percutaneous coronary intervention. N Engl J Med. (2009) 360(3):213–24. doi: 10.1056/NEJMoa0807611

PubMed Abstract | CrossRef Full Text | Google Scholar

12. De Bruyne B, Fearon WF, Pijls NH, Barbato E, Tonino P, Piroth Z, et al. Fractional flow reserve-guided PCI for stable coronary artery disease. N Engl J Med. (2014) 371(13):1208–17. doi: 10.1056/NEJMoa1408758

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Dattilo PB, Prasad A, Honeycutt E, Wang TY, Messenger JC. Contemporary patterns of fractional flow reserve and intravascular ultrasound use among patients undergoing percutaneous coronary intervention in the United States: insights from the national cardiovascular data registry. J Am Coll Cardiol. (2012) 60(22):2337–9. doi: 10.1016/j.jacc.2012.08.990

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Zaleska M, Kołtowski Ł, Maksym J, Tomaniak M, Opolski M, Kochman J. Alternative methods for functional assessment of intermediate coronary lesions. Cardiol J. (2020) 27(6):825–35. doi: 10.5603/CJ.a2019.0027

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Gibson CM, Diaz L, Kandarpa K, Sacks FM, Pasternak RC, Sandor T, et al. Relation of vessel wall shear stress to atherosclerosis progression in human coronary arteries. Arterioscler Thromb. (1993) 13(2):310–5. doi: 10.1161/01.atv.13.2.310

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Eshtehardi P, McDaniel MC, Suo J, Dhawan SS, Timmins LH, Binongo JN, et al. Association of coronary wall shear stress with atherosclerotic plaque burden, composition, and distribution in patients with coronary artery disease. J Am Heart Assoc. (2012) 1(4):e002543. doi: 10.1161/jaha.112.002543

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Eshtehardi P, Brown AJ, Bhargava A, Costopoulos C, Hung OY, Corban MT, et al. High wall shear stress and high-risk plaque: an emerging concept. Int J Cardiovasc Imaging. (2017) 33(7):1089–99. Disclosures the authors have no other relevant affiliations or financial involvement with any organization or entity with a financial interest in or financial conflict with the subject matter or materials discussed in the manuscript apart from those disclosed. doi: 10.1007/s10554-016-1055-1

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Kumar A, Thompson EW, Lefieux A, Molony DS, Davis EL, Chand N, et al. High coronary shear stress in patients with coronary artery disease predicts myocardial infarction. J Am Coll Cardiol. (2018) 72(16):1926–35. doi: 10.1016/j.jacc.2018.07.075

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Xiu D, Karniadakis GE. The wiener–askey polynomial chaos for stochastic differential equations. SIAM J Sci Comp. (2002) 24(2):619–44. doi: 10.1137/s1064827501387826

CrossRef Full Text | Google Scholar

20. Xiu D. Numerical methods for stochastic computations: a spectral method approach. New Jersey, United States: Princeton University Press (2010). ISBN: 9781400835348. doi: 10.1515/9781400835348

21. Crestaux T, Le Maître O, Martinez J-M. Polynomial chaos expansion for sensitivity analysis. Reliab Eng Syst Saf. (2009) 94(7):1161–72. doi: 10.1016/j.ress.2008.10.008

CrossRef Full Text | Google Scholar

22. Hosder S, Walters R. Non-intrusive polynomial chaos methods for uncertainty quantification in fluid dynamics. Orlando, Florida: 48th AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition (2012).

23. Eck VG, Donders WP, Sturdy J, Feinberg J, Delhaas T, Hellevik LR, et al. A guide to uncertainty quantification and sensitivity analysis for cardiovascular applications. Int J Numer Method Biomed Eng. (2016) 32(8):e02755. doi: 10.1002/cnm.2755

CrossRef Full Text | Google Scholar

24. Sankaran S, Kim HJ, Choi G, Taylor CA. Uncertainty quantification in coronary blood flow simulations: impact of geometry, boundary conditions and blood viscosity. J Biomech. (2016) 49(12):2540–7. doi: 10.1016/j.jbiomech.2016.01.002

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Gashi K, Bosboom EMH, van de Vosse FN. The influence of model order reduction on the computed fractional flow reserve using parameterized coronary geometries. J Biomech. (2019) 82:313–23. doi: 10.1016/j.jbiomech.2018.11.008

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Sankardas MA, McEniery P, Aroney C, Holt G, Cameron J, Garrahy P, et al. Directional coronary atherectomy for lesions of the proximal left anterior descending artery: initial clinical results, complications and histopathological findings. Aust N Z J Med. (1995) 25(6):676–80. doi: 10.1111/j.1445-5994.1995.tb02852.x

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Marcus JT, Smeenk HG, Kuijer JP, Van der Geest RJ, Heethaar RM, Van Rossum AC. Flow profiles in the left anterior descending and the right coronary artery assessed by MR velocity quantification: effects of through-plane and in-plane motion of the heart. J Comput Assist Tomogr. (1999) 23(4):567–76. doi: 10.1097/00004728-199907000-00017

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Wieneke H, Haude M, Ge J, Altmann C, Kaiser S, Baumgart D, et al. Corrected coronary flow velocity reserve: a new concept for assessing coronary perfusion. J Am Coll Cardiol. (2000) 35(7):1713–20. doi: 10.1016/S0735-1097(00)00639-2

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Koo BK, Erglis A, Doh JH, Daniels DV, Jegere S, Kim HS, et al. Diagnosis of ischemia-causing coronary stenoses by noninvasive fractional flow reserve computed from coronary computed tomographic angiograms. Results from the prospective multicenter DISCOVER-FLOW (diagnosis of ischemia-causing stenoses obtained via noninvasive fractional flow reserve) study. J Am Coll Cardiol. (2011) 58(19):1989–97. doi: 10.1016/j.jacc.2011.06.066

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Klocke FJ. Coronary blood flow in man. Prog Cardiovasc Dis. (1976) 19(2):117–66. doi: 10.1016/0033-0620(76)90020-7

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Mintz GS, Popma JJ, Pichard AD, Kent KM, Satler LF, Chuang YC, et al. Limitations of angiography in the assessment of plaque distribution in coronary artery disease: a systematic study of target lesion eccentricity in 1446 lesions. Circulation. (1996) 93(5):924–31. doi: 10.1161/01.cir.93.5.924

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Yamagishi M, Terashima M, Awano K, Kijima M, Nakatani S, Daikoku S, et al. Morphology of vulnerable coronary plaque: insights from follow-up of patients examined by intravascular ultrasound before an acute coronary syndrome. J Am Coll Cardiol. (2000) 35(1):106–11. doi: 10.1016/S0735-1097(99)00533-1

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Tang B. Orthogonal array-based Latin hypercubes. J Am Stat Assoc. (1993) 88(424):1392–7. doi: 10.1080/01621459.1993.10476423

CrossRef Full Text | Google Scholar

34. Bertoglio C, Caiazzo A, Fernández MA. Fractional-step schemes for the coupling of distributed and lumped models in hemodynamics. SIAM J Sci Comp. (2013) 35(3):B551–75. doi: 10.1137/120874412

CrossRef Full Text | Google Scholar

35. Nguyen MT, Jeon BJ, Chang H-J, Lee S-W. Domain decomposition based parallel computing for multi-scale coronary blood flow simulations. Comput Fluids. (2019) 191:104254. doi: 10.1016/j.compfluid.2019.104254

CrossRef Full Text | Google Scholar

36. Choi HG, Choi H, Yoo JY. A fractional four-step finite element formulation of the unsteady incompressible navier-stokes equations using SUPG and linear equal-order element methods. Comput Methods Appl Mech Eng. (1997) 143(3):333–48. doi: 10.1016/S0045-7825(96)01156-5

CrossRef Full Text | Google Scholar

37. Sharma P, Itu L, Zheng X, Kamen A, Bernhardt D, Suciu C, et al. A framework for personalization of coronary flow computations during rest and hyperemia. Conf Proc IEEE Eng Med Biol Soc. (2012) 2012:6665–8. doi: 10.1109/embc.2012.6347523

CrossRef Full Text | Google Scholar

38. Feinberg J, Langtangen HP. Chaospy: an open source tool for designing methods of uncertainty quantification. J Comput Sci. (2015) 11:46–57. doi: 10.1016/j.jocs.2015.08.008

CrossRef Full Text | Google Scholar

39. Slager CJ, Wentzel JJ, Gijsen FJ, Schuurbiers JC, van der Wal AC, van der Steen AF, et al. The role of shear stress in the generation of rupture-prone vulnerable plaques. Nat Clin Pract Cardiovasc Med. (2005) 2(8):401–7. doi: 10.1038/ncpcardio0274

PubMed Abstract | CrossRef Full Text | Google Scholar