Jhon Buendia

Julien Gornay

Marion Lacoue-Negre

Silvia Mas

Jihane Er-Rmyly

Ryad Bendoula

Jean-Michel Roger

A novel analysis methodology for preprocessing methods effectiveness determination in reducing undesired spectral variability in near-infrared spectra acquisition

Keywords: Spectral variability, Hierarchical Clustering Analysis (HCA), Principal Components Analysis (PCA), Preprocessing effectiveness, outliers, Qresidual, Hotelling's T², reproducibility, repeatability, sample temperature

This study uses a novel analysis methodology based on the Hierarchical Clustering Analysis (HCA) to determine the effectiveness of different preprocessing methods in minimizing undesired spectral variability in near-infrared spectroscopy due to both the consecutive and repetitive acquisition of the spectrum and the sample temperature. Nine preprocessing methods and different combinations of them were evaluated in four case studies: reproducibility, repeatability, sample temperature, and combination of the before mentioned cases. Eighty-four spectra acquired on seven different hydrocarbon samples from catalytic conversion processes have been selected as the real case study to illustrate the potential of the mentioned methodology. The approach proposed allows a more detailed discriminatory analysis compared to the classical methods for comparing the between-class and the within-class variances, such as the Wilks' lambda criterion, and hence constitutes a powerful tool to determine adequate spectral preprocessing strategies. This study also proves the potential of the discrimination analysis methodology as a general scheme to identify atypical behaviors either in the spectrum acquisition or in the measured samples.

Introduction

In the past few decades, the use of Near-Infrared Spectroscopy (NIRS) in the development of nondestructive and rapid measurement applications has been significantly increasing in several industries, such as food [START_REF] Fernández-Novales | A feasibility study on the use of a miniature fiber optic NIR spectrometer for the prediction of volumic mass and reducing sugars in white wine fermentations[END_REF][START_REF] Saeys | Multivariate calibration of spectroscopic sensors for postharvest quality evaluation: A review[END_REF] , pharmaceuticals [START_REF] Blanco | Analysis of pharmaceuticals by NIR spectroscopy without a reference method[END_REF][START_REF] De Beer | Near infrared and Raman spectroscopy for the in-process monitoring of pharmaceutical production processes[END_REF] , and petroleum [START_REF] Balabin | Neural network (ANN) approach to biodiesel analysis: Analysis of biodiesel density, kinematic viscosity, methanol and water contents using near infrared (NIR) spectroscopy[END_REF][START_REF] Zanier-Szydlowski | Control of Refining Processes on Mid-Distillates by Near Infrared Spectroscopy[END_REF] . Due to the recent growth boom in using the NIR spectrum for real-time acquisition data [START_REF] Wahl | Continuous monitoring of API content, API distribution and crushing strength after tableting via near-infrared chemical imaging[END_REF][START_REF] De Oliveira | Process modeling and control applied to real-time monitoring of distillation processes by near-infrared spectroscopy[END_REF] , the need to determine, analyze, and minimize spectral variability that is not associated with the physicochemical characteristics of the sample has notably arisen. This need becomes particularly evident when spectral variability is mainly generated by factors associated with the spectrum acquisition, such as spectrometer system, operator, measurement conditions, and environmental factors such as temperature and humidity, rather than the physicochemical characteristics of the sample. An example of this is the possible generation of spectral variability in the consecutive and repetitive acquisition of NIR spectra on a sample whose physicochemical characteristics remain constant over the spectrum acquisition. A lack of minimization of this type of spectral variability can result in inaccurate analysis and interpretation of spectroscopic information, misleading conclusions and flawed decision making [START_REF] Chauchard | Correction of the temperature effect on near infrared calibration-application to soluble solid content prediction[END_REF][START_REF] Igne | Robustness Considerations and Effects of Moisture Variations on near Infrared Method Performance for Solid Dosage Form Assay[END_REF] .

Among the classical performance parameters needed to validate a measurement methodology, precision is the most affected by the aforementioned factors. Precision is defined as the closeness of agreement between measured values obtained by replicate measurements on the same or similar samples under conditions of repeatability or reproducibility [START_REF] Betz | Accuracy, precision, and reliability of chemical measurements in natural products research[END_REF] . Repeatability conditions include the same measurement procedure, the same operator, the same instrument and measurement conditions, the same location, and a short interval between repetitions [START_REF] Olivieri | Validation and Error[END_REF] . On the other hand, reproducibility implies successive measurements of the same sample under changing measurement conditions [START_REF] Iso | Measurement management systems -Requirements for measurement processes and measuring equipment[END_REF] , such as measurement principle, measurement method, operator, measurement instrument, reference standard, location, conditions of use, and time. The NIR spectra acquisition is very sensitive to any change in measurement; even ensuring that both the spectrum acquisition conditions and the physicochemical characteristics of the sample do not change in a repetitive NIR spectrum acquisition, the resulting spectra may have differences that can lead to random errors and deviations, which must be corrected or minimized.

Due to its high impact on NIR spectra acquisition accuracy, the temperature is the most studied influencing parameter [START_REF] Wülfert | Influence of Temperature on Vibrational Spectra and Consequences for the Predictive Ability of Multivariate Models[END_REF][START_REF] Abe | A Study on the Universality of a Calibration with Sample Temperature Compensation[END_REF] . Hansen et al. [START_REF] Hansen | Tolerance of near Infrared Calibrations to Temperature Variations; A Practical Evaluation[END_REF] showed that molecular bonds vibration intensity depends on temperature, leading to changes in the spectrum according to temperature variation. Furthermore, some physicochemical properties of samples, such as viscosity and density, are temperature-dependent, and many changes in the sample due to temperature are not permanent and do not reflect the intrinsic nature of the sample [START_REF] Jasem | The effect of temperature and pressure on the viscosities of crude oils and their mixtures[END_REF][START_REF] Luo | Effects of asphaltene content on the heavy oil viscosity at different temperatures[END_REF][START_REF] Payri | The effect of temperature and pressure on thermodynamic properties of diesel and biodiesel fuels[END_REF] . Nevertheless, these changes can significantly affect spectrum acquisition. As with the spectral variability generated by repetitive spectrum acquisition of a specific sample, the variability caused by sample temperature must be minimized to ensure the reliable description of the sample physicochemical behavior from the spectroscopic information extracted. Data preprocessing is a common step for reducing undesired effects and for minimizing spectral variability. There are different preprocessing algorithms for the correction of the undesired spectral variation; these can be divided into two main categories: scatter-correction methods, employed to correct the additive and multiplicative effects, and spectral derivatives, used to minimize the sources of unwanted and non-informative spectral variations [START_REF] Rinnan | Review of the most common pre-processing techniques for near-infrared spectra[END_REF] . Among the most commonly preprocessing methods used in NIR spectroscopy, Savitsky-Golay derivative (Sav-Gol) [START_REF] Savitzky | Smoothing and Differentiation of Data by Simplified Least Squares Procedures[END_REF] , Extended Multiple Signal Correction (EMSC) [START_REF] Harald | Extended multiplicative signal correction and spectral interference subtraction: new preprocessing methods for near infrared spectroscopy[END_REF] , Standard Normal Variate (SNV) [START_REF] Barnes | Standard Normal Variate Transformation and De-Trending of Near-Infrared Diffuse Reflectance Spectra[END_REF] , and recently, the Variable Sorting for Normalization (VSN) [START_REF] Rabatel | VSN: Variable sorting for normalization[END_REF] , can be highlighted. However, the effectiveness of the preprocessing methods is highly dependent on the type of spectroscopic information analyzed and the factors that are causing its variability [START_REF] Chauchard | Correction of the temperature effect on near infrared calibration-application to soluble solid content prediction[END_REF] .

The preprocessing method effectiveness evaluation is generally based on the performance of prediction models [START_REF] Gholizadeh | Comparing different data preprocessing methods for monitoring soil heavy metals based on soil spectral features[END_REF][START_REF] Jiao | Preprocessing methods for near-infrared spectrum calibration[END_REF][START_REF] Brown | Derivative Preprocessing and Optimal Corrections for Baseline Drift in Multivariate Calibration[END_REF] . Among the contributions reported in the literature, the work of Gerretzen et al., [START_REF] Gerretzen | Simple and Effective Way for Data Preprocessing Selection Based on Design of Experiments[END_REF] which presents a novel approach for the selection of the most appropriate preprocessing methods based on the design of experiments, is worth mentioning. Similarly, the studies of Devos et al., [START_REF] Devos | Parallel genetic algorithm co-optimization of spectral pre-processing and wavelength selection for PLS regression[END_REF] and Allegrini et al., [START_REF] Allegrini | An integrated approach to the simultaneous selection of variables, mathematical pre-processing and calibration samples in partial least-squares multivariate calibration[END_REF] which, by means of a parallel workflow approach of preprocessing and variable selection, present an interesting alternative to the optimization of the preprocessing method selection. Nonetheless, the application of these approaches may be limited when the variability of the physicochemical characteristics of the samples is negligible, but significant spectral variability exists as a result of the repetitive spectrum acquisition and the sample temperature. In that case, a different analysis approach may yield more detailed results, helping to improve understanding of the impact of these parameters. Another less common approach to assessing preprocessing methods effectiveness is analyzing the spectral variance [START_REF] Luthria | A comparison of analytical and data preprocessing methods for spectral fingerprinting[END_REF] . Different statistical tools are available to determine both within-class variance (multiple measurements of the same sample) and between-class variance (measurements of different samples). One of the most common criteria used to evaluate betweenclass and within-class variances is the Wilks' lambda [START_REF] Wilks | The large-sample distribution of the likehoold ratio for testing composite hypotheses[END_REF] .

In this study, a novel and general strategy based on the Hierarchical Clustering Analysis (HCA) [START_REF] Rokach | Clustering Methods[END_REF] was proposed for evaluating the effectiveness of preprocessing methods in reducing the spectral variability generated by parameters related to the continuous and dynamic spectrum acquisition. To this aim, the effectiveness of nine preprocessing methods and different combinations of them in minimizing undesired spectral variability due to repeatability, reproducibility, sample temperature and combination of these parameters were evaluated. Eighty-four spectra acquired on seven different hydrocarbon samples from catalytic conversion processes have been selected as the real case study to illustrate the potential of the mentioned methodology.

To obtain reliable conclusions and validate the results obtained by the analysis methodology proposed, the Wilks' lambda criterion [START_REF] Wilks | The large-sample distribution of the likehoold ratio for testing composite hypotheses[END_REF] was used as a reference method.

Material and methods

Samples

Twenty-four vacuum gasoil (VGO) samples were processed in the catalytic conversion pilot plant reactors at IFPEN (Solaize, France). From these reactors, ninety-three different hydrocarbon samples, known as total effluent, were obtained (see references [START_REF] Bricker | Hydrocracking in Petroleum Processing[END_REF][START_REF] Speight | Hydrocracking[END_REF] for a detailed description of catalytic conversion processes). From these 93 samples, 7 samples were selected, ensuring their representativeness and physicochemical diversity. Table 1 summarizes four relevant physicochemical properties of the selected samples: the density [START_REF] Astm D | Standard Test Method for Refractive Index and Refractive Dispersion of Hydrocarbon Liquids[END_REF] and the simulated initial boiling point and distillation temperatures range to obtain both 5% and 95% of sample distillate (Simulated Distillation IBP, T5 and T95) [START_REF]Standard Test Method for Boiling Range Distribution of Petroleum Fractions by Gas Chromatography[END_REF] . It can be observed that physicochemical variability between the selected samples is guaranteed.

Spectral acquisition

The spectra were recorded with a Fourier Transform Near-Infrared spectrometer (FT-NIR) MATRIX-F (Bruker, Optik GmbH, Ettligen -Germany) within the range of 9090 -4600 cm -1 (1100 -2160 nm) and a resolution of 4 cm -1 . 32 scans were used to obtain the final spectrum after each measurement. For acquiring absorbance spectra, the spectrometer system was equipped with an immersion transflectance probe with an optical path fixed at 2 mm withstanding temperatures ranging from -40

°C to 200 °C. The software used with the spectrometer was OVP (OPUS Validation Program -Bruker, Optik GmbH, Ettligen -Germany) which automatically performs a series of analyses of the instrument's performance, evaluates them and ensures that it is operating within specifications.

Besides, to ensure the spectrometer operation within specifications and that the spectral variability generated was due to the parameters evaluated and not to the instrument's inadequate functioning, the spectrometer performance was validated once a day using cyclohexane as an external reference sample. Before NIR analysis, the samples were heated in closed flasks at 60°C for 1 hour in a water bath and shaken manually to ensure their liquid state and homogeneity. The initial boiling point (IBP) reported in Table 1 guarantees no loss of volatiles.

Ensuring the integrity and stability of both the sample and the NIR spectrum acquisition conditions, spectral variability due to repeatability, reproducibility, and sample temperature was generated. A short description of the spectrum acquisition for the cases evaluated in this study is presented below.

• Case 1. Spectral variability due to reproducibility: Each of the seven samples was analyzed one time per day during five consecutive days at 60 °C. Thirty-five spectra were obtained.

• Case 2. Spectral variability due to repeatability: Each of the seven samples was analyzed three times on the same day at 60 °C. All samples were analyzed in less than 8 hours.

Twenty-one spectra were obtained.

• Case 3. Spectral variability due to the sample temperature: Each of the seven samples was analyzed at five different temperatures, ranging from 60 °C to 80 °C with a temperature increment of 5°C. The samples were heated in closed flasks at the desired temperature for 1 hour. Evaporation losses of volatiles were null or negligible (see IBP in Table 1). Twenty-eight spectra were obtained.

• Case 4. Spectral variability due to the combination of the aforementioned cases: In this case, all spectra acquired in the above-described cases were used.

For each case, a matrix was generated. Each analyzed sample was defined as a class; thus, seven classes were defined in all matrices.

Analysis methodology

The main steps of the data analysis workflow proposed in this study are schematized in Figure 1. A brief explanation of the procedure is given as follows.

The first step consisted in preprocessing each of the generated matrices. The nine most commonly preprocessing methods used in NIR data were divided into two categories; the filtering and normalization methods. The preprocessing methods from each category were analyzed individually.

If the total reduction or compensation of the studied spectral variability were not achieved, the evaluated preprocessing method were complemented with the methods belonging to the opposite category. This allowed the evaluation of all possible combinations and order of use of the preprocessing methods. The methods evaluated are described in Table 2. It should be emphasized that each preprocessing scenario evaluated includes the data centering by columns. Afterward, a preliminary inspection and dimension reduction of the corresponding dataset was performed by principal component analysis (PCA) [START_REF] Wold | Principal component analysis[END_REF] . The number of the chosen principal components (PCs) captured at least 99 % of the total variance in the dataset. The Q residual and Hotelling's T² tests were performed to determine the possible presence of anomalous data [START_REF] Chen | A regularized Hotelling's T2 test for pathway analysis in proteomic studies[END_REF] .

The chosen PCs scores were then used to perform the hierarchical clustering analysis (HCA) employing Ward's algorithm and the Mahalanobis distance. The HCA aims to group clusters to form a new one to either minimize a statistical distance between classes or maximize a measure of similarity between them [START_REF] Sasirekha | Agglomerative Hierarchical Clustering Algorithm-A Review[END_REF][START_REF] Murtagh | Ward's Hierarchical Agglomerative Clustering Method: Which Algorithms Implement Ward's Criterion?[END_REF] . The analysis starts with as many groups as individuals contained in the dataset.

From these initial groups, clusters are formed in an ascending manner until all cases treated are included in at least one of them. Ward's algorithm seeks to minimize each group variance by calculating all samples mean in each cluster. The algorithm then calculates each case distance and the cluster mean, adding up the distances between all cases. Finally, the clusters whose sum of distances is minimal are grouped. This procedure creates homogeneous groups of a similar amount of individuals. For achieving the grouping of classes, it is necessary to define a comparison parameter to calculate the variance of each class concerning the others. The most common is the Mahalanobis distance [START_REF] Martens | Multivariate calibration[END_REF] .

A widespread manner of displaying the cluster analysis results is constructing a tree diagram known as a dendrogram. The resulting diagram shows the different groups' clustering order and the association measure's value, also known as the fusion level. The fusion level was defined for obtaining seven clusters corresponding to the seven classes. Finally, the number of correctly grouped sample measurements in each cluster was determined, and the percentage of clustering was calculated as: Number of correctly grouped samples / Total number of samples * 100.

These steps were repeated for each preprocessing method scenario, and the results obtained were compared using the percentage of samples correctly grouped as a figure of merit to determine the effectiveness of the preprocessing methods evaluated.

All the analyses were conducted with the PLS_Toolbox version 8.8 (Eigenvector Research Inc., Wenatchee, WA, USA) for MATLAB version R2019b (MathWorks, Natick, MA, USA).

Results comparison

To validate the results obtained with the methodology proposed in this study, the Wilks' lambda was used as a comparative criterion. This criterion evaluates how well the data set classes are separated by calculating a ratio involving between-class and within-class variances. Several versions of the Wilks' lambda exist. In this article, the ratio of the between-class variance over the total variance was used. This ratio varies between 0 and 1, where 0 means that all the classes are superimposed, and 1 means that all the classes are perfectly separated.

Results and discussion

This section shows the results obtained from the application of the proposed methodology. A comprehensive description of its application for the case study of reproducibility (Case 1) is presented. However, only the main results from the case studies of repeatability, temperature effect, and the combination of all cases are showed. Finally, the validation and the advantages of the proposed methodology are offered.

Case 1 -Reproducibility

For the reproducibility case, 35 spectra were used (5 spectra for each sample, see section 2.2). Figure 2 shows the raw spectra over the entire spectral range used. At first glance it can be observed the variability of spectra due to the physicochemical nature of the sample, showing, with some exceptions, a trend consistent with the properties reported in Table 1; i.e., the spectra acquired on the sample with the lowest density (sample 1 -dark blue color) is at the bottom of the plot, and the spectra acquired on the sample with the highest density (sample 7 -red color) is at the top. This observation becomes more evident in the black square of Figure 2 where the spectra are magnified over a defined range (1340nm -1430nm). In this same figure it is possible to visualize that the reproducibility measurements made on a single sample generate a variability that has a similar behavior to the variability caused by multiplicative effects, which can impact the final results obtained. In this first case, the proposed methodology application shows the impact of the variability in NIR spectra acquisition caused by reproducibility measurements as well as the effectiveness of preprocessing methods in minimizing this spectral variability.

Firstly, an exploratory analysis using PCA was performed to visualize the similarity/dissimilarity among the 5 spectra (1 per day) acquired for each of the 7 hydrocarbon samples analyzed. The data set (35 spectra) was merely centered. The first two components (PC1 and PC2) explained 96.7% of the data set total variance. From the score plot of the first two components (see Figure 3a), it can be seen that all 5 measurements of samples 1 to 6 are grouped in a consistent pattern regarding to the variance within classes. However, there are measurements of different samples that intersect with each other (between-class variance) (see Figure 3a-b), and hence their clustering can be influenced.

Additionally, from these score plots, it can also be observed that the measurement of sample 7 on the fifth day is relatively distant from the other measurements of this sample. This could mean the presence of possible outliers in the dataset. Figure 3c shows the Hotelling's T² and Qresidual scores for the dataset. It can be observed that the same measurement identified previously (Sample 7, 5 th acquisition day) was found to be above the threshold of both tests. Therefore, this measurement can be confirmed as an outlier. The scores of the first 4 principal components yielded by the PCA analysis (99% explained variance)

were used to perform the HCA analysis. Figure 4 shows the dendrogram obtained from the HCA, where a fusion level (black line) for obtaining seven clusters corresponding to the seven hydrocarbon samples is defined. This Figure also shows the correct grouping percentage achieved for each class.

From the table embedded in the Figure, it can be observed that no sample achieves an accurate grouping of all its measurements. Samples 2 and 7 (orange and red classes) are the classes with the highest correct grouping percentage (4 out of 5 for 80%), while classes 3 and 4 (yellow and purple classes) do not have any correctly grouped measurements. Furthermore, measurement 35 that belongs to class 7 (red class) and identified in the previous steps as a potential outlier is the only measurement grouped in class 4 (purple class). Based on these results, it can be considered that this measurement has no similarity with any of the other 34 measurements, which confirms its outlier status.

Table: correct grouping percentage for each class

In order to compensate the spectral variability caused by the lack of reproducibility and hence achieve a better clustering of all the classes, the use of appropriate preprocessing methods is needed. Nine preprocessing methods and different combinations of them were evaluated. Figure

App 1, which can be found in the article's supplementary section, summarizes the correct grouping percentage results for each preprocessing scenario evaluated in the 4 cases. From this Figure, it can be seen that for Case 1 (diamond shape -blue color), none of the evaluated scenarios reaches a correct grouping percentage of 100%, implying that none of the scenarios achieved the total reduction of the spectral variability generated by the reproducibility measurements. With a correct grouping of 85%, the EMSC was the most effective preprocessing method scenario. Figure 5 shows for this scenario the score plot of the first two principal components of the PCA analysis and the

Qresidual & Hotelling's T² results achieved using 2 PCs. It can be seen that all measurements are consistently grouped and are within the threshold of the two tests, except for measurement 35

(potential outlier identified). From the HCA dendrogram and the table of its corresponding correct grouping percentage (see Figure 6), it could be observed that this measurement is still wrongly grouped as the unique measurement in class 2 (orange class). From these results, it can be assumed that all other measurements could be correctly grouped without this misgrouped measurement.

Therefore, measurement 35 was removed from the data set, and the preprocessing methods scenarios were re-evaluated to confirm this assumption.

Figure App 1 (diamond shape -red color) shows that by removing measurement 35 from the data set, correct grouping of all measurements is possible in 6 scenarios (EMSC, AWLS-B+MSC, AWLS-B+SNV, AWLS-B+VSN, SG-D+SNV, and MSC+AWLS-B). However, only one scenario uses a single method, which is the EMSC. In order to prevent loss of relevant information, the use of a minimum number of methods in the data preprocessing is generally recommended. Therefore, in this case, EMSC could be selected as the most efficient preprocessing method scenario to reduce the spectral

Case 2 -Repeatability

Figure 7 shows the raw spectra used in the analysis of the variability caused by repeatability measurements. Analogous to case 1, the spectral variability generated by the physicochemical nature of the sample can be observed, presenting the same relationship (trend) with the properties reported in Table 1. However, the black square in Figure 7 shows that the variability caused by repeatability measurements seems to present a similar behavior to the variability generated by the combination of two different effects (additive and multiplicative), making the spectral differences of a single sample more evident in comparison with case 1.

In this second case, a total of 21 spectra (7 samples performed in triplicate) were analyzed to demonstrate the proposed methodology's ability to evaluate the effectiveness of preprocessing methods in minimizing unwanted spectral variability due to the lack of repeatability. The correct grouping percentage achieved from all the preprocessing methods scenarios is shown in Figure App 1 (squared shape -orange color). In the scenario that is assumed that the variability generated in the repeatability measurements do not have a significant impact on the grouping of sample measurements, that is, that no additional preprocessing methods are needed besides the data centering, only 48% of the measurements were correctly grouped (see Figure App 1 (squared shape -orange color) mean center scenario).

Therefore, it can be preliminarily concluded that preprocessing methods to reduce spectral variations due to repeatability are needed. From the PCA analysis (data not shown) of this scenario (centered data), it could be determined that the first repetition of each sample differs significantly from its second and third repetitions. However, all of them were within the Qresidual and Hotelling's T² test thresholds (data not shown). Therefore, they cannot be considered as anomalous data (outliers).

From all the evaluated scenarios, fourteen achieved a correct grouping of 100% (AWLS-B+MSC, AWLS-B+SNV, Dtd+MSC, Dtd+PQN, Dtd+SNV, EMSC+PQN, EMSC+SNV, MSC+AWLS-B, MSC+Dtd, PQN+EMSC, SNV+AWLS-B, SNV+Dtd, SNV+EMSC, and VSN). Among these scenarios, only one uses a single method, which is VSN (see HCA dendrogram in Figure App 4). Therefore, VSN could be selected as the most effective preprocessing method to minimize the spectral variability due to the lack of repeatability. As in Case 1, Figure App 1 shows that the order of use of the preprocessing methods affects the correct grouping result. Comparing the same scenarios as in Case 1, it can be reaffirmed that the results are more promising when the filtering methods are applied after the normalization methods.

From these results, it can be concluded that spectral variability due to repeatability has a lesser impact than those generated by reproducibility. Nevertheless, the proposed methodology demonstrated that the two cases' variability could be entirely compensated using an appropriate data preprocessing strategy.

Case 3 -Temperature effect

As previously mentioned in the introduction of the manuscript, sample temperature is one of the factors having significant impact on the NIR spectra acquisition. Figure 8 shows that the spectral variability generated by the sample temperature presents a behavior similar to the multiplicative effect. However due to the absorbance shift caused by the temperature increase, which prevents having a direct relationship between this parameter and the height of the acquired spectra, the spectral difference presents a nonlinear growth [START_REF] Roger | EPO-PLS external parameter orthogonalisation of PLS application to temperature-independent measurement of sugar content of intact fruits[END_REF] . This can be corroborated in the black square of Figure 8, where it is observed that with a temperature variation greater than 15°C the spectral variability is more evident than when the delta in temperature is less than 15°C. This non-linear impact of the sample temperature on the spectrum acquisition could limit the performance of the different preprocessing methods evaluated. To sum up the results obtained by using the analysis methodology proposed, it could be concluded that for cases 1 and 2 (reproducibility and repeatability), the variability affecting the measurement clustering was fully compensated using a single preprocessing method. On the other hand, cases 3 and 4 (sample temperature and combination of cases) needed the combination of two methods, and still, the correct grouping of all the measurements was not achieved. The sample temperature has a high impact on the spectrum acquisition. Therefore, it is recommended to use a methodology that evaluates this variable more thoroughly for a more efficient solution [START_REF] Chauchard | Correction of the temperature effect on near infrared calibration-application to soluble solid content prediction[END_REF][START_REF] Roger | EPO-PLS external parameter orthogonalisation of PLS application to temperature-independent measurement of sugar content of intact fruits[END_REF] .

Results comparison

In order to validate the consistency and reliability of the proposed approach, the analysis methodology results were compared with those obtained by Wilks' lambda criterion.

Table 3 summarizes each case's most relevant results achieved by both the proposed approach and Wilks' Lambda. The results shown in Table 3 validate the approach proposed in this study. It can be seen that the results between the two methodologies are comparable, except for Case 1, including all measurements, when the data have been only mean-centered. In this case, Wilks' Lambda value is close to 1, while the value obtained by the proposed methodology is 0.46. The difference may be attributable to the presence of the outlier identified in Case 1 and how each approach handles this type of data. While the Wilks' lambda criterion assumes that there is no presence of outliers in the analyzed dataset, the methodology used in this study provides a preliminary analysis of the dataset for the identification and removal of possible anomalous data. This premise can be supported by observing that the two approaches' results are comparable when the identified outlier is removed from the dataset (see Table 3 -Case 1 (Measurement removed)).

Moreover, the proposed methodology provides a more detailed discrimination analysis in comparison with Wilks' lambda criterion. As a way of example, both the proposed approach and the Wilks' lambda results obtained from the individual evaluation of the nine preprocessing methods in Case 1 were compared. From Table 4, it can be seen that Wilks' lambda criterion presents no significant differences in 4 preprocessing methods (SG-D, EMSC, MSC, SNV), which could lead to the conclusion that the 4 methods have equal effectiveness in minimizing the unwanted spectral variability. On the contrary, the proposed methodology allows to identify that out of these 4 methods, two are equally effective in minimizing the unwanted spectral variability (MSC, SNV), the least effective is the SG-D, and the most effective preprocessing method, which achieves the maximum compensation of the unwanted spectral variability, is the EMSC. Therefore, the methodology used in this work could help to select the preprocessing method in a more precise and reliable way. Finally, the proposed methodology allows identifying measurements and samples that have been properly or poorly discriminated, information that the Wilks' lambda does not provide as it is a global measurement.

Conclusions

The results obtained in this study shows the capacity of the analysis methodology used to assess the effectiveness of preprocessing methods in reducing the undesired spectral variability of near-infrared spectroscopy (NIRS) measurements in a more thoughtfully and detailed manner than other approaches based on the dataset variance analysis such as the Wilks' lambda criterion.

In this study, an original strategy not previously reported in literature was proposed to evaluate and determine the effectiveness of different preprocessing methods in minimizing the unwanted spectral variability due to parameters related to the continuous and repetitive NIR spectra acquisition such as repeatability, reproducibility, sample temperature, and the combination of these three parameters.

It is essential to stress the twofold benefit of using the proposed methodology. On the one hand, the detailed discrimination analysis provides a significant aid in determining the most effective data preprocessing scheme. On the other hand, the methodology provides a preliminary data analysis

Figure 1 -

 1 Figure 1 -Methodology flow diagram

Figure 2 -

 2 Figure 2 -35 raw spectra used in Case 1 evaluation over the entire spectral range. Legend: Color → samples analyzed. Black square → 35 raw spectra used in Case 1 evaluation magnified over the 1340nm -1430nm spectral range

Figure 3 -

 3 Figure 3 -a) Score plot of PC1 & PC2 for Case 1 with centered spectra, b) Score plot of PC3 & PC4 for Case 1 with centered spectra, c) Reduced Qresidual & Hotelling's T² for Case 1 using 4 PCs with centered spectra. Legend: Shapes → acquisition day. Color → samples analyzed (classes)

Figure 4 -

 4 Figure 4 -Dendrogram for Case 1 (all measurements) with centered spectra. Legend: Color samples measured (classes).

  variability due to the lack of reproducibility (see HCA dendrogram in Figure App 2 in the supplementary section).

Figure App 1 Figure 5 -

 15 FigureApp 1 also shows that both the order of use and the combination of preprocessing methods can influence the final result. An example of this statement can be seen by comparing the scenarios from Case 1 without the measurement 35, where the NW-D method was combined with the MSC, PQN, SNV, and VSN methods. When the NW-D method (filtering category) is used before applying the other preprocessing methods (normalization category), the correct clustering is lower (about 30%) compared to when using the normalization preprocessing methods before the NW-D method. It can also be found by comparing the same scenarios that using complementary preprocessing methods does not always yield better correct clustering results than using a single method (NW-D+PQN → 38% Vs NW-D → 47%). This is an important consideration when using more than one preprocessing method.

Figure 6 -

 6 Figure 6 -Dendrogram for Case 1 (all measures) with EMSC preprocessing. Legend: Color samples measured (classes).

Figure 7 -

 7 Figure 7 -21 raw spectra used in Case 2 evaluation over the entire spectral range. Legend: Color → samples analyzed. Black square → 21 raw spectra used in Case 2 evaluation magnified over the 1340nm -1430nm spectral range

Figure 8 -Figure 9 -

 89 Figure 8 -a) 35 raw spectra used in Case 3 evaluation over the entire spectral range. Legend: Color → samples analyzed. Black square → 35 raw spectra used in Case 3 evaluation magnified over the 1340nm -1430nm spectral range

Figure App 2 Figure App 3

 23 Figure App 2 Dendrogram for Case 1 (all measurements) with EMSC preprocessing. Legend: Color samples measured (classes).Table: correct grouping percentage for each class

Table 1 Samples physicochemical properties. SimDis IBP, T5 & T95 description: Simulated distillation to determine the temperatures to start the sample evaporation and to recover both 5% and 95% of sample distillate Sample Physicochemical Properties Distillation Temperatures Sample ID Density (gr/ml) IBP (°C) SimDis T5 (°C) SimDis T95 (°C)

 1

	Sample 1	0.8049	84.7	121.6	425.8
	Sample 2	0.8186	79.8	111.5	474.5
	Sample 3	0.8219	80.3	117.8	516
	Sample 4	0.8411	83.8	136.1	503.5
	Sample 5	0.847	83.7	129.1	512.6
	Sample 6	0.8962	148.9	227.6	504.6
	Sample 7	0.9181	159.8	231.6	502.9

Table 2 Preprocessing methods description

 2

	#	Category	Name	Acronym	Parameters
	1		Automatic Weighted Least Squares Baseline 20	AWLS-B
	2		Norris-Williams Derivation 38	NW-D	15-point window, gap size = 7, First order derivation
	3	Filtering	Savitsky-Golay Derivative 21	SG-D	15-point window, polynomial order = 2 First order derivative
	4		Detrend 23	Dtd	Polynomial order = 1
					Reference spectrum (basis to remove
	5		Extended Multiplicative Scatter/Signal Correction 22	EMSC	the scatter) = mean of each matrix generated, polynomial order = 2,
					whole spectral range
	6		Multiplicative Signal Correction 39 MSC	Reference data = mean of data, whole spectral range
	7		Standard Normal Variate 23	SNV
	8	Normalization	Probabilistic Quotient Normalization 40	PQN
	9		Variable Sorting for Normalization 24	VSN	Automatic calculation

Table : correct grouping percentage for each class

 :

Table 3 % Results comparison using Wilk's Lambda algorithm

 3

	Case	Preprocessing	Methodology	Wilk's
		method	(% Grouped/100)	Lambda
	Case 1	Mean Center	0.46	0.95
	(All Measurements)	EMSC	0.85	0.96
	Case 1	Mean Center	0.53	0.68
	(Measurement 35 removed)	EMSC	1.00	0.99
	Case 2	Mean Center VSN	0.48 1.00	0.69 0.99
	Case 3	Mean Center SavGol + SNV	0.49 0.82	0.55 0.95
	Case 4	Mean Center VSN+EMSC	0.51 0.96	0.66 0.99

Table 4 Case 1 detailed results comparison

 4

	Preprocessing Method	Explored Methodology (% Grouped/100)	Wilk's Lambda
	Mean Center	0.53	0.68
	AWLS-B	0.65	0.88
	NW-D	0.56	0.73
	SG-D	0.65	0.98
	Dtd	0.65	0.90
	EMSC	1.00	0.99
	MSC	0.82	0.99
	PQN	0.68	0.87
	SNV	0.82	0.99
	VSN	0.94	0.85

Table: correct grouping percentage for each class

Acknowledgments

The authors would like to thank IFP Energie Nouvelles for providing the hydrocarbon samples obtained in their HCK pilot plant reactors and the facilities for spectra acquisition and data analysis.

Thanks also go to Axel One for providing the spectrometer for the NIR spectra acquisition.

grouping in all their measurements (5 out of 5 = 100%). The class affecting the overall measurement clustering is sample 4 (purple), which does not have any measurements grouped correctly. As in Case 1, it could be assumed that the whole misgrouping of sample 4 is due to the presence of some atypical data. However, no measurement was found above the Qresidual and Hotelling's T² tests thresholds in any scenario (data not shown).

The results analyzed in this case show that sample temperature is a very influential parameter on the spectrum acquisition. Therefore, it is recommended to use a strategy that evaluates this variable's impact more thoroughly for a more efficient solution [START_REF] Chauchard | Correction of the temperature effect on near infrared calibration-application to soluble solid content prediction[END_REF][START_REF] Roger | EPO-PLS external parameter orthogonalisation of PLS application to temperature-independent measurement of sugar content of intact fruits[END_REF] .

Case 4 -cases 1, 2 & 3 combined

The parameters causing unwanted spectral variability evaluated in the 3 cases previously described are likely to occur simultaneously, mainly when online NIR measurement is used for real-time data analysis. For this reason, a fourth case was evaluated where the spectral variability generated by these three cases was combined.

The dataset used in this case comprises 84 spectra (35 from reproducibility Case, 21 from repeatability Case, and 28 from sample temperature Case). As in the 3 cases already studied, the Qresidual and Hotelling's T² analyses were applied to the centered dataset to determine the possible presence of atypical data (data not shown). The measurement identified as an outlier in Case 1 (Sample 7, 5 th acquisition day) once again exceeds the thresholds of the two tests. Thus, this measurement was removed, and the methodology proposed in this study was applied to the other 83 measurements.

In this case, 5 scenarios (EMSC, EMSC+PQN, PQN+EMSC, SNV+EMSC, and VSN+EMSC) had the best performance, where 80 out of 83 measurements were correctly grouped, that is, 96% of correct grouping. All 5 scenarios identified involve the use of EMSC, which, when evaluated individually, yields the same percentage of correct grouping (96%). It could be preliminarily inferred that for this case, the normalization methods do not give any additional improvement, and thus EMSC would be selected as the most effective preprocessing method scheme. Nevertheless, the use of the proposed methodology allows a more detailed analysis to determine at what level each type of variability is minimized in each of the 5 scenarios mentioned. In this way, the most effective preprocessing scheme's determination can be done more conveniently according to the researcher's objective.

Figure 9 shows the PCA score plot and the Qresidual & Hotelling's T² results for the VSN+EMSC preprocessing method scenario. The three measurements that were not correctly grouped correspond to the measurement at 80°C of samples 3 (highest value of residual Q test) and 4, and the step for identifying and removing the potential anomalous data from the dataset, thus improving the reliability of the final results.

The results obtained using the proposed analysis methodology suggest that the variability caused by repeatability and reproducibility can be fully corrected when using the adequate preprocessing scheme; however, no preprocessing scenario could entirely compensate the unwanted spectral variability caused by the sample temperature. Similarly, the detailed discriminant analysis employed in this study showed that the EMSC preprocessing method presents interesting and promising results in all cases.

The preprocessing scheme's ultimate selection should be conducted in a careful manner considering the researcher's objective. The proposed methodology offers an analysis strategy that could help determine the most effective preprocessing scheme more reliably.

The conclusions reached in this work promote further optimization and automation of the proposed methodology to improve its implementation in large datasets.

The strategy proposed was shown to work for a case study including seven different hydrocarbon samples but can be generally applicable in any study involving spectroscopic information analysis 6 Declaration of conflicting interests The Author(s) declare(s) that there is no conflict of interest

Appendix