STUDYING EVOLUTION OF HYDROTHERMAL ALTERATION MATERIALS IN THE TURRIALBA VOLCANO TROUGH MULTISPECTRAL AND HYPERSPECTRAL IMAGES

The aim of this work is to develop a geospatial methodology for the analysis of the time evolution of The Turrialba volcano using different automatic imaging techniques compared to expert-based remote sensing techniques. Change detection of hydrothermal alteration materials in relation with time series from multisensor data acquired in spectral ranges of the visible (VIS) and short wave infrared (SWIR) have been calculated. We used for this purpose multispectral and hyperspectral scenes of the Sentinel 2, ALI and Hyperion sensors, respectively, on four dates from 2013 and 2018. This work adopts a multi-source approach, applied to the analysis of the correlations between hydrothermal materials and spectral anomalies in The Turrialba volcano complex, located in The Central Volcanic Range (Costa Rica). An expert-based technique called Crosta’s technique for detecting hydrothermal materials have been applied. We have chosen four variables for generating a different Principal Component Analysis (PCA) for groups of channels, two highly reflective and two highly absorptive for each mineral. We have tested another technique to detect hydrothermal materials based on a discrete spectral profile analysis and an unsupervised data mining approach. In other sense, we have applied an automatic technique called anomaly detection to compare with the hydrothermal materials results. Results are presented as an approach based on a comparison of two different strategies whose main future interest lies in the automated identification of patterns of hydrothermally altered materials without prior knowledge or poor information about the area, which has relevant implications in image-based prospecting.


INTRODUCTION
Multispectral and hyperspectral imagery has been used in the last decades with satisfactory scientific results to study and monitor a huge number of different applications in the biosphere, geosphere, cryosphere and hydrosphere such as deforestation, fires, volcanic activities, landslides, earthquakes, environment, pollution, sea ice extent, soil moisture or flooding. Data and remote sensing technologies have been improving exponentially, in such a way that they provide information about the Earth with many details and high frequency. The European Earth Observation Program Copernicus has been specially developed for this purpose, one of its main tools being the Sentinel 2 twin satellites (Baillarin et al., 2012), forming a multispectral system of high resolution, both geometric/radiometric and temporal.
Imaging spectrometry, one of the most important imaging techniques, has proven in particular to be effective in the characterization of minerals based on statistical methods using specific reflection and absorption bands. Multispectral configurations in SWIR and VNIR (Visible-Near Infrared) have been successfully used to map hydrothermal alteration materials in different geological scenarios (Antón-Pacheco et al., 2001;Cipar et al., 2011;Crosta et al., 2003). There are some minerals that indicate hydrothermal alteration that we can detect remotely by their spectral responses in the diagnostic absorption and reflection bands. These specific spectral features can be extracted from properly calibrated remote sensing data Argillaceous minerals, such as kaolinite, illite and alunite, have a concrete spectral characteristic presenting a high reflectance between the wavelengths of 1.55 µm and 1.75 µm, they also present high absorption between 2,08 µm and 2,35 µm. Another characteristic of these minerals is that the rocks that have not suffered the hydrothermal process, usually present regular values for the previously mentioned wave lengths. In addition, the minerals with high content in Fe present a very high contrast between the wavelengths 0,63 µm and 0,69 µm and the wavelengths 0,45 µm to 0,52 µm.
In this research, we also rely on Anomaly Detection (AD), been developed in the last decades, allowing a better understanding of a domain that has the sub-pixel differentiation of the spectral mixture and its implications in anomalous responses (Duran and Petrou, 2007;Manolakis and Shaw, 2002;Malpica et al., 2008;Shaum, 2005).
Due to the increasing activity of this complex that requires special attention for its proximity to main urban areas, an important aim in this research is to establish relationships to link spectral patterns to surface phenomena in the area. Results are presented based on two different strategies. The main future interest of the approach in the automated identification of patterns of hydrothermally altered materials without prior knowledge or poor information about the area, which has relevant implications in image-based prospecting.

DATA AND PRE-PROCESSING
We used multispectral and hyperspectral scenes of the Sentinel 2, ALI and Hyperion sensors, respectively, in a period between 2013, 2017 and 2018. This work implies a multi-source approach, applied to the analysis of the correlations between hydrothermal materials and spectral anomalies in The Turrialba volcano complex, located in The Central Volcanic Range (Costa Rica). Data pre-processing techniques to extract information from all images have been applied. They include geometric correction, radiometric correction and a number of image arithmetic and statistical analysis. This previous step allowed to obtain a perfectly overlapped image data set, the identification of spectral diagnostic bands of hydrothermal alteration minerals and the spatial correlation accuracy.

Principal component analysis to detect hydrothermal alteration materials
There are some minerals that indicate hydrothermal alteration that we can detect remotely by their spectral responses in the diagnostic absorption and reflection bands. These specific spectral features can be extracted from the remote sensing data that are properly calibrated.
Argillaceous minerals, such as kaolinite, illite and alunite, have a concrete spectral characteristic presenting a high reflectance between the wavelengths of 1.55 µm and 1.75 µm, they also present high absorption between 2,08 µm and 2,35 µm. Another additional characteristic of those minerals is that the rocks that have not suffered the hydrothermal process, usually present regular values on the wave lengths previously mentioned. The minerals with high content in Fe present a very high contrast between the wavelengths 0,63 µm and 0,69 µm and the wavelengths 0,45 µm to 0,52 µm.
An expert-based technique called Crosta's technique for detecting hydrothermal materials have been applied. We have chosen four variables for generating a different Principal Component Analysis (PCA) for groups of channels, two highly reflective and two highly absorptive for each mineral (Crosta et al., 2004;Rejas et al., 2012).
For Sentinel 2 the channels used to analyse the clay were 2, 8, 11 and 12, while for ALI the channels were 3,7, 9 and 10, and for Hyperion after radiometric correction the channels were 14, 48, 150 and 206. We also defined diagnostic bands following the same criteria for a group of Fe minerals, and for this purpose the channels used of Sentinel 2 were 2, 4, 8 and 13, for ALI the channels were 3, 5, 6 and 9, and for Hyperion the channels 14, 31, 48 and 151. The subtraction of the clay bands plus iron oxides between the different dates, adopting 2017 as a reference, have been calculated. Once the statistics of the differences were calculated, we verified that the new variables fit the Gaussian bell, as expected.
Everything out of -3-sigma and + 3-sigma are the abnormal changes in the hydrothermal alteration material concentrations, the most interesting changes. In our case, it is observed that they appear in the crater of the Turrialba volcano and on the South-West slope, in addition to some cultivation plots in the South of the volcano (see Figure 3). Figure 3. Scatter plot between Ar+Fe between the differences for 2013 and 2018 taking 2017 as a reference (in green 3-sigma, red 2-sigma, blue 1-sigma, yellow 3+sigma).

Discrete spectral profile analysis
We have tested another technique to detect hydrothermal materials based on an unsupervised data mining approach. It uses the quantized spectra of the pixels of the multispectral/hyperspectral images, and in these discrete spectral profiles it finds local shapes that are shared by a significant number of pixels. In addition, it retains a shape only if it is coherent over space, meaning that if the shape is found in the spectral profile of a pixel, then it tends to occur also in the spectral profiles of the adjacent pixels. Let us consider the spectral profiles over the 13 channels of Sentinel 2, after quantization, a discrete spectral profile for a pixel can be, for instance, "1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 2, 2", where the intensity in each channel has been encoded using three symbols: 1 (low), 2 (medium) and 3 (high). A local shape is a sequence of consecutive symbols in this discrete spectral profile, as for example, "2, 2, 3, 3, 3, 3, 3, 3, 2". Among all possible shapes, a local shape is retained only if it occurs in at least S% of the pixels of the image (a minimum coverage constraint) and if when it occurs in a profile then, on average, it also occurs in K of the 8 neighbouring pixels (a minimum average connectivity constraint). Thus, a selected local shape denotes a regularity in the spectral profiles occurring in a coherent way over space.
The efficient extraction of these local shapes is based on a pattern mining technique (Méger et al. 2019) adapted here to handle spectral profiles. Examples of shapes that can be retrieved from these profiles, include (but are not limited to) plateau (flat zone), peak, valley. Indeed, the kinds of shapes are not predefined and can be of various lengths (i.e., shapes covering different numbers of channels), the only syntactic constraint is that a shape must range over channels that are adjacent in the spectra.
Each selected local shape is then visualized using an "occurrence map" where highlighted pixels represent the pixels containing the local shape in their spectral profile. For a given local shape, comparing its occurrence maps at two different dates provides insights about the evolution in the area. This is illustrated by Figure 4, showing the occurrence maps of two local shapes extracted in an unsupervised way from two Sentinel 2 images.
The first image is the one given Figure 1 (26/01/2017) and the second is an acquisition over the same area on 12/12/2018. The occurrence maps (left up, right up) are the ones obtained for the local shape 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3 ranging from channels 1 to 13, where the symbol 3 denotes high intensity in these channels.
These maps correspond to clouds in the images (26/01/2017 left up, 12/12/2018 right up). Another local shape obtained from the discrete spectral profiles in these images is 2, 3, 3, 3, 3, 2, 1, 2 ranging from channels 5 to 12, and its occurrence maps (26/01/2017 left down, 12/12/2018 right down) correspond to vegetation that do not change between the two dates (trees along brooks). These local shapes where obtained for S=1% and K=4.

Anomaly Detection
Detection of spectral anomalies aims at extracting automatically pixels that show significant responses in relation of their neighbourhood. Several methods have been developed in the last decades, allowing a better understanding of the relationships between spectral data dimensionality and the optimization of search procedures as well as the sub-pixel differentiation of the spectral mixture and its implications in anomalous responses (Duran and Petrou, 2007;Manolakis and Shaw, 2002;Malpica et al., 2008;Shaum, 2005).

26/01/2017 (left down Sentinel 2) and 12/12/2018 (right down
Sentinel 2). High anomalies in black, the big anomaly areas corresponding to clouds. A RX algorithm (Reed and Xaoli, 1996), widely accepted as a standard spectral anomaly detection, has been applied to all set of data. Higher concentrations of hydrothermal alteration minerals, in scenarios where the sources of error (mainly vegetation and clouds) are minimized, are correlated with the anomalies calculated in the spectral range reflective for Sentinel 2 and ALI. Spectral anomalies computed from Hyperion data mainly show the high signal-to-noise ratio presented in this particular scene.
Anomaly detection results for ALI, Hyperion and Sentinel 2 sensors are similar, but we have observed the effect of a pattern due to the acquisition process in the Hyperion case ( Figure 5, right up). High anomalies areas corresponding to clouds in this tropical spectral background. However, the main anomalies detected in relation with this study are those corresponding to the hydrothermal alteration pixels detected by the others image methods.
The correlation between RXD spectral anomalies and hydrothermal alteration, clays plus iron oxides have been carried out. In order to do that, it has been linearly adjusted, a sample space of 30 pairs of points, placing spectral anomalies in the Y axis and altered minerals in the X axis.
Independent samples of hydrothermal alteration and RX spectral anomalies have been statistically compared. Both samples have standardized kurtosis values outside the normal range. In the case of hydrothermal alteration, the check points corresponding to clouds are eliminated, the sample is adjusted to normality. While the non-normality of spectral anomalies is expected in a tropical scenario due to the anomalous response depends on the relationship between the spectral background and the target.
All regressions have been calculated at a confidence level of 95%, removing in each adjustment, the sample values that showed unusual residues and which correspond mainly with clouds and shadow. The results obtained are summarized in Table 1. For example, for 2017 (see Figure 6), since the P-value in the ANOVA table is less than 0.05, there is a statistically relationship significant between anomaly detection and Ar + Fe_hydrothermal alteration with a confidence level of 95.0%. The R-Square statistic indicates that the adjusted model explains 52.0811% of the variability in Anomaly Detection.

Models AD vs HT P-value Correlation R 2 (%)
The mean absolute error (MAE) of 20.7276 is the average value of the residuals. The residuals were examined to determine if there is any significant correlation based on the order in which it is presented in the data file. Since the P-value is greater than 0.05, there is no indication of a serial autocorrelation in the residuals with a confidence level of 95.0%.

VALIDATION AND DISCUSSION
A check of radiometric and geometric corrections of images have been performed from data measured in supervision campaigns. The selected training areas have been visited and for each one of them, the type of material has been identified. In August 2017 and July 2018 has been carried out a field survey, which allowed us samples of 29 and 30 points respectively, measured in field and laboratory with USB400 and ASD FieldSpec 4 Hi-Res spectroradiometers. The spectra from the laboratory and the ground truth data have been directly compared with the reflection apparent data on the sensor. This process has allowed to detect the spectral pure pixel in the training areas. The spectra have been used to characterize hydrothermal alteration materials and to check Sentinel 2 reflectance images by empirical linear regression.
Hydrothermal alteration for all images have been classified (see Figure 7) with an unsupervised K-Means algorithm in six clusters. It should be noted that the maps were obtained by unsupervised classification of hydrothermal alteration Crosta's images. Therefore, the result should be understood as a map with six classes supposed to correspond to six degrees of hydrothermal alteration.
The next step we have taken is to make a correspondence between the information of these groups with the real information of the use of the land for 2017 and 2018. For this validation we used data from the two years, in addition to our fieldwork experience in the field.
Well, therefore, the conclusion is that in Sentinel 2 maps of hydrothermal alteration in 2017 and 2018, group 4 corresponds in both cases to the highest concentration of hydrothermal materials, clays and iron oxides. The other groups are correlated with a mixture of vegetation cover (vegetation with different phenological states, water, clouds, shadows and certain artificial materials).
Note that these are not land cover maps. In this sense, the hydrothermal alteration maps would correspond (with an inherent error obviously) to a certain ground cover as follows: green-lush vegetation; light blue-vegetation; dark bluevegetation affected by acid rain; yellow-crop fields in medium or high phenological state, corpses; red-hydrothermal alteration materials, clays and iron oxides, and sum of clays and iron oxides; white-clouds and fumaroles, that are the amount of gases released into the atmosphere as a result of the volcano's activity.
The most relevant thing is to identify the group in red as the one with the highest concentration of hydrothermal materials. The rest of the colors and groups contains a lot of mismatches with respect to the ground cover. Hydrothermal alteration detected from Sentinel 2 in Turrialba and environment, is linked to crops and bare areas, as well as to artificial surfaces (roofs of buildings mainly). This does not mean that in the vicinity of these areas there are no hydrothermal alteration minerals. It can be observed that the detection is strongly influenced by the contribution of the broad, dense and homogeneous natural vegetation cover existing on the spectral information recorded by space sensors.
The evaluation of the hydrothermal alteration detection for Sentinel 2 in 2017 and 2018 images was performed with the method of SNAP programme. The Crosta's technique applied on the classification of this area for 6 cluster presented satisfactory results. Nine of the twelve validation points of hydrothermal alteration for 2017 are satisfactorily labelled such as hydrothermal alteration (75%). Eight of the nine validation points of hydrothermal alteration for 2018 are satisfactorily labelled such as hydrothermal alteration (88%). Dense vegetation is confusing with low vegetation and crops, and gullies are confusing with vegetation affected by acid rain.
The discrete spectral profile analysis presented Section 3.2, led to occurrence maps coherent with the clusters obtained on the Sentinel 2 images in 2017 (see Figure 8) and 2018. With parameters S=1% (minimum surface coverage) and K=4 (minimum average connectivity constraint) the technique selected 50 (resp. 48) local shapes from the spectral profiles in the 26/01/2017 (resp.12/12/2018) image. Figure 9. Check and test points for hydrothermal alteration (in red of K-means, right) for Sentinel 2 2017.
The occurrence maps provide additional evidences of the pertinence of the clusters. For instance, the maps given Figure  10 (2017 left up and 2018 right up) correspond to the red cluster (hydrothermal alteration, see Figure 9 for 2017) in the bottom left corner and is coherent with the evolution of this cluster in this area. In the more central part of the image (south-west flank of the volcano), the maps highlight parts of the red and yellow clusters and again follow their evolution from 2017 to 2018. These maps are those of the local shape 3, 3, 2, 2, 2, 2, 2, 3, 3, 3 ranging from channels 4 to 13, where symbol 3 (resp. 2) denotes high (resp. medium) intensity in the channels. This local shape, in the spectral profiles, is consistent with the spectral signatures of Ar + Fe hydrothermal alterations. A second pair of maps is given Figure 10 (2017 left down and 2018 right down). They correspond to local shape 3, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 ranging from channels 1 to 13, and are coherent with the evolution of the upper right part of the dark blue cluster (vegetation affected by acid rain and gullies).

CONCLUSIONS
The characteristic of high resolution data, both spatial and spectral, for natural covers has been studied by different remote detection methods using ALI, Hyperion and Sentinel 2 test data sets in 2013, 2017 and 2018 dates. This paper evaluates the performance of detection methods of hydrothermal alteration materials in scenes with similar backgrounds in the Turrialba volcano (Costa Rica) and it contains a comparison of two different strategies whose main future interest lies in the automated identification of patterns of hydrothermally altered materials without prior knowledge or poor information about the area. The two strategies lead to convergent results that exhibit hydrothermal alteration patterns that are coherent in space and also in their evolution over time.
Higher concentrations of some hydrothermal alteration materials in bared areas where the sources of error are minimized, mainly clouds and vegetation, are correlated with the spectral anomalies in the VNIR range. That is an interesting point because it might be possible to identify some particular spectral bands to carry on monitoring and also detect their environmental impact at the same time.
The spectral mixing associated directly with the spatial resolution has an impact in a significant way on the characterization of the tropical volcanic backgrounds, and thus in the detection of hydrothermal alteration materials used to monitor the dangerous consequence of the volcano activity.