The comparison of this study with our previously published work [20, 31] can be found in Table S1. The reproducibility of image evaluation results between the two radiologists was excellent for the DSA score (ICC, 0.978; 95% CI, 0.961–0.988; p < 0.001) and CTA score (ICC, 0.992; 95% CI, 0.990–0.994; p < 0.001).

A total of 56 patients (38 men; median age, 69 years (IQR, 64–78 years)) with 56 lower limbs were included in this study. The patients were divided into two groups according to the DSA runoff score. The group with a low DSA score (≤ 7) included 36 patients (26 men; median age, 68 years (IQR, 63–78 years); mean DSA score, 3.3 ± 2.6; mean CTA score, 6.1 ± 4.0), and the group with a high DSA score (> 7) included 20 patients (12 men; median age, 70 years [IQR, 64–79 years]; mean DSA score, 12.3 ± 3.0; mean CTA score, 11.8 ± 5.1).

The clinical characteristics of the enrolled patients are shown in Table 1. The height of the patients was 1.7 m (IQR, 1.6–1.7 m), the weight was 68 ± 12 kg, and the BMI was 25.2 kg/m2 (IQR, 22.2–27.4 kg/m2). In the study cohort, 15 patients (27%, 15/56) had coronary heart disease, 41 (73%, 41/56) had hypertension, 31 (55%, 31/56) had diabetes, 7 (13%, 7/56) had hyperlipidemia, 31 (55%, 31/56) had a smoking history, and 18 (32%, 18/56) had a history of drinking alcohol. As shown in Table 1, there was no statistical difference in the clinical features between the two groups.

Table 2 presents the results of the univariate analysis of the 20 histogram features extracted from the lower leg muscles. Seven features, which measured the intensity of CT values from different aspects, showed significant differences between the two groups (DSA score ≤ 7 vs. DSA score > 7): 10th percentile (24.2 vs. 19.0, p < 0.001), 90th percentile (64 vs. 61, p = 0.006), energy (6.8 × 108 vs. 5.2 × 108, p = 0.04), mean (44.6 vs. 39.5, p < 0.001), median (45 vs. 39, p < 0.001), mode (47 vs. 40, p = 0.001), and root mean squared (47.2 vs. 43.1, p < 0.001). All these features showed higher values in the mild PAD group.

In addition, there were significant differences between the groups in the following seven features that indicate the dispersion or deviation of CT values (DSA score ≤ 7 vs. DSA score > 7): CV (35.6 vs. 41.0, p < 0.001), IQR (20 vs. 21, p = 0.009), mean absolute deviation (12.4 vs. 12.9, p = 0.04), robust mean absolute deviation (8.5 vs. 8.9, p = 0.01), skewness (0.0 vs. 0.3, p = 0.002), standard deviation (16.0 vs. 16.7, p = 0.03), and variance (256 vs. 280, p = 0.03). These features showed higher dispersion in the severe PAD group. For the other six histogram features, no differences were observed between the two groups. Details of the statistical results are shown in Table 2.

Figure 3 shows the analysis results of the 75 texture features. Figure 3a presents the p values of these features. The values in the red, yellow, and blue areas represent p < 0.01, 0.01 < p < 0.05, and p > 0.05, respectively. Figure 3b was the feature names corresponding to Fig. 3a. As shown in Fig. 3, there were significant differences between the two groups in 31 texture features. Results of the shape features can be found in Table S2 and Supplementary Materials. Compared to the severe PAD group, the mild PAD group had larger muscle volume and average area (volume (cm3), 1084 vs. 876, p = 0.01; area (cm2), 32 vs. 25, p = 0.04).

Univariable analysis of texture features of the lower leg muscles. Texture features contain 75 statistics in 5 categories: 24 gray level co-occurrence matrix (GLCM) features, 14 gray level dependence matrix (GLDM) features, 16 gray level run length matrix (GLRLM) features, 16 gray level size zone matrix (GLSZM) features, and 5 neighboring gray tone difference matrix (NGTDM) features. a The p values of the statistical analysis of the 75 texture features. The values in the red, yellow, and blue areas represent p < 0.01, 0.01 < p < 0.05, and p > 0.05, respectively. b Feature names corresponding to (a). The numbers before the feature names in b corresponding to the leftmost line number in a

Considering the collinearity between muscle features, we further selected 45 features (14 histogram features and 31 texture features) with p < 0.05 using LASSO regression. As shown in Fig. 4, five representative image features were finally identified. Figure 4a–c describes the process of feature selection, and Fig. 4d shows the feature names and weight coefficients of the final selected features.

Feature selection based on least absolute shrinkage and selection operator (LASSO) regression. a Forty-five features (14 histogram features and 31 texture features) with p < 0.05 selected by univariable analysis. b The trend graph of the mean square error (MSE) with different λ (Lamda) during cross-validation. λ is an important parameter of LASSO regression that is usually adjusted by cross-validation to find the optimal value. The red dots represent the average values of the MSE. The blue error bars represent the standard deviation of the MSE. The black dotted line indicates the best value of λ. c The convergence graph of the weight coefficients of the features under different λ values. Each convergence line corresponds to a feature, and the color of the line matches the color before the feature name in a. As shown in b and c, the MSE is minimized (0.21 ± 0.07) at λ = 0.044 (the black dotted line), where five representative features were finally identified (weight coefficient ≠ 0). d Feature names and weight coefficients of the five selected features. GLCM gray level co-occurrence matrix, GLDM gray level dependence matrix, GLRLM gray level run length matrix, GLSZM gray level size zone matrix

Table 3 presents the multivariable analysis results of the two LRMs (LRM-I and LRM-II). As shown in Table 3, histogram_10Percentile (OR, 0.28; 95% CI, 0.13–0.61; p = 0.001) and gldm_DependenceNonUniformityNormalized (OR, 0.51; 95% CI, 0.27–0.98; p = 0.04) retained significant difference in LRM-I. CTA runoff score (OR, 3.27; 95% CI, 1.42–7.53; p = 0.006), histogram_10Percentile (OR, 0.33; 95% CI, 0.14–0.77; p = 0.01), and gldm_DependenceNonUniformityNormalized (OR, 0.51; 95% CI, 0.25–1.03; p = 0.06) were eventually selected in LRM-II.

According to the coefficients of variables and the constant in Table 3, we establish logistic regression equations of the two models: LRM-I, \(Logit\left(P\right)=-1.29\times }1-0.67\times }2-0.85\); LRM-II, \(Logit\left(P\right)=1.18\times \mathrm-1.11\times }1-0.68\times }2-0.90\). In the two equations, Logit (·) represents logit transformation, P represents the output probability of the LRMs, Feature1 represents the value of histogram_10Percentile, and Feature2 represents the value of gldm_DependenceNonUniformityNormalized.

Table 4 shows the performance of the CTA runoff score and the LRMs constructed with lower leg muscle features. Both LRMs were significant (Omnibus test, p < 0.001) with a high goodness of fit (Hosmer–Lemeshow test, p > 0.05). Compared to the CTA score (AUC, 0.81; sensitivity, 75%; specificity, 75%; accuracy, 75%), LRM-I achieved a better predictive performance (AUC, 0.84; sensitivity, 80%; specificity, 81%; accuracy, 80%). LRM-II achieved the best results (AUC, 0.89; sensitivity, 80%; specificity, 83%; accuracy, 82%).

Figure 5 shows two examples of CT values of lower leg muscles with lower (≤ 7) or higher (> 7) DSA runoff score. Comparing the two sets of heatmaps, the lower leg muscles of patient with mild PAD is shown as red regions with higher CT values. And the patient with severe PAD is shown as blue regions with lower CT values.

Examples of CT values of lower leg muscles with lower (≤ 7) or higher (> 7) digital subtraction angiography (DSA) runoff score. a A 67-year-old man with mild peripheral arterial disease (PAD) (DSA score = 1). The first image shows a coronal view of lower extremity vessels. The second image is a localizer. The third column are axial computed tomography angiography (CTA) images, corresponding to the white lines in the localizer image. The fourth column shows the heatmaps of CT values of the muscles in the third column images. The red area of the color bar reflects high CT values and the blue reflects low values. b A 63-year-old woman with severe PAD (DSA score = 16). Images of b represent the same meaning as a

To validate whether the severity of PAD can be predicted using only partial lower extremity images, we extracted muscle features from sub-datasets 1 to 5 and input them into the logistic regression equations of LRM-I and LRM-II to evaluate the predictive performance of each subset. The results (Table 5) showed that the AUC of LRM-I on sub-datasets 1 to 5 were 0.79, 0.79, 0.82, 0.83, and 0.78, respectively. The AUC of LRM-II on the subsets were 0.86, 0.86, 0.88, 0.90, and 0.86, respectively. The trend of the AUC values for the five subsets was consistent for both LRMs, and the largest AUC were obtained in sub-dataset 4 (LRM-I, 0.83; LRM-II, 0.90). Details of the other evaluation indices can be found in Table 5.