J.
Vrábel
*ab,
E.
Képeš
ab,
P.
Nedělník
a,
J.
Buday
ab,
J.
Cempírek
c,
P.
Pořízka
*ab and
J.
Kaiser
ab
aCEITEC, Brno University of Technology, Purkyňova 123, Brno, 612 00, Czech Republic. E-mail: jakub.vrabel@ceitec.vutbr.cz; pavel.porizka@ceitec.vutbr.cz
bInstitute of Physical Engineering, Brno University of Technology, Technická 2, Brno, 616 00, Czech Republic
cDepartment of Geological Sciences, Faculty of Science, Masaryk University, Kotlářská 2, Brno, 602 00, Czech Republic
First published on 8th February 2023
The mutual incompatibility of distinct spectroscopic systems is among the most limiting factors in laser-induced breakdown spectroscopy (LIBS). The cost related to setting up a new LIBS system is increased, as its extensive calibration is required. Solving this problem would enable inter-laboratory reference measurements and shared spectral libraries, which are fundamental for other spectroscopic techniques. We study a simplified version of this challenge where LIBS systems differ only in the spectrometers used and collection optics but share all other parts of the apparatus and collect spectra simultaneously from the same plasma plume. Extensive datasets measured as hyperspectral images of a heterogeneous rock sample are used to train machine learning models that can transfer spectra between systems. The transfer is realized using a composed model that consists of a variational autoencoder (VAE) and a multilayer perceptron (MLP). The VAE is used to create a latent representation of spectra from a primary system. Subsequently, spectra from a secondary system are mapped to corresponding locations in the latent space by the MLP. The transfer is evaluated using several figures of merit (Euclidean and cosine distances, both spatially resolved; k-means clustering of transferred spectra). We demonstrate the viability of the method and compare it to several baseline approaches of varying complexities.
As a trade-off for its instrumental robustness, LIBS exhibits a high sensitivity to changes in the analyzed target's topography9 and matrix.10 More importantly, LIBS is significantly affected by the changes in experimental conditions (such as ablation11,12 and collection geometry13). A complete list of physics-related factors responsible for the structure of spectra is beyond the scope of this work and was described elsewhere.14 The most prominent and easily recognizable changes in the spectrum's structure can be attributed to the spectrograph and camera. Namely, spectrographs commonly differ in their spectral range and resolving power. Similarly, pixel size and the related quantum efficiency of a detector affect the resolution and overall structure of detected spectra. While the detector's response can be addressed by intensity calibrations using standard light sources, corrections for the spectrograph's impact remain elusive.
Problems with non-trivial signal dependence on experimental conditions and instrumentation could be partially neglected for a simple qualitative analysis. Tasks such as element detection and sample classification require only a wavelength position-calibrated spectrograph and a sufficient resolution. Therefore, even cost-efficient spectrographs could provide adequate performance for qualitative applications.15 In contrast, a reliable quantitative analysis is significantly more demanding, with the necessity of spectral intensity correction, compensation for the matrix effects, and instrumental limitations. A common practice is to create a calibration curve from a set of known standards, which has limited extrapolation reliability. An intriguing alternative to calibration curves is calibration-free (CF) LIBS,16 where plasma emission models are used. Despite many efforts, the applicability of CF-LIBS to real-world tasks and general experimental conditions is disputable.
The challenges described up to this point concerned measurements on a single instrument. However, many applications would benefit from the cooperation of several systems and the combination of their outputs, i.e. the transfer of data libraries. A representative example is the ChemCam device mounted on the Curiosity Mars rover that has an exact copy operating on Earth.17 The compatibility of measurements from two distinct systems is the missing element for a trustworthy inter-laboratory comparison and creation of shared spectral libraries in the LIBS community. Both technical and physics-related limitations rule out generally valid transfer (as spectral ranges and resolving power of detectors may vary significantly). However, even a reasonable approximation would have great potential for the practice.
A method we propose forspectral transfer is based on a client-server relationship between a preferred primary system and one or possibly many secondary systems. We use a variational autoencoder (VAE) to find a latent representation of the primary system spectra, which is expected to abstract from the specifics of the primary system. Then a multilayer perceptron (MLP) is trained to map the corresponding spectral pairs from the secondary system to the latent representation of the primary system. The final model for the spectral transfer between the two systems is obtained by joining the MLP with the decoder part of the autoencoder. The transfer methodology is illustrated in Fig. 1. It should be noted that the latent space is identical for all participating secondary systems. As a result, the amount of data that needs to be exchanged to include a new secondary system or to enrich a possible shared library of spectra is greatly reduced (as opposed to an approach where the full spectral space of the primary system is shared and not only the latent space), leading to multiple practical advantages. Additionally, a shared, low-dimensional representation can be directly used to create shared calibration, classification, and clustering models or to enhance any efforts directed toward the interpretability of the data and utilized models. In future work, we plan to expand on this point further by introducing physics-relevant constraints on the latent space. Finally, the latent space formed by the VAE can be used to generate new spectra.
We demonstrate that our methodology can efficiently transfer spectra (secondary → primary) even when the secondary system covers a considerably narrower spectral region. The performance is evaluated on out-of-sample data that were measured on a different day and location on the sample. Using multiple figures of merit, we compare the model to several baselines and provide a detailed analysis of the results.
Note that to establish transfer between the primary and the secondary systems, a one-to-one correspondence is required for training spectra. A common solution for this requirement is to use a set of shared standardized samples. A number of shared samples and their spectral variability determine the efficiency of the transfer. In practice, only a small number of shared samples are being used (max. tens) due to cost and unavailability of suitable standards. Even more, considering system fluctuations between the measurements of the standards and experimental error, the one-to-one correspondence is only approximate. We propose a novel approach for obtaining matched pairs of spectra, where we measure on both systems simultaneously from the same plasma plume. While this setup imposes certain limitations (e.g. a necessity of two systems being physically present together during the training phase, tedious synchronization, etc.), it creates the possibility to obtain an unprecedented amount of unique spectral pairs. This, in combination with the large-scale mapping of heterogeneous samples, allows us to precisely study the performance and limitations of the transfer between the given two systems.
Considering mainly the CT approaches that utilize some form of the spectral transfer as a part of the process, the main contrasts between our work and prior efforts are the possibility of a higher discrepancy between studied systems and the use of non-linear techniques (up to the exceptions mentioned below, linear models were predominantly used in previous studies).
Similarly to our work, the authors of (ref. 22) transferred NIR spectra from a system with a narrow wavelength range to a broader one using the direct standardization (DS) method (with various regularizers). Major drawbacks of the DS approach are the impossibility to express a complex non-linear dependency between systems and high computational space complexity (scaling N × M, where N is the number of shared data samples and M is the dimension of spectra). The latter rules out its use for large shared datasets.
A non-linear approach (pseudo-autoencoder) was used in ref. 23 to transfer NIR spectra. However, the transfer was performed only between systems with matching wavelength ranges and similar spectral responses. Additionally, the term autoencoder is used incorrectly in the mentioned work. The model's architecture only resembles that of a standard autoencoder (by the presence of a bottleneck), but it was not trained or used as an autoencoder (i.e. for dimension reduction and reconstruction). The utilized model is an extreme learning machine (architecture and forward pass are analogical to the MLP but weights are not trainable), which was shown to be significantly outperformed by MLPs (trained by the stochastic gradient descent) for large datasets.24
In LIBS, CT was studied in ref. 19, as an example of a general manifold alignment problem. The authors proposed a novel method called low rank alignment (LRA) that creates a shared low-dimensional embedding space, where the distance between the corresponding spectral pair embeddings is minimized. LRA was used for CT on ChemCam spectra measured at two different energies, but the dimensions were the same for both setups. The calibration model was built in the embedding space on high-energy spectra and used for predictions on low-energy spectra. The reconstruction from the embedding space to the spectral space of each system was not provided. The performance of the transferred calibration model is only compared to another manifold alignment technique (affine matching), but the absolute value of RMSE is missing. Also, as the authors state, the method has cubic time complexity w.r.t. a number of data examples, which makes it unsuitable for large shared datasets.
Last to mention is our previous work,25 where we used an MLP to transfer spectra from the ChemCam to the SuperCam to improve calibration models for studied oxides. We were able to improve the RMSE over models that were trained solely on the SuperCam system. The focus was placed on calibration models, and the spectral transfer was not examined separately. Furthermore, we experienced certain limitations in the spectral transfer performance that are addressed in the present work by adding an intermediate step to the process (the VAE part).
The contribution of the present work is a novel spectral transfer methodology based on deep artificial neural networks (ANNs) that can process large shared datasets, works between systems with reasonably different dimensions (from lower to higher), and allows to connect more secondary systems to the same primary system. In addition, we present a unique way for obtaining unprecedentedly large datasets of corresponding spectral pairs, which results in more robust models for the spectral transfer.
Fig. 2 Illustration of the experimental setup. All devices (laser, motorized stage, and both spectrometers) are synchronized using a digital delay generator (DDG). |
To ensure a large variability of training data, we selected a heterogeneous rock sample (rare-element granitic pegmatite, further described below) and measured it in a mapping regime (i.e. raster over the sample surface, obtaining one spectrum from each spot). Each row (spectrum) xi,primary corresponds to the row xi,secondary, i.e., there is a one-to-one correspondence between the measurements, which were obtained from a slightly different spot of the same laser-induced plasma.
The test sample was cut from a strongly heterogeneous rock, classified as granitic pegmatite (locality Maršíkov D6e). Pegmatites are vein rocks specific by an extreme degree of chemical and textural fractionation, resulting in large contents of chemical elements that are otherwise rare in the Earth's continental crust. These include especially e.g. Li, Be, F, Rb, Cs, Ta, and Nb. The studied rock sample is mineralogically dominated by three minerals: quartz, albite, and muscovite. In subordinate amounts, it contains Fe–Mn garnet (almandine–spessartine), beryl, and Nb and Ta-oxides (columbite, fermite, and microlite).
All LIBS measurements were performed on a Firefly LIBS system (Lightigo, Czech Republic). The Firefly was equipped with a diode-pumped solid-state laser (266 nm, 5 ns, 50 Hz, 5 mJ) used for sample ablation and plasma generation. The laser beam was focused onto the sample (30 μm spot size) and plasma emission was then collected using a wide-angle objective and through an optical fiber bundle (400 μm core diameter), where each of the fibers collected radiation from a slightly different part of the plasma. Plasma emission was transferred to both Czerny–Turner spectrometers. The primary system ranged from 245 to 407 nm with a resolution varying from 0.035 to 0.046 nm, and for the secondary system from 260 to 370 nm with a resolution varying from 0.023 to 0.032 nm (see Fig. 4). The slit size was 25 μm for both spectrometers. Samples were placed onto a computer-controlled motorized stage. Spectra were acquired for each sampling position resulting in a raster with 40 μm spacing for training and 20 μm for the test dataset, in both X and Y axes.
The collected data from both systems were separated into three datasets; the training dataset (used for model training), the validation dataset (used for hyperparameter optimization and model validation during the training process), and the test dataset (used for final one-shot testing). To ensure that the results are representative of the performance on unseen data, each dataset corresponds to a separate measurement (out-of-sample evaluation) of a hyperspectral image. Their respective dimensions are 560 × 560 (training), 266 × 500 (validation), and 500 × 500 (test). Note that the test dataset was measured on a different day than the training and validation datasets. The sample photo with highlighted locations for each of the datasets can be seen in Fig. 3.
Fig. 3 Photograph of the measured mineral sample. The highlighted regions correspond to the training, validation, and test datasets. |
Fig. 4 Mean spectra of the primary and secondary systems from the test dataset. Significant lines that are relevant to the sample composition are labeled. |
The training/validation/test datasets consist of 313600/133000/250000 spectra (about a 45–20–35% ratio) with 4062 and 3872 features for primary and secondary systems, respectively. For example, this means that the dimensions of the primary training dataset are and for secondary. The combined size of all the datasets is roughly 22 GB.
All spectra were baseline corrected by performing a sliding minimum operation followed by a sliding Gaussian mean smooth (window size 100 and smooth parameter 50).26 In addition, each dataset was individually mean-centered and scaled to feature-wise unit variance (standard scaling). This was performed to compensate for the significant shift in the baseline of the spectra from the training and the test datasets (see Fig. 5). Since the dataset's original mean and variance can be saved, this transformation is reversible and lossless. All results are presented in the unscaled form. Selected examples of raw and processed spectra are available in the Appendix.
Fig. 5 Mean spectra from the primary system before preprocessing. Comparison of test and training datasets. There is a notable shift in the baseline. |
We propose a two-step approach: using a VAE (see Section 4.2), we obtain a (low-dimensional) latent representation of the data (denoted as Lprimary). In the second step, which is repeated for every secondary system, we train an MLP (see Section 4.1) to map the Xsecondary spectra to the latent space Lprimary. The second step requires the one-to-one correspondence of training spectra between the two systems. The ground truth vector for the MLP loss function is obtained from coordinates in the Lprimary (as illustrated in Fig. 6). Finally, by combining the shared decoder part of the VAE with the newly trained MLP we get the desired mapping f′. Parts of the proposed model are further explained below, followed by the evaluation methodology and baseline models for comparison.
In this work, the model hyperparameters (architecture and learning) were optimized on the validation data (see Section 3) with a hyperband algorithm.33 The considered domains, as well as other fixed hyperparameters, were selected heuristically, using prior experience with processing of spectroscopic data, along with non-exhaustive manual experimentation. The MLP is composed of two layers. For the first (second) layer 2048, 1024, 512, and 128 (1024, 512, and 128) neurons were considered, and 128 (512) was chosen by the optimization algorithm as the best. The leaky ReLU activation function was used in every layer except the last one that utilized a linear activation. Mean squared error with L2 regularization was used as the loss function. For training, the Adam optimizer was selected with a learning rate of 1 × 10−4, which was optimized from the following options: 1 × 10−3, 1 × 10−4, and 1 × 10−5. The model was trained for 100 epochs with early stopping set on validation loss and a batch size of 128. To reiterate, the inputs to the model are the Xsecondary spectra and the outputs are mapped as closely as possible to the corresponding Lprimary embeddings.
We used a variational autoencoder (VAE)38 that differs from the standard AE by bottleneck regularization and sampling from the bottleneck, in order to achieve desired characteristics of the latent space. The aforementioned characteristics are continuity (two “close” points should give similar results) and completeness (for a chosen distribution of the latent space, any point sampled from it should give meaningful results). In practice, the VAE is trained to predict a distribution (defined by its expected value and its variance), rather than predicting the latent encoding directly. The distribution of the latent space is regularized to fit the normal distribution centered at zero with unit variance by using the Kullback–Leibler divergence.39
Similarly to the MLP, the model parameters were optimized on the validation data with a hyperband algorithm. The autoencoder is composed of five layers, mirrored around the bottleneck. For the first and last layer we considered the following options: 2048, 1024, 512, and 128, and out of these 1024 was chosen. For the second and fourth (second to last) layers we considered 1024, 512, and 64 neurons, from which 512 was selected. The bottleneck was optimized to 64 neurons out of 3, 8, 32, and 64. The leaky ReLU activation function was used in every layer except the last and the bottleneck, which utilized a linear activation. Mean squared error with L2 regularization was used as the loss function. The Kullback–Leibler divergence was scaled by 4-cycle linear annealing40 going from 0 to 0.5. For training, the Adam optimizer was selected. The learning rate was optimized to 1 × 10−4 from the following options: 1 × 10−3, 1 × 10−4, and 1 × 10−5. The model was trained for 50 epochs with a batch size of 128. Both, the input and the output of the model are Xprimary. Note that for the proposed model (MLP + VAE), the output (prediction) is primary, which is obtained from the input Xsecondary.
(1) |
Euclidean distance of the i-th spectral pair:
(2) |
Cosine distance of the i-th spectral pair:
(3) |
Additionally, we also calculate the relative squared error (RSE) as:
(4) |
Lastly, for the impact on the qualitative analysis, we train a k-means clustering algorithm on the Xprimary train dataset and compare its predictions on the Xprimary test and the corresponding primary test measurements (i.e. spectra transferred from the secondary to the primary system). A simple accuracy metric was used:
(5) |
(6) |
〈Euclid〉 | 〈cosine〉 | RSE | k-score | |
---|---|---|---|---|
VAE + MLP | 4658.72 | −0.9919 | 0.1219 | 0.9758 |
One-step MLP | 4201.44 | −0.9931 | 0.1021 | 0.9798 |
KNN baseline | 17437.70 | −0.9075 | 1.8272 | 0.4495 |
PLS baseline | 5141.14 | −0.9917 | 0.1346 | 0.9670 |
A comparison between the total emissivity map of the Xprimary dataset and the primary can be seen in Fig. 7. It should be noted that the spectra outside the interquantile range Q0.99–Q0.01 are aggregated together and displayed with the same color as the corresponding quantile value. This is performed to filter out the outliers present in the data (the same approach is repeated for each subsequent hyperspectral image).
The spatial distribution of the transfer error can be seen in Fig. 8. Ideally, the error would be spatially invariant; however, it is clear that this is not the case. The spectra with the highest error are assembled on the borders of distinct matrices. This follows from the under-representation of boundary spectra in the dataset and contamination of the emission signal by the previous measurements.
Qualitatively, the transfer was highly successful, and the sample topology is well-preserved in the predicted map (see Fig. 7b). However, some important features (spectral lines) are predicted imperfectly. For the applications related to the qualitative analysis, these errors are not significant, but for quantitative analysis they could be considerable. To further investigate the prediction error and its wavelength dependency, we selected three representative spots, each from a different matrix present in the sample (see Fig. 7) and compared the original and predicted spectra from these spots (see Fig. 9–11). Predictions of representative spectra depicted in Fig. 9 and 10 show a good performance of the transfer methodology, where only a minor discrepancy between specific line intensities is present. More significant reconstruction error can be seen in Fig. 11, especially between 390 nm and 410 nm, where some lines are missed in the prediction, and a couple of background wavelengths are predicted with negative intensity. This error is representative of most of the spectra of the same matrix. This is likely a consequence of the background signal change (specifically for the test spectra of this matrix), since the aforementioned line intensities are below the background in the training dataset (see Fig. 5).
Fig. 9 Blue spot, from Fig. 5 (coordinates 408, 167). Muscovite (ideally KAl3Si3O10(OH)2; composition: K 9.8 wt%, Al 20.3 wt%, Si 21.2 wt%, O 48.2 wt%, and H 0.5 wt%). The original spectrum is from the primary system and prediction is from the secondary system. |
Fig. 10 Green spot from Fig. 5 (coordinates 46, 221). Albite (ideally NaAlSi3O8; element contents: Na 8.8 wt%, Al 10.3 wt%, Si 32.1 wt%, and O 48.8 wt%). The original spectrum is from the primary system and prediction is from the secondary system. |
Fig. 11 Orange spot from Fig. 5 (coordinates 408, 365). Quartz (ideally SiO2; element contents: Si 46.7 wt% and O 53.3 wt%). The original spectrum is from the primary system and prediction is from the secondary system. |
Lastly, we present the results from the k-means experiment described in the methodology section. The hyperspectral images of the test dataset were clustered into four clusters given by the k-means model built on the training dataset. We compared predictions on the original primary test and predictions from secondary test spectra. Misclustered spots after the prediction are plotted in Fig. 12. The k-score was 97.578, which demonstrates great potential for qualitative applications of the methodology (including classification tasks). Furthermore, the results show the robustness of the proposed method and the ability to generalize when faced with previously unseen data (even if they are coming from different measurement and samples with different matrix/composition ratios). However, they also reaffirm the suspected limitation of predicting the spectra on the borders of distinct matrices.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d2ja00406b |
This journal is © The Royal Society of Chemistry 2023 |