Deep learning models typically need many samples for training. This is because a model will need to learn the weights for numerous parameters and uses gradient-based optimization to converge to a global optimum. When the model attempts to fit data patterns, it may learn random noise if a dataset is small. This causes what is known as overfitting and has as a consequence the inability of the model to generalise well. We have conducted a practical experiment to show that as the training data size increases, model loss decreases (refer to Supplementary Figure S9).The need for larger realistic datasets, and particularly those that are more private/secure (obfuscate real patient data), have driven studies into the generation of synthetic data. It is important that the synthetic data follows the underlying distribution of the real data and captures correlations between features in order to be plausible. As yet, there are no public datasets of arterial CFD that can be used for training. Since in our study we have only 50 patients in scope, it is essential to augment the training set by generating synthetic data.

In machine learning, and particularly in the case of images, there are many simple approaches that can be utilised to create a larger sample. These methods are referred to as data augmentation and involve transforming an image via shifting (horizontally or vertically), flipping (horizontally or vertically), rotating (clockwise or anti-clockwise), zooming in or out, and changing the brightness (Chlap et al., 2021). These methods can create more variation in the data and artificially expand the dataset.

The use of flipping as an approach of data augmentation is commonly seen in medical image processing (Cossio, 2023; Perez and Wang, 2017). Nishio et al. applied random 10° rotations, translations, and horizontal flipping to chest X-ray radiography (CXR) images in order to improve the accuracy of their CNN model for COVID-19 diagnosis (Nishio et al., 2020). Schmid et al. (2023) propose the use of statistical shape and intensity models (SSIM) to generate augmented CT images of hip bones including flipping, shifting and rotation. Although the reversal of inlet and outlet may not be physiologically representative, we examined the ratio of inlet and outlet radii, which is 0.86 ± 0.12 and suggesting the tapering is not substantial. Flipping remains a crucial step to avoid overfitting and creating a more general model.

Inspired by these approaches, we modified the geometry of the coronary artery and running the CFD analysis in ANSYS. In that way, we were able to create realistic synthetic data by making 3 modifications to the original patient data. After implementing the following modifications listed, 550 synthetic patients were created that were used for training.

• Reverse inlet and outlet boundaries: The inlet and outlet boundaries were reversed (i.e., flipping was applied) for the original patient data whilst maintaining the original patient inlet velocity. Thus, this created 1 new dataset per patient and overall 50 synthetic datasets.

In order to train the model, we split the synthetic dataset generated. Of the total 600 patient dataset (50 real + 550 synthetic), 80% was used for training (40 real + 440 synthetic generated from those) and 20% was used as a test set (10 real patients). The validation set, still used as part of the training process for hyperparameter tuning, is automatically created once the training begins and comprises 20% of the training set (96 real/synthetic data mixture). It is important to note that the train/test split is completed at the real patient level, such that only synthetic data associated with the 40 real patients (440 patients) can be used for training. This is done in order to prevent leakage of information from the training set into the test set. In other words, the model has not seen any information similar to the test patients, not even their synthetic data. Once all data creation and pre-processing was complete, the neural network architecture was designed.

We adopted and modified the neural network architecture used in Gharleghi et al. (2020) by removing the concatenation with global features such as bifurcation angle, which is not relevant for our use case. We also added a spatial dropout layer in order to reduce overfitting of the model. The final neural network can be categorised as a U-Net and had an architecture shown in Figure 2. This network uses average pooling to scale the input data down to half of its resolution. This is done twice, leading to a quarter representation of the original dataset. The data is then passed through 2 convolutional layers where a 3 × 3 convolution was applied, followed by up-sampling and concatenation with the higher resolution data. Adam optimizer was then used in addition to a 10% spatial dropout as a regularization technique. The activation function selected was ReLU such that the output is constrained to a positive value.

As the loss function, we selected the mean squared error (MSE) although other similar studies in the literature, including Gharleghi et al. (2020) and Jordanski et al. (2018), typically adopt the mean absolute error (MAE). The MSE was calculated as per Eq. 3, where n represents the number of vessels, τ i represents the observed values (i.e., CFD-based WSS), and τ ̂ i represents the predicted values (i.e., DL-based WSS). M S E = 1 n ∑ i = 1 n τ i − τ ̂ i 2 (3)

In the case of our patient dataset, we are investigating patients with stenosis where WSS is significantly elevated. Our choice of MSE is to better predict the peak WSS, as it incurs a higher penalty when mismatched. To highlight the importance of this, we also trained the model using MAE for comparison. The model was trained for 1,000 epochs. The computational resource utilised was an NVIDIA T4 Tensore Core GPU with 52GB RAM. Training time is approximately 40 min on average. Increasing infrastructure resources can lead to faster model training and inference.

The deep learning model was assessed for WSS prediction accuracy. We implement k-fold cross-validation with 5 folds, where the DL-based prediction of WSS for the 10 patient data left aside in Section 2.3 are compared against their original CFD-based WSS, over 5 repeated training ∼ validation steps.

DL models are typically considered as a black box models, and the associations between particular inputs and outputs have not been formally investigated. However, as our work also aims to understand the mechanistic insight behind WSS profiles, explainability (xAI) techniques were used to quantitatively assess feature importance on a trained model. Many techniques exist that have been used in DL models such as Shapley values (Lundberg et al., 2017), however as this approach is computationally expensive, and the model internals are available, we opted for the integrated gradients approach.

Integrated Gradients is a common technique for explaining differentiable models such as neural networks. It is based on two key properties: implementation invariance and sensitivity. It works by creating a straight path between a reference point (typically zeros) and the inputs to the model. By partitioning and interpolating the path, the model can compute predictions at the different partition points. The gradient information at the points of this path are calculated, making it computationally efficient (Holzinger et al., 2022). The intuition is that the gradient in the points where a model’s predictions have flattened out are zero and therefore do not contain information that contributes to the explanation. With this process, the significance of each input feature in the prediction of output (i.e., WSS) can be estimated. In order to apply this method on our deep learning model, we utilise the Innvestigate library1 (Alber et al., 2018).

The results of the MSE on the test set are shown in Table 3. The mean absolute error (MAE) is also displayed to show how different the absolute value of the prediction is from the gold-standard, CFD-based WSS for the same patients. In the literature, the normalised mean absolute error (NMAE) is typically reported, thus we include this metric in order to provide comparison with the other state-of-the-art methods. The NMAE can be calculated by either dividing the MAE by the mean or by the difference between the maximum and minimum WSS. For this paper, we divided by max-min, where the difference in WSS is 33.33 Pa.

As previously mentioned, the model training with our data took on average 40 min. This means that it took less than 2.5 s per epoch; not impractically long in the current form, yet it can be accelerated even further in the future with additional infrastructure resources. Model training results are presented in the Supplementary Table S4. Prediction time was approximately 35 milliseconds, which is many orders of magnitude faster than the CFD processing time, which ranges from 20 min to approximately 3 h in our cases. Furthermore, the model was able to predict the WSS on the test set to a NMAE of 6.03% with a standard deviation of 0.47% when normalised by the difference between the maximum and the minimum WSS.

In order to show correlation between the CFD-based WSS and the DL model prediction, a Bland-Altman plot was generated for the test set and shown in Figure 3. Here, the minimum, maximum, and mean predominant WSS in every 3-mm vascular segments were used. Briefly, predominant WSS is defined as a moving average of WSS within a window around a point of interest, which has a size of 90 ° circumferentially and 3 mm longitudinally. These metrics were originally proposed to compensate image data uncertainty, associated with motion of an intravascular imaging catheter during a heartbeat, in the creation of a 3D vascular model. They have been used in previous studies examining the value of WSS in predicting adverse clinical events (Stone et al., 2018), including ours (Tufaro et al., 2021), and shown to be practical and effective. The means and differences are calculated between the CFD-based and DL-predicted WSS and plotted along with the ± 1.96 SD lines. Figure 3 shows the values for all the patients. The majority of points lie between 95% CI, demonstrating a correlation between the CFD-based and DL-based WSS values, with the bias of maximum 0.38 Pa.

By scanning along the slices of data in the direction depicted by Figure 1, we can plot the WSS profile for a particular patient. Figure 4 illustrates the WSS profile for a patient in more detail, the DL-predicted WSS values are plotted over the wall surface from proximal (left) to distal (right) and across the circumference of each cross-section.

The plotted data display the output of the model trained with MAE as the loss (red) and the MSE as the loss (green), along with the ground truth CFD-based WSS that are displayed in blue. It is clear that the model trained by MSE predicts WSS that fits better to the ground truth. Although the MAE result had a low overall error, it did not capture the local variations of WSS as much as the MSE-based training, especially in the areas of stenosis. This justifies the choice of MSE loss for model training.

Among the 10 test cases, we performed a qualitative assessment using 3D WSS surface plot in order to illustrate in more detail the patient-level WSS prediction results. Two of those vessels with the most accurate and inaccurate WSS predictions are presented in Figure 5, and similar plots for the full 10 test cases are presented in the Supplementary Figure S10. It is apparent that in both cases the WSS patterns are generally captured by the DL model, but there are differences after looking the data in detail, e.g., in the distal region of Patient 7 where the vessel is more tortuous. This trend is observed in the other models shown in the Supplementary Material. Patient 7 appears to be the most curved/tortuous and having multiple stenosis sites, which could have given the difficulty in the prediction of WSS by the DL model. Additionally, we observed that the prediction by the DL model underestimates WSS in the stenosis region and tends to overestimate otherwise.

Figure 6 shows the value of the 17 features based on the integrated gradients estimations. The integrated gradients are calculated per each data point on the vascular wall surface and give an evaluation on how the input features contribute to the prediction. In order to conduct an overall assessment of the model, the samples need to be aggregated; the mean feature importance across all points along a vessel surface, and the mean across all test patients are summed to have the final results. Our findings indicate that the radius and inlet velocity are the most important features while the circumferential coordinate θ , the components of centreline tangential and normal vectors t x and n y were found to be the less contributing features. The positive and negative signs indicate whether the relationship between the input feature and the output feature (WSS) is positively correlated or inversely proportional.

The top six important features were then selected (velocity, radius, theta, curvature, tangential vector in Y direction and normal vector in Z direction) and the model was retrained with MSE as the loss function. This resulted in an MSE of 8.89 P a 2 , an MAE of 2.28 Pa, NMAE of 6.84%. As expected, the results are comparable to using the full set of features available. This is because the model is placing a low weight on the remaining features, as such they do not have a strong impact on the final prediction. Nonetheless, utilizing the full feature set (with MSE loss) results in an average MSE of 8.6 P a 2 , an MAE of 2.01 Pa, and NMAE of 6.03%; thus still outperforming the experiment with the subset of features.

As mentioned earlier, training of DL models tend to require large data sets. For example, Cho et al. (2016) investigated the impact of increasing the number of training samples of axial Computed Tomography (CT) images for classification into six anatomical classes. The results showed that an increase in the number of samples from 5 to 200 improved model accuracy from 8% to 95% (Cho et al., 2016). Similar findings have been also reported in other studies that highlighted the relationship between a larger data set and improved model performance, not only for machine learning but especially for deep learning (Sarker, 2021; Halevy et al., 2009; Zhang et al., 2019).

Gharleghi et al. (2020) generate a synthetic dataset out of 127 patients by modifying the bifurcation angle geometry and, as a result, obtain 3,302 synthetic patients. In our paper, to augment our dataset, we took a simple approach by flipping the geometry in combination with changing the input velocity. This results in 550 synthetic patients that were then used to train a DL model that can predict WSS with an accuracy that is similar to previous reports. An advantage of the proposed method is that it does not require creation of vessel geometry and computational mesh generation. These are time-consuming processes requiring up to 3 h to estimate the WSS patterns of a single new synthetic “model.”

The training dataset size in this study, is deemed effective in terms of the trade-off between accuracy and training time, based on our preliminary investigation using idealised curved tubes with a stenosis (Supplementary Material). However, it is likely that a larger training datasets would provide more accurate estimation of WSS patterns as real coronary anatomy varies between patients and vessels and this critically determines WSS. To mitigate this, acquisition of real patient data is ideal, but considering the largest number of patients reported in a CFD analysis is approximately 500 (Stone et al., 2012), synthetic data will still be necessary. Additionally, patients’ data are often associated with confidentiality concerns, which may prevent the dataset from being shared with the research community. Thus, future studies should try to combine the inclusion of larger clinical datasets and the creation of synthetic data from these data.

A common method for synthetic data generation of geometries in biomedical engineering is statistical shape modelling (SSM). Three-dimensional (3D) SSM facilitates the study of shape variability and allow the creation of new geometries with a wide range of variability (Alemany et al., 2019). Such statistically varied models offer an opportunity for experiments by exploring how changes to the shape geometry impact other factors (Sarker, 2021), which has been utilised in combination with CFD (Lamecker and Zachow, 2016). For example, Bruse et al. (2017) show that SSMs can aid to improve cardiac device development by modelling hemodynamic and geometric boundary conditions in cardiac anatomy. These models have been shown to improve both the efficiency robustness and value of synthetic patient data (King et al., 2019).

As previously noted, the MSE was used as a loss function as it facilitates a better fitting of the model prediction to WSS peaks due to the larger penalty in these areas. This is ideal for our vessels as there are sharp increases of WSS in the areas of a tight stenosis, reaching up to 70 Pa, whereas normal WSS range is reported to be 1–7 Pa (Malek, 1999). The WSS prediction error from the model trained on MAE turned out to be 2.10 ± 0.28 Pa, with a NMAE of 2.51 ± 0.33% (Table 3). This value on average is comparable to that calculated when MSE is the loss (2.01 ± 0.18 Pa MAE, p = 0.60 vs. MAE of the models trained on MAE loss). Therefore, as far as global error metrics are used to evaluate model accuracy, a difference is not evident and that may be the reason why most of the previous works utilised the MAE for simplicity. Additionally, when models are trained on MAE, the standard deviation of error (MAE) tends to be larger than MSE-based training (0.28 Pa vs. 0.18 Pa) though that is not statistically significant (p = 0.19 for F-test of training outcome MAE variances in Table 3).

Bland-Altman plots were generated for the minimum, maximum and mean predominant WSS to assess the difference between the CFD-based and DL-predicted estimations (refer to Figure 2). The limits of agreement range from −3.59 to 3.37 Pa. These ranges may have implications when WSS is low, but in the case of prediction of WSS in stenotic vessels where WSS in our dataset reached over 70 Pa, this result shows that the DL performed well and may have a value for real time computation of the WSS. The heterogeneous distribution of the WSS data in the test set underscores the challenge of the prediction task. Of note at the lesion site, which is the most clinically relevant segment, the range of difference between DL and reference standard was −4.04–3.96 Pa and this is comparable to that reported by Tufaro et al. (2022) (−4.1–5.7 Pa) who compared the estimations of two CFD-based approaches (one conventional CFD analysis performed using the ANSYS software a CFD analysis performed by a dedicated software CAAS Workstation WWS) in models reconstructed by two different software (Medis vs. Pie Medical).

An example of patient-level WSS prediction in Figure 4 shows that the error values are comparable for both the MAE trained and MSE trained models (2.85 Pa vs. 2.55 Pa respectively) however the WSS profile of the MSE approached better the reference standard in areas of stenosis. Nonetheless, in order to identify the areas that may be more difficult for prediction, a difference 3D plot is generated in addition to 3D visualization for DL-predicted and CFD-based coronaries (refer to Figure 5). Although the DL prediction is able to capture the overall WSS distribution, it can be seen that in diseased segments, the prediction is typically lower than the CFD-based results and the error is higher. Nevertheless, our prediction accuracy is similar to other existing models that range from 1.6% to 10.1% NMAE (normalised by the maximum) (Su et al., 2020; Suk et al., 2022), despite the fact that we included distinct stenoses.

Explainability is one of the current key advancements in machine learning models, moving forward from the use of black-box models (Belle and Papantonis, 2020). It is important to explain why a model has made a particular prediction before applying this in clinical practice. In typical supervised learning, e.g., linear regression models, model coefficients associating the input features to prediction output can be used to explain the importance of a specific input feature. However, this is more challenging with DL approaches, as they are typically black-boxes. In this paper, we introduced engineered features based on domain knowledge in order to improve model prediction. Feature engineering plays a key role in model prediction (Heaton, 2016). Previous works have demonstrated a relationship between key features and WSS, such as curvature and velocity (van Oijen, 2005).

Our results indeed indicate that velocity carries significant importance, following the vessel radius. This makes a mechanistic sense since WSS on a tube wall is theoretically determined by the flow rate and radius assuming Newtonian fluid and parabolic velocity profile ( τ = 4 μ Q π R 3 ) . Curvature as Menger curvature is not shown as significant, but the components of centreline tangential and curvature vectors t x ⃗ and n Y ⃗ were found to be relatively significant features, along with the circumferential position θ , indicating that the curved vascular geometry and relative location of the wall play a role in WSS distributions. Further research is needed towards this direction and better understand WSS distribution. For instance, different features may play a different role in disease-free and stenotic segments.

This study was an initial attempt to predict WSS in patient-specific stenotic vessel geometries, hence it has limitations. First, the number of patient-specific vessels was limited to 50, of which only 40 were used for training. This also limited the range of stenosis degree and inflow velocity as presented in Table 1. Although the results indicated comparable predictive capability to the literature, this can be improved by taking advantage of the full set of patients ( n = 293 ) included in our previous study (Tufaro et al., 2021). A larger number of patients with clinical outcome data will facilitate further assessment of DL-based WSS prediction more towards its clinical utility. Furthermore, the current output of the model is the magnitude of the WSS; not a WSS vector which is of high interest in hemodynamics modeling in combination with their transient behaviour. Steady-state flow is also an assumption that was made; by extending the work to unsteady flow, time-resolved WSS behaviour can be learned. Addressing these points will allow the study to include multidirectional WSS metrics, such as Oscillatory Shear Index, transverse WSS and WSS topology, etc. (Morbiducci et al., 2020; Hoogendoorn et al., 2019; Kok et al., 2019).

Future work is also needed to improve model performance and optimise DL architecture. The U-net approach used in this study is a well established approach, but it has been developed in 2015 (Ronneberger et al., 2015) and since then several other DL approaches have been introduced in cardiovascular research such as physics-informed neural networks (PINNs) (Raissi et al., 2019) and graph-neural networks (GNNs) (Scarselli et al., 2009). These DL architectures may enhance model performance.

Improvements to data augmentation in terms of quality and quantity can also be introduced. For instance, additional features may be engineered such as torsion in order to enrich the data set and, as a result, the data quality. Moreover, a finer level of granularity such as more slices per patient may be beneficial for achieving more realistic results, in addition to generating a larger cohort of synthetic data. Finally, a further deep dive into explainability should be considered, in particular in the areas of stenosis compared to the rest of the artery.