## 4 Discussion  
In this study, we proposed a method to compute representative SHAP values in a repeated nested CV procedure. Regarding the standard performance of machine learning models, in which average performance metrics are usually given, representative explainable values acting as a final assessment of the behavior of the entire model are also essential. Whereas the current literature mainly focuses on SHAP values computed on the entire dataset, we propose separate representative SHAP values for the training and test sets to allow a rigorous assessment of the generalization abilities of the SHAP explanations of a trained model. Based on traditional Shapley values [10], SHAP uses a game-theoretic framework to reframe the task of explaining the contribution of different features to the model output for a particular instance. However, since the SHAP values depend on the model predictions, variability in the performance of re-trained models leads to variability of the model's explainability through SHAP values. An example of this effect is shown in Figure 6, in which we report the frequency with which each feature was identified as the most important (highest impact as measured by the absolute value of the SHAP value over all the given samples) across all outer CV test folds in 100 repetitions for each regression/classification task. It is apparent that the most impactful feature changes over the different folds and repetitions. Still, the order of the most impactful features in a single iteration may differ between training and test sets--see, for example, the summary bar plot representing global feature importance for the regression task using the NKI dataset in Figure 7.  
Figure 6: Relative frequency of neuroimaging features. For each regression/classification task, the frequency with which each feature was identified as the most important (highest impact as measured by the absolute value of the SHAP value over all the given samples across all outer CV test folds in 100 repetitions) is shown. Figure 7: Summary bar plot representing global feature importance for a single iteration of the regression task using the NKI dataset, separately, for the training (on the **left**) and test sets (on the **right**). The order of the most impactful features in one iteration differs between training and test sets. Previously, SHAP values were estimated for the test set of a single nCV repetition, thus generating potentially unstable explanations [36; 66; 67]. [PERSON] et al. proposed using SHAP values in the test sets of a repeated CV without detailing the procedure adopted and considering only the average impact on model output magnitude, thus losing information about the signs of the SHAP values, which inform the positive/negative association with each feature [68]. In two recent works [11; 37], the average of SHAP values of samples in test sets among 100 repetitions of an ML model has been applied. In [11], the authors trained a deep neural network (DNN) model for age prediction using MRI data from the Autism Brain Imaging Data Exchange (ABIDE I) dataset collected from 17 international sites. For hyperparameter tuning, they used a leave-one-site-out CV, where the data from one site was adopted as a test set to evaluate the model's performance, while the data from all other sites was used as a training set. After each CV, they randomly undersampled the training set 100 times by removing a percentage of the samples in each iteration to produce small variations of the composition of the set and trained the DNN model to predict the subjects' age. They tested the DNN models on each test set sample, collecting 100 MAEs and SHAP values and averaging them for each sample. In [37], the authors adopted a leave-one-subject-out CV strategy using MRI data from the Alzheimer's Disease Neuroimaging Initiative (ADNI) dataset, splitting the dataset into as many sets as the number of subjects. One subject was randomly selected for testing, while the others were used to train the model.