Co-Orbital Sentinel 1 and 2 for LULC Mapping with Emphasis on Wetlands in a Mediterranean Setting Based on Machine Learning

: This study aimed at evaluating the synergistic use of Sentinel-1 and Sentinel-2 data combined with the Support Vector Machines (SVMs) machine learning classiﬁer for mapping land use and land cover (LULC) with emphasis on wetlands. In this context, the added value of spectral information derived from the Principal Component Analysis (PCA), Minimum Noise Fraction (MNF) and Grey Level Co-occurrence Matrix (GLCM) to the classiﬁcation accuracy was also evaluated. As a case study, the National Park of Koronia and Volvi Lakes (NPKV) located in Greece was selected. LULC accuracy assessment was based on the computation of the classiﬁcation error statistics and kappa coefﬁcient. Findings of our study exempliﬁed the appropriateness of the spatial and spectral resolution of Sentinel data in obtaining a rapid and cost-effective LULC cartography, and for wetlands in particular. The most accurate classiﬁcation results were obtained when the additional spectral information was included to assist the classiﬁcation implementation, increasing overall accuracy from 90.83% to 93.85% and kappa from 0.894 to 0.928. A post-classiﬁcation correction (PCC) using knowledge-based logic rules further improved the overall accuracy to 94.82% and kappa to 0.936. This study provides further supporting evidence on the suitability of the Sentinels 1 and 2 data for improving our ability to map a complex area containing wetland and non-wetland LULC classes. of images designed to implement decision rules. The tree is made up of a series of decisions that are used to determine the correct class for each pixel. The rules can be based on any available characteristic of the dataset. For example, additional elevation information was used in this study to assist the correct classification of vegetation types (i.e., forest–marshes). The classified images from previously implemented steps, along with the ancillary information (e.g., cloud mask, NDVI, NDWI, slope, elevation, and crop features), were used herein as variables.


Introduction
Land use and land cover (LULC) consists of fundamental characteristics of the Earth's system intimately connected with many human activities and the physical environment [1].Information on LULC is of key importance in environmentally or ecologically protected ecosystems or native habitat mapping and restoration (Council Directive, 92/43/EEC, 1992).Wetlands in particular represent one of the world's most important and productive ecosystems, having a critical role in climate change, biodiversity, hydrology, and human health [2,3].Wetlands include permanent water bodies, lands that remain completely dry over several months, and areas where water is below a dense vegetation cover, such as peat bogs or mangroves [4].Those also include important natural complex habitat types such as fresh water marsh and riverine forests, scrublands, as well as agricultural landscapes [5].Although freshwater wetlands cover only 1% of the Earth's surface, these areas provide shelter to over 40% of the world's flora and fauna species [6].As such, wetlands are internationally recognized as an indispensable resource for humans [2] providing a wide range of services that are dependent on water, such as freshwater, agricultural production, fisheries and tourism [7].
Despite their importance, wetlands are one of the most threatened ecosystems due to anthropogenic factors such as intensive agricultural production, irrigation, water extraction for domestic and industrial use, urbanization, infrastructure, industrial development and pollution [7][8][9].Many wetlands are under pressure due to the natural and anthropogenic climate change (namely, changes in rainfall patterns and temperature and extreme events [9], as well as changes in land use brought by increasing populations and urban expansion).Environmental concerns on the degradation of wetlands came to the fore during the Ramsar Convention [8].Over the last century, it is estimated that 50% of the world wetlands have disappeared, with an increased rate of 3.7 times that during the 20th and 21st centuries [10].Thus, mapping and monitoring their dynamics over time is of crucial importance, as well as in the broader context of quantifying the temporal and spatial patterns of land use/land cover (LULC) and of its changes [11].
Earth Observation (EO) offers a repeated and frequent coverage of the Earth's surface over long time periods, which is ideal for monitoring wetlands [8].This has resulted in EO becoming the preferred method for natural resource managers and researchers [12].Identification and characterization of key resource attributes allows resource managers to monitor landscape dynamics over large areas, including those where access is difficult or hazardous, and also facilitates extrapolation of expensive ground measurements for monitoring and management [13].LULC mapping using satellite or airborne images allows for short-or long-term change detection and monitoring in such vulnerable habits [14][15][16].
Use of Synthetic Aperture Radar (SAR) imagery has been highly effective in wetland mapping and, although they are less frequently used in land-cover classification studies than optical data, they can be an important alternative or complementary data source [17].Synergistic use of optical data with SAR imagery may enhance the wetland-related information.This is because SAR provides data associated with the inundation level, biomass, and soil moisture, complimentary to optical sensors' information [18].In this context, Sentinel-1 offers dual-polarimetric C-band data and Sentinel-2 offers a wide range of high resolution spectral bands, thus testing their synergistic capability is of crucial importance.
Although spectral information about ground objects is important in information extraction from remote sensing data, previous studies have indicated that texture features are also important for differentiating wetland classes [15,16,[31][32][33][34][35].Spectral features describe the tonal variations in different bands, whereas texture features describe the spatial distribution [36].The most commonly used texture features are first order measurements (minimum, maximum, mean, range, standard deviation, skewness, and kurtosis) and second order measurements (mean, angular second moment, contrast, correlation, homogeneity, dissimilarity, variance, and entropy).First order statistics are used to quantify the distribution properties of the images' spectral tone for a given neighborhood, while second-order statistics contain the frequency of co-occurring gray scale values and are calculated from the gray-level co-occurrence matrix (GLCM) [15].
The recent growth of EO technology has resulted to the launch of sophisticated instruments such as that of the Sentinels series from the European Space Agency (ESA).This has opened up opportunities for new techniques development which aim at improving our ability to map wetland ecosystems.The Sentinel Mission is part of Copernicus Programme for monitoring climate change, natural resources management, civil protection and disaster management.Copernicus open data are available at Copernicus Open Access Hub (https://scihub.copernicus.eu/).Sentinel-1 (S1) is a Synthetic Aperture Radar (SAR) mission, providing data regardless of weather conditions and cloud coverage.Sentinel-2 (S2) is a land monitoring mission part of the Copernicus Program that provides high-resolution optical imagery to perform terrestrial observations in support of land services.The mission provides a global coverage of the Earth's land surface every five days and at a spatial resolution of 10, 20 and 60 m, making the data of great use in studies related to land use/cover mapping and quantification of its changes.The Sentinel satellites can play a vital role in future land surface monitoring programs.Thus, exploring and evaluating the use of Sentinel data for wetland mapping is a thematic area of key interest and priority to develop.The improved discrimination capabilities offered by S1 and S2 and their effectiveness for wetland mapping combined with contemporary classification algorithms (e.g., Support Vector Machines (SVMs)) would therefore be important to be investigated.
This study aims at exploring the synergistic use between a range of spectral information products derived from Sentinel 1 and 2 and the SVMs classifier in evaluating their ability to map a complex area containing wetland and non-wetland LULC classes.In particular, the study objectives are: (1) to analyze a number of secondary derivatives produced from S2 data with the SVMs to evaluate their added value in mapping LULC and specifically wetlands; and (2) to investigate the suitability of S2 data and their synergistic use with S1 data with contemporary LULC mapping techniques (SVMs and knowledge rules) for LULC mapping with emphasis on mapping wetlands.

Study Site
The National Park of Koronia and Volvi Lakes (NPKV) lies in the Mygdonia basin, a semi-urban area located in a tectonic depression in northern Greece (Figure 1).The basin includes a large lowland area around the lakes that offers abundant soil for cultivation, as well as some mountain ranges at its borders.The rugged terrain creates a dense hydrographic network that ends up in the lakes.Climate is typical Mediterranean and annual rainfall ranges from 400 to 450 mm, distributed almost entirely during the winter season [37].The wide rainfall range causes floods during the winter and serious droughts during the summer.The area is characterized by relatively low temperature during winter, ranging between −10 • C and 17 • C, whereas summer is warm with temperature ranging between 12 • C and 42 • C. NPKV is one of the most important Ramsar wetlands of Greece [38].According to the Greek Biotope/Wetland Centre (http://www.ekby.gr/ekby/en/),many flora and fauna reproduce, nest, feed and rest in the wetland habitat.Two perennial plane trees between the lakes are characterized as "Monument of Nature" and provide shelter to numerous bird species.It is a significant habitat of structural and species diversity, also providing an important nesting and roosting site for many endangered bird species (e.g., Milvus migrans, Haliaeetus albicilla and Hieraaetus pennatus) [5].In addition, 19 amphibian and reptile species, 34 species of mammals (some of them are under protection such as Myotis bechsteini, Myotis blyth, Lutra lutra, Canis aureus, Lutra lutra, Canis lupus and Capreolus capreolus) and 24 species of fish (among them, the rare "Aspius aspius") live and reproduce within the area.
The wetland is protected by numerous national and international conventions due to its high environmental interest.It is also included in the European ecological network of protected sites "NATURA 2000".Koronia Lake is a very important ecosystem of Greece which is almost vanishing because of non-sustainable water management practice in the region.It faces serious environmental issues, such as decreasing water levels, deterioration of water quality and water salinization [38].The ongoing unsystematic economic growth of the area has resulted in water depletion and environmental degradation with serious social and economic impacts.

Datasets
Sentinel data were acquired from Copernicus Open Access Hub (https://scihub.copernicus.eu/).The selection was based on the following criteria: (1) low cloud coverage for Sentinel 2 images; (2) coincidence of at least one Sentinel 1 and Sentinel 2 image; and (3) seasonal coverage.According to these, the images selected were a Single Look Complex (SLC) Sentinel-1 (C-band) image captured on 2 August 2016 in Interferometric Wide Swath Mode (IW), a Sentinel 2 image also captured on 2 August 2016 and a Sentinel 2 image captured on 28 January 2016 (Table 1).
SLC product produces a 250 km swath at approximately 5 × 20 m resolution.The acquired imagery was obtained in dual-polarization mode at VV + VH.Sentinel 2 imagery was acquired in processing Level-1C which includes radiometric and geometric corrections, ortho-rectification and spatial registration on a global reference system with sub-pixel accuracy.

Datasets
Sentinel data were acquired from Copernicus Open Access Hub (https://scihub.copernicus.eu/).The selection was based on the following criteria: (1) low cloud coverage for Sentinel 2 images; (2) coincidence of at least one Sentinel 1 and Sentinel 2 image; and (3) seasonal coverage.According to these, the images selected were a Single Look Complex (SLC) Sentinel-1 (C-band) image captured on 2 August 2016 in Interferometric Wide Swath Mode (IW), a Sentinel 2 image also captured on 2 August 2016 and a Sentinel 2 image captured on 28 January 2016 (Table 1).
SLC product produces a 250 km swath at approximately 5 × 20 m resolution.The acquired imagery was obtained in dual-polarization mode at VV + VH.Sentinel 2 imagery was acquired in processing Level-1C which includes radiometric and geometric corrections, ortho-rectification and spatial registration on a global reference system with sub-pixel accuracy.

Pre-Processing
The acquired S1 and S2 data were pre-processed using the SNAP open source software and georeferenced in the WGS84 coordinate system and UTM projection, zone 34.
S1 data first had to split to the study site extend, de-burst the sub-swaths and apply the precise orbit file to offer the highest geometric precision.A Refined Lee 7 × 7 speckle filter was applied after suggestions [35] and the σ 0 outputs were terrain corrected using SNAP's "Range Doppler Terrain Correction" algorithm with the SRTM 1 arc-sec DEM.H-Alpha (H-a) decomposition [41] was included, allowing for entropy and alpha derivatives to be extracted from the data.Finally, all data were resampled at 10 m using a bilinear method.
S2 data pre-processing included atmospheric correction to convert the Top-of-Atmosphere reflectance values (TOA) to corrected Bottom-of-Atmosphere reflectance values (BOA).For the atmospheric correction, ESA's Sen2Cor plug-in was used.Sen2Cor is a processor for Sentinel-2 Level 2A product generation and formatting; it performs the atmospheric, terrain and cirrus correction of Top-Of-Atmosphere Level 1C input data.Sen2Cor creates Bottom-Of-Atmosphere, optionally terrain and cirrus corrected reflectance images, as well as Aerosol Optical Thickness, Water Vapor, snf Scene Classification Maps, and Quality Indicators for cloud and snow probabilities.Its output product format is equivalent to the Level 1C User Product.Subsequently, bands 1, 9 and 10 were removed from the dataset.Then, the image was resampled at 10 m using a bilinear method and was also subset to the study site extent.
The SRTM DEM (30 m) was processed using ArcGIS software.Topographic information, namely slope, aspect and elevation, were derived and then resampled at 10 m utilizing the cubic convolution resampling method to match the Sentinel-2 spatial resolution.The accuracy of the derived DEM and the co-registration accuracy with the Sentinel images was assessed by comparing on a cell-by-cell basis with some same scale reference vector data of the area.A positional accuracy within the sensor pixel range (i.e., <10 m) was achieved which was considered satisfactory.

Analysis of Spectral Features
Spectral bands of S2 from 2 to 8A, 11 and 12 were used for this study (Figure 2).The Principal Component Analysis (PCA) and Minimum Noise Fraction (MNF) transformation methods were implemented to decrease the high correlation between the spectral bands providing independent information.PCA is a classical statistical method for transforming attributes of a dataset into a new set of uncorrelated attributes called principal components (PCs) and it is used to reduce the dimensionality of a dataset, while still retaining as much of the variability of the dataset as possible [42].The MNF is used to determine the inherent dimensionality of image data, to segregate noise in the data and to reduce the computational requirements for subsequent processing [43].PCA components 1 to 3 and MNF components 1 to 5 were used in this study to compose the new image for classification (Figure 3).The used components contained over 98% of the original information.
In addition, the commonly used Normalized Difference Vegetation Index (NDVI) (Equation ( 1)) and the Normalized Difference Water Index (NDWI) (Equation ( 2)) were used in this study to help discriminate vegetation types and water surfaces (Figure 2).Bands 8 (NIR) and 4 (Red) were used to calculate NDVI according to the following equation: Likewise, bands 8 (NIR) and 3 (Green) were used to calculate NDWI according to the following equation:

Analysis of Texture Features
Apart from the tone (spectral variation), satellite images also consist of texture (spatial variation) [36].While spectral information is relatively easy to quantify, texture is more difficult as it involves measurements of variation in pattern, shape and size [44].Texture is described mainly by histograms, the gray-level co-occurrence matrix (GLCM), local statistics, etc., with GLCM being the most often used [32].In this study, four statistical indicators of texture information, namely homogeneity (HOM), dissimilarity (DIS), entropy (ENT) and angular second moment (ASM), were implemented (Figure 4).These indicators were selected as effective indicators for the texture description of different land cover types [35].Texture measures are influenced by window size since the scale of the spatial patterns measured, is dependent on window size [45].After several window size implementations, it was found that a window size of 7 × 7 pixels is suitable for the particular study site.The texture measures were applied on the NDVI index, PC3 and MNF-C5.

Analysis of Texture Features
Apart from the tone (spectral variation), satellite images also consist of texture (spatial variation) [36].While spectral information is relatively easy to quantify, texture is more difficult as it involves measurements of variation in pattern, shape and size [44].Texture is described mainly by histograms, the gray-level co-occurrence matrix (GLCM), local statistics, etc., with GLCM being the most often used [32].In this study, four statistical indicators of texture information, namely homogeneity (HOM), dissimilarity (DIS), entropy (ENT) and angular second moment (ASM), were implemented (Figure 4).These indicators were selected as effective indicators for the texture description of different land cover types [35].Texture measures are influenced by window size since the scale of the spatial patterns measured, is dependent on window size [45].After several window size implementations, it was found that a window size of 7 × 7 pixels is suitable for the particular study site.The texture measures were applied on the NDVI index, PC3 and MNF-C5.

Analysis of Texture Features
Apart from the tone (spectral variation), satellite images also consist of texture (spatial variation) [36].While spectral information is relatively easy to quantify, texture is more difficult as it involves measurements of variation in pattern, shape and size [44].Texture is described mainly by histograms, the gray-level co-occurrence matrix (GLCM), local statistics, etc., with GLCM being the most often used [32].In this study, four statistical indicators of texture information, namely homogeneity (HOM), dissimilarity (DIS), entropy (ENT) and angular second moment (ASM), were implemented (Figure 4).These indicators were selected as effective indicators for the texture description of different land cover types [35].Texture measures are influenced by window size since the scale of the spatial patterns measured, is dependent on window size [45].After several window size implementations, it was found that a window size of 7 × 7 pixels is suitable for the particular study site.The texture measures were applied on the NDVI index, PC3 and MNF-C5.

Analysis of Shape Features
Shape is one of the most useful features in a satellite image and is an important aid for the extraction of information about ground objects [32].The methods for discovering shape information about ground objects include perimeter, area and other shape measurements [46].From the visual interpretation of the images, it was found that the shapes of the crops were relatively small and rectangular, while the shape of natural vegetation had more irregular patterns.Therefore, crops can be discriminated from natural vegetation, besides their spectral differentiation, based additionally on the area and the near-rectangle shape.Thus, two shape indicators were used, namely rectangle fit and compactness.An image segmentation was preceded using the "edge method" with scale factor 40% and merge factor 80% to ensure that the segments would be right-side to avoid the confusion in the separation of land cover types, as well as the limitation the shape attributes of the features [47,48].The segmentation was applied to red (R), green (G), blue (B), near-infrared (NIR) and NDVI bands to both summer and winter images to identify seasonal changes in crops (Figure 5).

Analysis of Shape Features
Shape is one of the most useful features in a satellite image and is an important aid for the extraction of information about ground objects [32].The methods for discovering shape information about ground objects include perimeter, area and other shape measurements [46].From the visual interpretation of the images, it was found that the shapes of the crops were relatively small and rectangular, while the shape of natural vegetation had more irregular patterns.Therefore, crops can be discriminated from natural vegetation, besides their spectral differentiation, based additionally on the area and the near-rectangle shape.Thus, two shape indicators were used, namely rectangle fit and compactness.An image segmentation was preceded using the "edge method" with scale factor 40% and merge factor 80% to ensure that the segments would be right-side to avoid the confusion in the separation of land cover types, as well as the limitation the shape attributes of the features [47,48].The segmentation was applied to red (R), green (G), blue (B), near-infrared (NIR) and NDVI bands to both summer and winter images to identify seasonal changes in crops (Figure 5).

Analysis of Shape Features
Shape is one of the most useful features in a satellite image and is an important aid for the extraction of information about ground objects [32].The methods for discovering shape information about ground objects include perimeter, area and other shape measurements [46].From the visual interpretation of the images, it was found that the shapes of the crops were relatively small and rectangular, while the shape of natural vegetation had more irregular patterns.Therefore, crops can be discriminated from natural vegetation, besides their spectral differentiation, based additionally on the area and the near-rectangle shape.Thus, two shape indicators were used, namely rectangle fit and compactness.An image segmentation was preceded using the "edge method" with scale factor 40% and merge factor 80% to ensure that the segments would be right-side to avoid the confusion in the separation of land cover types, as well as the limitation the shape attributes of the features [47,48].The segmentation was applied to red (R), green (G), blue (B), near-infrared (NIR) and NDVI bands to both summer and winter images to identify seasonal changes in crops (Figure 5).

Crop Features Extraction
To extract the crop features a combination of the previously mentioned methods was adopted.The segmented image (Figure 5) was used and a range of permitted values for both spectral and shape features was defined (Table 2).The values used for the extraction are shown in Table 2.

LULC Classes
The classification scheme used in this study was developed in two stages.In the first stage, land use land cover (LULC) classes were defined which were representative of the scene's attributes for a particular image acquisition date (2 August 2016) (Table 3).In the second stage, a multi-seasonal approach was attempted to classify all important area's attributes, thus both summer and winter multi-spectral images were used.LULC classes were selected based on knowledge about the specific area.For this purpose, "Crops" were divided into "Summer Crops", "Winter Crops" and "Permanent Crops", and a new class named "Grassland" was implemented.A set of training points were selected representative of each class.

Support Vector Machines
Support Vector Machines (SVMs) is a nonlinear and non-parametric large margin supervised machine learning classifier implementing Vapnik's structural risk minimization principle [49].SVMs have several advantages in comparison to other hard classification approaches (for an overview of SVMs uses, see [50]).SVMs separate the samples of different classes by finding the separating hyperplane related to maximal margin minimizing the hinge loss function.Such solution guarantees a minimal generalization error.Essentially, the hyperplane is the decision surface on which the optimal class separation takes place.Intuitively, a good separation is achieved by the hyperplane that has the largest distance to the neighboring data points of both classes.Each training example is represented by a feature vector.To avoid computational overload, this is not done by evaluating all training points, but only a subset, called the "support vectors" of the algorithm.By using nonlinear kernel functions (e.g., Gaussian RBF and polynomial), the SVMs implicitly work linearly in a higher dimensional space, corresponding to a nonlinear solution in the input space, where the data naturally exist.Such mapping into the higher dimensional kernel space is implicitly performed by applying a kernel function k(•,•), evaluating the dot product between samples mapped in some higher dimensional space as φ(x i ) φ(x j ), where denotes vector transpose.For the standard binary SVMs formulation implemented in this paper, the hyperplane f (x) = w x + b optimally separating the two classes is found by minimizing: The slack variables ξ allow some training errors, guaranteeing robustness to noise and outliers.C corresponds to a user selected parameter to control the complexity of the model, acting as a trade-off parameter between nonlinearity and number of training errors.This quadratic optimization is solved by introducing Lagrange multipliers α to obtain the following dual form: When the optimal solution of the latter optimization is found, i.e., the α, labels of unknown samples x t are predicted by the side of the margin in which they lie, i.e., by the following expression: Note that standard SVMs are sparse in the α coefficients, so the final solution may be equivalently expressed only by the samples having a non-zero α.These samples are the ones lying on the separating margins f (x) = 1 and f (x) = −1.
To represent more complex hyperplane shapes than the linear methods, the techniques can be extended by using kernel functions.In this case, the problem transforms into an equivalent linear hyperplane problem of higher dimensionality.The use of the kernel function essentially allows the data points to be classified to spread in a way that allows the fitting of a linear hyperplane.Commonly used SVMs kernels include polynomial and radial basis function (RBF), and sigmoid kernels.In addition, SVMs also introduces a cost parameter C to quantify the penalty of misclassification errors in order to handle non-separable classification problems.In this study, RBF kernel was used for performing the pair-wise SVMs classification due to its promising capabilities compared to linear and polynomial [23].RBF kernel is defined from the following equation: Radial Basis Function: This kernel requires the definition of only a small number of input parameters (i.e., the C and γ parameters) and has also already been shown to produce generally good results in a range of classification studies (e.g., [29,30,51]).SVM optimum parameters of the classification implementation were established.The γ value was kept as suggested, 1/number of features.After several tests, the optimum C value was found at 2000 (Figure 6).From a C value of 1, overall accuracy rapidly rose up to 500 before beginning to plateau at 2000.After parameter selection, the entire scene was classified with the entire dataset tested.For the classification, different scenarios were implemented aiming at assessing the added value of the different derivatives generated in the previous processing steps.All scenarios included NDVI, NDWI and elevation, as their contribution is considered to be effective in recent studies [52][53][54].Initially (SB), only the original S2 channels were used.The first scenario (T) examined the performance of the transformed images without the original channels.In the second scenario (T + SB), the contribution of the transformed channels to the overall classification accuracy was examined using them along with the initials.The third scenario (SAR) examined the contribution of the S1 data derivatives.The fourth scenario (GLCM) examined the contribution of the texture characteristics resulting from the GLCM.Finally, in the fifth scenario (MS), a seasonal approach was examined, using features from both the summer and the winter images.The process was repeated until all components were examined and in each case, the channels that appeared to contribute to the increase in accuracy were also used in the following scenarios.

Post-Classification Corrections
Although LULC maps reached high accuracies, each scenario was found promising for different purposes.For example, the multi-seasonal approach (MS) contributed in the improvement of the classification accuracy in vegetation classes more, while sand/soil and impervious surfaces were better differentiated when the texture features were used (GLCM) (Table 4).A post-classification refinement was developed and applied using the LULC maps and ancillary information built in a hypothesis framework of Knowledge Engineer (ERDAS Imagine, 2016) to reduce classification errors.The process involves the definition of hypotheses and variables, as well as the implementation of true/false rules in a detailed decision tree (Figure 7).The variables in the decision tree classifier refer to a band of data, while the terminal hypotheses represent the final classes of interest.A decision tree is a type of multistage classifier that can be applied to a stack of images designed to implement decision rules.The tree is made up of a series of decisions that are used to determine the correct class for each pixel.The rules can be based on any available characteristic of the dataset.For example, additional elevation information was used in this study to assist the correct classification of vegetation types (i.e., forest-marshes).The classified images from previously implemented steps, along with the ancillary information (e.g., cloud mask, NDVI, NDWI, slope, elevation, and crop features), were used herein as variables.After parameter selection, the entire scene was classified with the entire dataset tested.For the classification, different scenarios were implemented aiming at assessing the added value of the different derivatives generated in the previous processing steps.All scenarios included NDVI, NDWI and elevation, as their contribution is considered to be effective in recent studies [52][53][54].Initially (SB), only the original S2 channels were used.The first scenario (T) examined the performance of the transformed images without the original channels.In the second scenario (T + SB), the contribution of the transformed channels to the overall classification accuracy was examined using them along with the initials.The third scenario (SAR) examined the contribution of the S1 data derivatives.The fourth scenario (GLCM) examined the contribution of the texture characteristics resulting from the GLCM.Finally, in the fifth scenario (MS), a seasonal approach was examined, using features from both the summer and the winter images.The process was repeated until all components were examined and in each case, the channels that appeared to contribute to the increase in accuracy were also used in the following scenarios.

Post-Classification Corrections
Although LULC maps reached high accuracies, each scenario was found promising for different purposes.For example, the multi-seasonal approach (MS) contributed in the improvement of the classification accuracy in vegetation classes more, while sand/soil and impervious surfaces were better differentiated when the texture features were used (GLCM) (Table 4).A post-classification refinement was developed and applied using the LULC maps and ancillary information built in a hypothesis framework of Knowledge Engineer (ERDAS Imagine, 2016) to reduce classification errors.The process involves the definition of hypotheses and variables, as well as the implementation of true/false rules in a detailed decision tree (Figure 7).The variables in the decision tree classifier refer to a band of data, while the terminal hypotheses represent the final classes of interest.A decision tree is a type of multistage classifier that can be applied to a stack of images designed to implement decision rules.The tree is made up of a series of decisions that are used to determine the correct class for each pixel.The rules can be based on any available characteristic of the dataset.For example, additional elevation information was used in this study to assist the correct classification of vegetation types (i.e., forest-marshes).The classified images from previously implemented steps, along with the ancillary information (e.g., cloud mask, NDVI, NDWI, slope, elevation, and crop features), were used herein as variables.

Support Vector Machines
Regarding the overall classification accuracy, the curve ascended from 90.83% to 93.85% and kappa from 0.894 to 0.928, which indicates a strong agreement with reality (Figure 8).Some exceptions were observed when the initial S2 bands were removed.S2 bands seem to have the greatest contribution to the classification results.The most significant contribution to the accuracy of the classification appeared to be the texture features derived from the GLCM analysis and the seasonal approach for the vegetation classes.Overall, the highest accuracy (93.85%) was reached when the NDVI and MNF texture bands (entropy, homogeneity, dissimilarity, and angular second moment) were used.On the other hand, the lowest accuracies were observed when the PCA components were used alone (82.50%) and when the H-a dual decomposition derivatives were implemented (91.37%).

Support Vector Machines
Regarding the overall classification accuracy, the curve ascended from 90.83% to 93.85% and kappa from 0.894 to 0.928, which indicates a strong agreement with reality (Figure 8).Some exceptions were observed when the initial S2 bands were removed.S2 bands seem to have the greatest contribution to the classification results.The most significant contribution to the accuracy of the classification appeared to be the texture features derived from the GLCM analysis and the seasonal approach for the vegetation classes.Overall, the highest accuracy (93.85%) was reached when the NDVI and MNF texture bands (entropy, homogeneity, dissimilarity, and angular second moment) were used.On the other hand, the lowest accuracies were observed when the PCA components were used alone (82.50%) and when the H-a dual decomposition derivatives were implemented (91.37%).4).
As mentioned above, each scenario's results were promising for different LULC classes.Individual class accuracies in vegetation classes decreased (up to ~20%) when S2 red bands (B5, B6, B7, and B8A) were removed.S1 derivatives contributed mostly to distinguish bare ground and impervious surfaces.Texture derivatives increased class accuracies for both vegetation (2-5%) and bare ground classes (4-7%).
As shown in Table 6, the main misclassified classes are "marshes" and "crops", followed by "swamps" and "urban".Producer's and User's accuracies (PA and UA) indicate that, while "swamps" ranged between ~70% and ~80% and "urban" between ~77% and ~93%, "crops" did not exceed 75% in any scenario tested."Marshes" appears to be the most unreliable class, as many pixels misclassified as marshes (mainly crops).UA in "marshes" class was lower than 70% in all scenarios while PA was almost 95%.This means that, even though almost 95% of the reference marsh areas have been correctly identified as "marshes", less than 65% of the areas identified as "marshes" in the classification were actually marshes.On the other hand, classes "water" and "forest" had the best results (>98.00%), as far as User's and Producer's accuracy (UA and PA) is concerned.

Post Classification
After the completion of the post-classification step, overall classification accuracy reached 94.82% and kappa coefficient 0.9362 (Table 7).All LULC classes were also classified with fairly high accuracy (for all classes, above 75%).For all vegetation classes, an accuracy above 80% was also reported (Table 7).The classes with the highest misclassification error are "urban", "sand" and "soil", while "forest" and "shrubs" appear to have some misclassification errors as well.Specifically, "shrubs", "urban" and "sand" classes appear more unreliable, with relatively high commission error, which means a relatively low percentage of each class's pixels actually represents the particular class on the ground.For "forest" and "soil", higher omission error was reported, which indicates that a relatively low percentage of each class's ground pixels also appears the same in the classified image.

Discussion
Accuracy assessment results suggest that the synergistic use of S1 and S2 data can achieve high classification accuracy, especially when combined with information on elevation, texture and shape.In a wetland area, the classes that prevail are water, aquatic vegetation and agriculture.In our case, water covers almost 15% of the area.The above, combined with the seasonality of both the lakes and the different types of vegetation, constitute a highly complex ecosystem in which the classification of LULC types may prove difficult.Moreover, weather conditions affect the quality of optical data (e.g., cloud cover) and the ground conditions (e.g., drought) affect the quality of SAR data, since it reduces the backscatter variations of the different classes due to their moisture difference.This is confirmed by both OA, UA and PA for specific classes.Wetland classes (especially marshes) were the most misclassified and the results were more accurate when texture and seasonality were considered.
Although vegetation classes have a similar spectral response, spectral signatures between some of the different LULC types can be distinguished in the red-edge area of the spectrum.That may explain the significant contribution of S2 red bands in the classification accuracy improvement.UA and PA for vegetation classes (marshes, swamps and crops) decreased significantly when S2 bands were removed (v2.1-T).Moreover, vegetation's phenological cycle helped the separation of different types in the multi-seasonal approach.A difficulty in distinguishing the aquatic and the inland vegetation was also observed.Texture features helped in this respect, as indicated by several previous studies [35,[56][57][58].Crops were successfully separated from natural vegetation using shape features, which can be easily explained as the cultivated land has more regular patterns and shape than natural vegetation areas.
Spectral features showed to be less effective in the separation of artificial, soil and sand classes, which was possibly due to their low spectral separability (refer to Figure 2).While the first image was acquired in mid-summer, the presence of moisture is significantly low on the ground which can lead to high reflectance values and, thus, confusion between artificial surfaces and natural surfaces with no vegetation (soil, sand).Specifically, according to data collected from the Meteorological Station of Lagkadas (http://www.meteo.gr/meteoplus/index.cfm), the cumulative rainfall height for the previous month was only 6.0 mm, with temperature ranging between 14 • C (nighttime) and 38.8 • C (daytime).These weather conditions cause severe droughts to the area during the summer (see Section 2.1).
Findings reported herein also align with the other studies conducted independently using different multispectral satellite data underlining as well the promising capabilities of SVMs [23,24,[59][60][61][62][63][64][65].Evidently, a proper parameterisation can highly affect the SVMs performance [66,67].The technique we chose to use in this study to parameterize the SVMs has been successfully implemented in other LULC investigations [27,68,69].RBF kernel was used due to its promising capabilities as suggested by previous studies [23,50,68,70].It is possible that another kernel function might have been more appropriate to the particular study site and data.SVMs algorithm achieved high classification accuracies, especially in wetland classes.
A number of studies have also indicated that SAR data may improve the classification accuracy, especially in extracting artificial classes (roads, buildings, etc.) [71,72] and change detection (wetland delineation, flood, etc.) [67,[71][72][73].However, in our study, S1 data slightly improved the classification accuracy, since the accurateness enhancement was about 1%.As mentioned above, weather conditions were not ideal at the specific time of the images acquisition, with high temperature and very low rainfall which had caused severe drought in the area.Although SAR data are only affected by wind speed, which in our case did not exceed 1.8 km/h (WSW direction), the severe drought conditions may have affected the data quality.Finally, the differences in the spatial resolution between the combined datasets and the image processing implemented to unify this spatial resolution could be another source of misclassifications reported herein.
Previous studies have indicated the improvement of accuracy in classifying land cover when knowledge rules are used [31][32][33].In our study, post classification results achieved very high accuracy and indicated the efficiency of this method in producing a high-accuracy LULC map (Figure 9).
Remote Sens. 2017, 9, 1259 15 of 19 accuracy, since the accurateness enhancement was about 1%.As mentioned above, weather conditions were not ideal at the specific time of the images acquisition, with high temperature and very low rainfall which had caused severe drought in the area.Although SAR data are only affected by wind speed, which in our case did not exceed 1.8 km/h (WSW direction), the severe drought conditions may have affected the data quality.Finally, the differences in the spatial resolution between the combined datasets and the image processing implemented to unify this spatial resolution could be another source of misclassifications reported herein.Previous studies have indicated the improvement of accuracy in classifying land cover when knowledge rules are used [31][32][33].In our study, post classification results achieved very high accuracy and indicated the efficiency of this method in producing a high-accuracy LULC map (Figure 9).More broadly, the proposed methodology includes image analysis, extraction of key characteristics (e.g., texture) and usage of additional information (e.g., elevation) to compose a new dataset that will be used in the classification process.In this respect, findings of our study provide important assistance towards the development of relevant EO-based products not only in Greece (where the method was tested herein), but potentially also allowing a cost-effective and consistent long term monitoring of wetland/non-wetlands at different scales based on the latest EO sensing technology.With the provision of continuous, reliable spectral information from the ESA's Sentinels 1 and 2 missions, the protraction of the proposed herein methodology has to some extend been assured.The proposed methodology may be easily transferred to be implemented at another location in terms of geographical scalability as well as in different wetland types (e.g., coastal wetlands, river estuaries, etc.).For this to be done successfully, the post-classification rules should be appropriately adapted when the method is implemented in different settings.The inclusion of this adjustment may both further improve the accuracy of the current technique and lead to a more detailed classification hierarchy.Finally, the present study also meliorates existing EU international figures, assisting efforts towards the development of operational products via the EU's Copernicus service.

Conclusions
This study is, to our knowledge, among the first to investigate the combined use of Sentinel with advanced machine learning algorithms for LULC mapping with an emphasis on wetlands, particularly in a Mediterranean setting.Results showed that spatial, spectral and temporal resolutions of S2 data are suitable for classifying wetland areas.Overall accuracy using S2 bands alone reached 90.83% and kappa coefficient 0.894, which indicates a strong agreement with reality.The inclusion of S1 data did not significantly improve classification overall accuracy (<1%), at least this was the case in our study.However, the inclusion of additional spectral (NDVI, NDWI, PCA, and MNF), texture (GLCM) and shape indicators improved classification results, in some cases not significantly.Highest accuracies were achieved when the texture features (3%) where included, while More broadly, the proposed methodology includes image analysis, extraction of key characteristics (e.g., texture) and usage of additional information (e.g., elevation) to compose a new dataset that will be used in the classification process.In this respect, findings of our study provide important assistance towards the development of relevant EO-based products not only in Greece (where the method was tested herein), but potentially also allowing a cost-effective and consistent long term monitoring of wetland/non-wetlands at different scales based on the latest EO sensing technology.With the provision of continuous, reliable spectral information from the ESA's Sentinels 1 and 2 missions, the protraction of the proposed herein methodology has to some extend been assured.The proposed methodology may be easily transferred to be implemented at another location in terms of geographical scalability as well as in different wetland types (e.g., coastal wetlands, river estuaries, etc.).For this to be done successfully, the post-classification rules should be appropriately adapted when the method is implemented in different settings.The inclusion of this adjustment may both further improve the accuracy of the current technique and lead to a more detailed classification hierarchy.Finally, the present study also meliorates existing EU international figures, assisting efforts towards the development of operational products via the EU's Copernicus service.

Conclusions
This study is, to our knowledge, among the first to investigate the combined use of Sentinel with advanced machine learning algorithms for LULC mapping with an emphasis on wetlands, particularly in a Mediterranean setting.Results showed that spatial, spectral and temporal resolutions of S2 data are suitable for classifying wetland areas.Overall accuracy using S2 bands alone reached 90.83% and kappa coefficient 0.894, which indicates a strong agreement with reality.The inclusion of S1 data did not significantly improve classification overall accuracy (<1%), at least this was the case in our study.However, the inclusion of additional spectral (NDVI, NDWI, PCA, and MNF), texture (GLCM) and shape indicators improved classification results, in some cases not significantly.Highest accuracies were achieved when the texture features (3%) where included, while the lowest when only selected transformed components were used (<0.5%).Cultivated land was successfully separated from natural using shape features.Finally, post classification techniques achieved very high overall accuracy (94.82%) with individual class accuracies more than 75% for all classes, thus providing an accurate map of the specific study site.
Wetlands are complex systems with a high presence of water and both natural and artificial surfaces.Thus, mapping wetlands and monitoring their changes over time is a challenging procedure, especially using multispectral imagery for this purpose.Sentinel missions can provide researchers and decision-makers invaluable data, with a satisfactory spatial and temporal resolution and, most importantly, free of charge.S2 proved to be appropriate for mapping such complex areas, thanks to the high spectral and spatial resolution of the data.Further work is required to verify the results obtained by testing the methods on other sites.The expansion of the techniques investigated herein using multi-seasonal as well as multi-annual datasets would also be another direction to be perhaps investigated next.Besides this, the use of other classifiers (e.g., random forests, artificial neural networks, etc.) would be an interesting avenue to be explored, as it may also provide information that can result in improving thematic information extraction accuracy from S1 and S2 satellites to map wetlands worldwide.

Figure 1 .
Figure 1.The National Park of Koronia and Volvi Lakes (NPKV) lies in the Mygdonia basin, northern Greece.The detail shows the land uses according to Corine Land Cover (CLC) 2012.

Figure 1 .
Figure 1.The National Park of Koronia and Volvi Lakes (NPKV) lies in the Mygdonia basin, northern Greece.The detail shows the land uses according to Corine Land Cover (CLC) 2012.

Figure 2 .
Figure 2. Spectral signatures of selected land use/land cover (LULC) classes (left).Vertical and horizontal axis represents the reflectance values (multiplied with scale factor of 10,000) of the Sentinel 2 summer image and the spectral bands used in this study in nanometers, respectively.Spectral mean values for Normalized Difference Vegetation Index (NDVI) and Normalized Difference Water Index (NDWI).Vertical axis represents the index values (right).The chart shows the values for the summer image (2 August 2016).

Figure 3 .
Figure 3. Spectral mean values of: Principal Component Analysis (PCA) (left); and Minimum Noise Fraction (MNF) (right) components.Vertical axis represents the transformed values in both charts.

Figure 2 . 19 Figure 2 .
Figure 2. Spectral signatures of selected land use/land cover (LULC) classes (left).Vertical and horizontal axis represents the reflectance values (multiplied with scale factor of 10,000) of the Sentinel 2 summer image and the spectral bands used in this study in nanometers, respectively.Spectral mean values for Normalized Difference Vegetation Index (NDVI) and Normalized Difference Water Index (NDWI).Vertical axis represents the index values (right).The chart shows the values for the summer image (2 August 2016).

Figure 3 .
Figure 3. Spectral mean values of: Principal Component Analysis (PCA) (left); and Minimum Noise Fraction (MNF) (right) components.Vertical axis represents the transformed values in both charts.

Figure 3 .
Figure 3. Spectral mean values of: Principal Component Analysis (PCA) (left); and Minimum Noise Fraction (MNF) (right) components.Vertical axis represents the transformed values in both charts.

Figure 4 .
Figure 4. Mean values of texture indicators.Vertical axis represents the texture values for Normalized Difference Vegetation Index (NDVI).

Figure 5 .
Figure 5.A representative scene from the segmented image with both summer and winter crops using R, G, B, Near Infrared (NIR) and Normalized Difference Vegetation Index (NDVI) channels.

Figure 4 .
Figure 4. Mean values of texture indicators.Vertical axis represents the texture values for Normalized Difference Vegetation Index (NDVI).

Figure 4 .
Figure 4. Mean values of texture indicators.Vertical axis represents the texture values for Normalized Difference Vegetation Index (NDVI).

Figure 5 .
Figure 5.A representative scene from the segmented image with both summer and winter crops using R, G, B, Near Infrared (NIR) and Normalized Difference Vegetation Index (NDVI) channels.

Figure 5 .
Figure 5.A representative scene from the segmented image with both summer and winter crops using R, G, B, Near Infrared (NIR) and Normalized Difference Vegetation Index (NDVI) channels.

Figure 6 .
Figure 6.C value selection graph for Support Vector Machines (SVM) plotted against the overall accuracy of the scene.

Figure 8 .
Figure 8.The figure above represents the overall accuracy.All Support Vector Machines (SVM)-based classifications were grouped into six groups according to the legend: SB, Spectral Bands; T, Transformations; MS, Multi Seasonal; GLCM, Texture Analysis and SAR, Synthetic Aperture Radar (as also shown in Table4).

Figure 8 .
Figure 8.The figure above represents the overall accuracy.All Support Vector Machines (SVM)-based classifications were grouped into six groups according to the legend: SB, Spectral Bands; T, Transformations; MS, Multi Seasonal; GLCM, Texture Analysis and SAR, Synthetic Aperture Radar (as also shown in Table4).

Figure 9 .
Figure 9.The final classified image after the post-classification corrections (PCC).

Figure 9 .
Figure 9.The final classified image after the post-classification corrections (PCC).

Table 1 .
Summary of the characteristics of the remotely sensed datasets used in this study.

Table 1 .
Summary of the characteristics of the remotely sensed datasets used in this study.

Table 2 .
Value range used for crop extraction.

Table 5 .
Training and validation samples for each land use/land cover LULC class collected as group of pixels.

Table 5 .
Training and validation samples for each land use/land cover LULC class collected as group of pixels.

Table 7 .
Overall accuracy (OA) and kappa coefficient after the post-classification and User's and Producer's accuracies (UA and PA) for all the classes.