The study was conducted at Beijing Anzhen Hospital of the Capital Medical University and Shengjing Hospital of China Medical University. The inclusion criteria for patients in this study had complete clinical data and underwent CTA, quantitative coronary arteriography (QCA), transthoracic echocardiography, and invasive FFR measurement within 30 days. The exclusion criteria for this study were as follows: 1) Poor quality of CT images; 2) Unstable angina; 3) Prior myocardial infarction; 4) Prior percutaneous coronary intervention or coronary artery bypass grafting; 5) Diffuse coronary stenosis; 6) Left ventricular ejection fraction (LVEF) < 50%; 7) Severe valve disease; 8) Atrial fibrillation; 9) Severe microcirculation disturbance; 10) Allergy to contrast agents and vasodilators. Finally, 73 patients with stable angina were enrolled between August 2013 and April 2019 in this study. The quality of CT images in all patients was examined and assessed by two experienced radiologists. Cardiac output and left ventricular ejection fraction in all patients were measured and calculated by two experienced echocardiographers based on the structural characteristics of the heart (Li G. Y. et al., 2021).

The measurement of invasive FFR relied on three steps: 1) Adenosine (140 μg/kg/min) was administrated through intravenous infusion to induce maximum hyperemia of the coronary artery; 2) We obtained pressure waveforms of aortic pressure and distal arterial pressure using pressure wire measurement; 3) We further calculated the invasive FFR for the ratio of the mean pressure at a cross-Section 3 cm downstream of the stenosis (P_d) to the aortic pressure (P_a) at least three cardiac cycles (Taylor et al., 2013; Zhang et al., 2021a). Patient informed consent was waived due to the retrospective nature of the study.

Based on the patient’s CT image information, the 3D model of the coronary artery was reconstructed using the commercial software MIMICS (Materialise, Leuven, Belgium). The coronary arteries with diameters larger than 1 mm were reconstructed. Subsequently, the coronary reconstructed noise was removed using the Freeform tool (Artec 3d, Luxembourg). And the reconstructed coronary surfaces were further repaired and smoothed using the commercial software GEOMAGIC (Geomagic, Research Triangle Park, North Carolina). Then, the coronary centerline was identified to calculate the diameter and length of the vessel using the MIMICS software. GEOMAGIC was used to divide the reconstructed coronary surfaces into curved surfaces for hemodynamics simulation. Later, the inlet and outlet of the reconstructed coronary artery were cut into planes, and boundary conditions for pressure or mass flow were loaded using the SOLIDWORKS software (Dassault Systemes, Waltham, Massachusetts). Finally, the reconstructed coronary artery was imported into ANSYS CFX (ANSYS Corporation, Canonsburg, Pennsylvania) for CFD simulation. The process of 3D reconstruction of the coronary artery was shown in Figure 1.

Based on the available literature and clinical reports, coronary blood flow was 4% of cardiac output (Kim et al., 2010), and we can calculate coronary blood flow by transthoracic echocardiography measurement. Based on Poiseuille’s and the law of minimum energy dissipation, the allometric scaling law between coronary morphological and functional parameters was quantified using in vitro experiments (Taylor et al., 2013; Zhang et al., 2021b). We further established the distribution law of blood flow of the coronary branch on the basis of the form-follow-function scaling law described by Huo, in which the blood flow of the parent and daughter branches of the coronary artery was related to their respective effective diameters (Huo et al., 2012). These findings had a good consistency with the earlier in vitro experimental results of Zhou (Zhou et al., 1999). The allometric scaling law between blood flow and the diameter of the coronary branch followed a power-law relationship as shown in Eq. 1. Q s Q max = D s D max 7 3 (1)

In Eq 1, Q _max and D _max represent the blood flow and diameter of the parent branch, and Q _s and D _s represent the blood flow and diameter of the daughter branch.

Broadly speaking, vessel resistance was calculated using morphological parameters (cross-sectional area and length) according to Poiseuille’s law. The equation followed a power-law relationship, as shown in Eq 2. R = 8 π μ L A 2 (2) A = A i n l e t + A o u t l e t 2 (3)

In Eq 2, R, L, and A represent the resistance, length, and cross-section area of the coronary branch, and μ represents the dynamic viscosity of blood flow. In Eq 3, A _inlet and A _outlet represent the cross-section area at the inlet and outlet of the coronary branch.

Ohm’s law describes the relationship between current, resistance, and voltage in a circuit. The blood flow, resistance, and pressure drop along the coronary arteries represent the corresponding parameters, respectively (Li B. et al., 2021). So the pressure drop ( ∆ P ) along the coronary arteries was described by Eq 4. ∆ P = R Q s (4)

Systolic and diastolic blood pressure measurements were performed during the statistical process of clinical data. The calculation of mean blood pressure was described by Eq 5. P M B P = P s p + 2 P d p 3 (5)

In Eq 5, P M B P , P s p , and P d p represent mean blood pressure, systolic blood pressure, and diastolic blood pressure, respectively.

Based on the mean blood pressure and pressure drop along the coronary arteries, the inlet pressure of the microcirculation was calculated stepwise from the proximal to the distal end of coronary arteries using Eqs 4, 5. The microcirculation resistance was calculated using Eq 6. R i = P i Q s (6)

In Eq 6, R i and P i represent the resistance and inlet pressure of the microcirculation, and i represents the number of coronary branches, i = 1, 2, 3…16.

Extensive literature and clinical cases have shown that mean blood pressure and microcirculation resistance decrease as the coronary circulation changes from a resting to a hyperemia state (Wilson et al., 1990; Tang et al., 2020). Many pieces of literature have reported that the severity of coronary epicardial stenosis has no affected on minimal microvascular resistance (Asrnoudse et al., 2004; Fearon et al., 2004; Zhang et al., 2016). In this study, mean blood pressure was reduced by 12% (Wilson et al., 1990), and microcirculation resistance was taken to be 0.23 times the resting state to mimic the hyperemia state caused by the effect of adenosine (Wilson et al., 1990). P 0 = P M B P − 0.12 P M B P (7) R j = 0.23 R i (8)

In Eqs 7, 8, P 0 and R j represent mean blood pressure and microcirculation resistance at hyperemia state, respectively, i = j = 1, 2, 3…16.

Eqs 9, 10 were coded in Fortran by using the user-defined function (UDF) of ANSYS CFX, running on an HP Z8 workstation. Then, the fluid dynamics analysis of the coronary artery was updated with the under-relaxation scheme as formulated in Eqs 9, 10. Eq. 9 until the sum of the pressure gradient of the epicardial coronary and microcirculation resistance matched the inlet pressure at hyperemia, and Eq. 10 until the target residual of the pressure gradient at the outlet of the coronary branch was determined to be 1e-4. Finally, we identified the blood flow at the outlet of the coronary branch under the hyperemia state. Q j , n e w = 1 − α Q j , o l d + α P j , o l d R j (9) P j , n e w = P j , o l d + α P 0 − P j , o l d + Q j , n e w R j (10)

In Eqs 9, 10, P j , o l d , P j , n e w , Q j , o l d , and Q j , n e w represent the pressure and blood flow before and after iterative calculation at the outlet of each coronary branch, respectively, and α represents the under-relaxation factor.

Blood flow was modeled as a Newtonian fluid. The blood flow state was steady, and the properties of the arterial walls were set to non-slip rigid (Zhang et al., 2014). The density and viscosity of blood flow were set at 1,050 kg/m³ and 0.0035 Pa s, respectively (Zhang et al., 2021a). The mesh of the coronary geometry model discretized the computational domain into tetrahedral elements. In this study, the CPU of an HP Z8 workstation is a Dual Intel Xeon Silver 4,210 processor, and the memory of the workstation is 128 GB. The mesh of the geometries was generated by using nonstructural tetrahedron elements. The maximum grid size was 0.23 mm based on the grid independence test. Then, we used the 3D N-S function to calculate the pressure and blood flow field of the coronary arteries.

We determined the pressure at the coronary inlet based on the effect of adenosine on mean blood pressure. Then, we can quickly identify the blood flow of each coronary branch by compiling a user-defined function at the outlet of the coronary branch. Later, combining the mean blood pressure as the inlet boundary condition and the blood flow at the outlet of each coronary branch to compile a user-defined function as the outlet boundary condition, we implemented a CFD simulation of the coronary arteries using ANSYS CFX. Finally, we extracted the pressure field of the coronary arteries to calculate FFR_U, the ratio of the pressure at a cross-Section 3 cm downstream of the stenosis to the aortic pressure, as shown in Figure 2.

Based on the mean values of aortic diastolic pressure, myocardial mass, and heart rate of the patient, we calculated a patient-specific coronary blood flow rate according to the empirical formula of total coronary blood flow rate. Based on the volume or diameter from CCTA, we achieved patient-specific distribution of the blood flow of the left and right coronary artery. Based on the blood flow of the left and right coronary artery and the distribution rule, we achieved a patient-specific coronary blood flow rate at each terminal branch. Combining the mean value of aortic diastolic pressure as the inlet boundary condition and the coronary blood flow rate at each terminal branch as the outlet boundary condition, we obtained the patient-specific boundary conditions of the fluid dynamics analysis. Based on the fluid dynamics analysis of the coronary artery, we extracted the pressure field of the coronary artery, and we calculated the FFR_D as the ratio of the mean pressure at a cross-Section 3 cm downstream of the stenosis to the mean arterial pressure.

In this study, to identify the blood flow at the outlet of the coronary branch and accelerate the calculation convergence, a user-defined function (UDF) was compiled to identify the blood flow and pressure at the outlet of each coronary branch. Then, a user-defined function was used to integrate the interaction of stenotic resistance, microcirculation resistance, and inlet aortic pressure to identify coronary blood flow to quantify FFR_U calculation based on boundary conditions of clinical statistics of hyperemia state. In contrast, as for FFR_D calculation, the blood flow at the outlet of each coronary branch was distributed step by step along the proximal to the distal blood flow direction based on the hyperemia state, and it ignored the impact of coronary stenosis on the distribution of coronary blood flow.

Clinical data analysis included clinical statistics, CTA, QCA, transthoracic echocardiography, and invasive FFR. Continuous and categorical were shown as mean, frequency, and/or percentage, respectively.

To evaluate the diagnostic accuracy of the novel method for calculating FFR_U, we used the novel method to calculate FFR_U for 86 vessels in 73 patients and then compared these data with those derived from existing methods and invasive FFR.

Bland-Altman plot with 95% confidence intervals (CI) was used to evaluate the consistency of the novel method for FFR_U calculation and the existing method (FFR_D) with invasive FFR. Based on the reference value of invasive FFR≤0.8 for the diagnosis of lesion-specific myocardial ischemia, we adopted the metrics of receiver operating characteristic curves (AUC) with 95% (CI), sensitivity, specificity, positive predictive value (PPV), and negative predictive value (NPV), to assess the diagnostic accuracy of FFR_U calculation and existing methods (CTA≥50% stenosis, QCA≥50% stenosis and FFR_D).

The baseline demographics of the patients were shown in Table 1. Of these, more than half of the patients were men (69.86%), and the mean patient age was 59 ± 16 years.

QCA and invasive FFR measurements were successfully performed in 86 vessels. The coronary arteries of most patients (79.45%) were right-dominant patterns by QCA. Among the 86 vessels, more than half of the lesions (68.6%) occurred in the left anterior descending artery (LAD), and 67 vessels (77.91%) had luminal stenosis≥50%. Of these, only 32 vessels (37.21%) showed significant ischemia (FFR≤0.8), as shown in Table 2.

The calculation of FFR_U and FFR_D was successfully performed on 86 vessels using a HP Z8 workstation. The average computational time for FFR_U was significantly reduced, taking only 5 min per simulation, which was 2-fold less than the FFR_D method. Table 3 shows the distribution of FFR_D and FFR_U calculated by the two approaches in 86 vessels. Regarding the calculation of FFR_D, the numbers 27 and 5 represent that 27 vessels were calculated as having FFR_D≤0.8 and 5 vessels were FFR_D>0.8 among the vessels with invasive FFR≤0.8 in 32 vessels. Concerning FFR_CT calculated by the proposed approach, the numbers 4 and 50 represent that 4 vessels were calculated as having FFR_U≤0.8, and 50 vessels were FFR_U>0.8 among those with invasive FFR>0.8.

Bland-Altman analysis was used to test the consistency of FFR_U and FFR_D with invasive FFR, respectively. The mean difference in FFR_U-FFR (0.006) for all vessels was less than FFR_D-FFR (0.018). The 95% CI of FFR_U and FFR_D was [-0.073, 0.085] and [-0.191, 0.228], respectively, and most of the data fell within the interval, indicating that FFR_U and FFR_D were in good agreement with invasive FFR, as shown in Figure 3.

The AUC of the receiver operating characteristics curve analysis at the lesion level for CTA, QCA, FFR_D and FFR_U were 0.62 (95% CI: 0.51–0.74), 0.67 (95% CI: 0.56–0.79), 0.85 (95% CI: 0.76–0.94), and 0.93 (95% CI: 0.87–0.98). At the patient level, the AUC was 0.61 (95% CI: 0.48–0.74) for CTA, 0.65 (95% CI: 0.53–0.77) for QCA, 0.83 (95% CI: 0.74–0.92) for FFR_D, and 0.92 (95% CI: 0.89–0.96) for FFR_U, as shown in Figure 4. The AUC at the lesion level demonstrated the diagnostic accuracy of FFR_U was higher than that of CTA, QCA, and FFR_D. Similar results were observed at the patient level with FFR_U (AUC, 0.92) compared with CTA, QCA, and FFR_D.

Table 4 shows that the per-vessel level sensitivity analysis of CTA, QCA, FFR_D, and FFR_U were 71.87%, 93.75%, 84.38%, and 93.75%; specificity of 51.85%, 40.74%, 85.19%, and 92.59%; PPV of 46.94%, 48.39%, 77.14%, and 88.24%; and NPV of 75.68%, 91.67%, 90.2%, and 96.15%. The diagnostic accuracy of the four metrics of FFR_U was higher than that of CTA, QCA, and FFR_D. Similar phenomena were observed at the patient level in FFR_U compared with CTA, QCA, and FFR_D.

Although the valuable information derived from our novel method, improved the diagnostic accuracy and reduced the computational time for FFR_U calculation, several limitations are notable. First, 73 patients undergoing CTA, QCA, and invasive FFR were enrolled from two central databases. The diversity and number of patients were relatively small. So the diversity and number of patients are enrolled from multiple centers in our future work. Second, the compliance of the epicardial coronary artery was neglected from resting to hyperemia state, and some studies had reported that the compliance of the epicardial coronary artery had almost no difference in the changes of coronary blood flow and pressure (Zeng et al., 2008; Zhang et al., 2014). Third, mean blood pressure was reduced by 12% to mimic the mean blood pressure changes from resting to hyperemia state in all patients. It ignored the effect of patient-special mean blood pressure changes on FFR_U calculation, and Zhang et al. reported that the mean blood pressure changes had almost no difference for FFR_CT calculation (Zhang et al., 2020). Finally, the steady-state numerical simulation employed in this study calculated the FFR_U, ignoring the effect of the pulsatile flow characteristics on FFR_U calculation. And the pulsatile flow was reported to be less important in FFR_CT calculation (Zhang et al., 2014).