#### 2.1.4 Data Processing  
The matrix X, containing the spectra of the standard samples, was pre-processed first with the Weighted Least Squares (WLS) algorithm based on the Whittaker filter (\(\lambda\) = 100, \(p\) = 0.001) to establish a baseline, followed by a Standard Normal Variate (SNV) normalization. The adaptive iteratively reweighted Penalized Least Squares (airPLS) algorithm [72] was also tested for the baseline treatment. The matrix containing the concentration of nineteen elements, including major and minor components, in the seventeen standard samples (Y17,19) was column-wise autoscaled. The instrument was calibrated through the sequential selection of individual elements (i.e., one column vector at a time). The Kennard-Stone algorithm [73] was used to split the X spectral matrix and the Y concentration matrix into calibration (X\({}_{\rm cal}\) and Y\({}_{\rm cal}\), 9 samples) and validation (X\({}_{\rm val}\) and Y\({}_{\rm val}\), 8 samples) sets. According to [PERSON] et al. [74], the validation step performed with a subset never used during the calibration and cross-validation steps can be regarded as an external test. However, to enable a more general assessment of the model's predictive capabilities, another validation step obtained using a strict external test set, namely, a distinct dataset differing not only in composition but also originating from a separate sampling campaign, would be desirable. The initial calibrated PLS regression model (X\({}_{\rm cal}\) vs. Y\({}_{\rm cal}\)) was input to a VIP analysis [75]. Since the average of squared VIP scores equals 1, a spectral variable with a VIP score equal or greater than 1 (one) was considered important in the model, whereas variables with VIP scores significantly lower than 1 (one) were excluded from the model. The performance of this step was assessed using the leave-one-out cross-validation approach. The output of the VIP analysis was a reduced spectral matrix (i.e., X\({}_{\rm cal}\)?) that was further independently fitted against each Y\({}_{\rm cal}\) column vector, corresponding to an element, to build the final PLS-VIP regression model. The optimal number of latent variables for the PLS-VIP analysis was selected by minimizing the Root Mean Squared Error (RMSE) in leave-one-out cross-validation. Then, the model was validated using the X\({}_{\rm val}\) and Y\({}_{\rm val}\) subsets. After the validation step, the model was applied to the spectral matrix of the archaeological samples (i.e., A) to obtain unknown element concentration, together with the uncertainty associated with the prediction. In all steps, namely, calibration (cal), cross-validation (CV), and validation (val), the quality of the results was assessed by evaluating the RMSE (see Equation (1)), bias, and \begin{table}
\begin{tabular}{c c c} \hline
**Sample Name** & **Material** & **Replicates** \\ \hline TY\_MA\_7 & Colorless glass & 1 \\ TY\_MA\_8 & Green glass & 1 \\ TY\_MA\_9 & Green glass & 3 \\ TY\_MA\_15 & Glass & 3 \\ TY\_MA\_20 & Glass & 3 \\ TY\_MA\_51 & Obsidian & 2 \\ TY\_MA\_52 & Volcanic rock & 7 \\ TY\_OU\_14 & Brick & 1 \\ TY\_OU\_15 & Brick & 1 \\ \hline \end{tabular}
\end{table} 
Table 1: Number of replicates collected for each investigated archaeological sample. correlation coefficient (R\({}^{2}\), see Equation (2)). The subscripts \"cal\", \"CV\", and \"val\" were added to RMSE and R\({}^{2}\) to define to which phase it refers (e.g., RMSE\({}_{\text{cal}}\)). The Limit of Detection (LOD) and Limit of Quantification (LOQ) were assessed in the calibration phase by Equations (3) and (4).