## 2 Materials and Methods
### Shap  
A common measure used to gauge the strength of each feature's effect on the outcome, which is the focus of the present study, is the so-called SHapley Additive exPlanation (SHAP) [19]. This is the equivalent in the machine learning literature to the Shapley values in cooperative game theory introduced by [PERSON] in 1951 [44]. SHAP values provide interpretability by constructing a simpler explainable model in the local neighborhood of each point in the dataset. Thus, given a learned Machine Learning model \(f\), a local approximation \(g\) can be formulated as follows [19]:  
\[g(\mathbf{u})=\phi_{0}+\sum\limits_{i=1}^{n}\phi_{i}u_{i}, \tag{2}\]  
where \(n\) is the number of features, \(\mathbf{u}\in\mathbb{R}^{n}\) is a binary vector whose value in the \(i^{\text{th}}\) position denotes whether or not the corresponding feature was used in the prediction, and \(\phi_{i}\) denotes the SHAP value of that feature, i.e., the strength of its contribution to the model's output.  
The values of the \(\phi_{i}\)s, following the notation of [PERSON] et al. [45], are computed as follows: let \(N=1,2,\ldots,n\) be the set of features used and let \(S\subseteq N\) be a subset of \(N\); then, we have [45, 46]  
\[\phi_{i}=\sum_{S\subseteq N\setminus\{i\}}\frac{|S|!(n-|S|-1)!}{n!}[f(S \cup\{i\})-f(S)]. \tag{3}\]  
Intuitively, this corresponds to the weighted average over all feature combinations (coalitions) of trip inhe difference in the model prediction with and without the inclusion of the \(i^{\text{th}}\) feature.  
As has been mentioned, the above \(\phi_{i}\) values pertain to a specific point. For example, considering the pair (Red, Black), there were 102 Red structures in the undersampled dataset and 90 Black ones, yielding \(90\times 102=9180\) pairs, i.e., samples in the transformed space, as shown in Table 1.  
Thus, we have a matrix \(\mathbf{\Phi}\in\mathbb{R}^{9180\times 13}\) in which each value \(\phi_{ij}\) is the SHAP value of the \(j^{\text{th}}\) feature calculated at the \(i^{\text{th}}\) sample. Thus, in order to obtain an aggregated value for the whole dataset, we used a normalized norm of each column in the matrix. We compared the results obtained using the \(\mathcal{L}_{1}\) norm (sum of absolute values), which is the most commonly used in the literature, and the well known Euclidean norm \(\mathcal{L}_{2}\), which increases the contribution of larger values while simultaneously reducing the effect of smaller noisy components. Thus, for each feature \(j=1,2,\ldots,13\) we considered the alternatives as obtained by Equation (4):  
\[\overline{\phi}_{ij}=\left\{\begin{array}{ll}\frac{1}{m}\sum_{i=1}^{m}|\phi _{ij}|,&\text{using }\mathcal{L}_{1}\\ \frac{1}{m}\sqrt{\sum_{i=1}^{m}\phi_{ij}^{2}},&\text{using }\mathcal{L}_{2} \end{array}\right. \tag{4}\]  
where \(m\) is the number of samples in the transformed space for each pair, as shown in Table 1. The computation of the SHAP values was carried out using the dedicated Python library by [PERSON] et al. [47].  
Thus, our overall proposed methodology comprises the following steps. For each damage threshold: (1) obtain the transformed inputs \(\mathbf{X}^{new}\) and outputs \(y^{new}\) as described in Section 2.2; (2) train the corresponding binary classification ML model using data for the particular damage threshold; and (3) obtain the feature importance metrics of the trained ML model using SHAP values via Equation (4).