VERY HIGH RESOLUTION LAND USE AND LAND COVER MAPPING USING PLEIADES-1 STEREO IMAGERY AND MACHINE LEARNING

Anthropocene is featured with increasing human population and global changes that strongly affect landscapes at an unprecedented pace. As a flagship, the coastal fringe is subject to an accelerated conversion of natural areas into agricultural ones, in turn, into urban ones, generating hazardous soil artificialization. Very high resolution (VHR) technologies such as airborne LiDAR or UAV imageries are good assets to model the topography and classify the land use/land cover (LULC), helping local management. Even if their spatial resolution suits with the management scale, their extent covers a few km, making large-scale monitoring complex and time-consuming. VHR spaceborne imagery has a great potential to address this spatial challenge given its regional acquisition. This research proposes to evaluate the capabilities of a Pleiades-1 stereo-satellite multispectral imagery (blue, green, red, BGR, and nearinfrared, NIR) to both model the surface topography and classify LULC. Horizontal and vertical accuracies of the photogrammetrydriven digital surface model (DSM) attain 0.53 m and 0.65 m, respectively. Nine LULC generic classes are studied using the maximum likelihood (ML) and support vector machine (SVM) algorithms. The classification accuracy of the basic BGR (reaching 84.64% and 76.13% with ML and SVM, respectively) is improved by the DSM contribution (5.49 % and 2.91 % for ML and SVM, respectively), and the NIR contribution (6.78 % and 3.89 % for ML and SVM, respectively). The gain of the DSM-NIR combination totals 8.91 % and 8.40 % for ML and SVM, respectively, making the ML-based full combination the best performance (93.55 %). * Corresponding author


Global and Local Change Monitoring
Facing worldwide changes at an unprecedented pace, the Earth system is today driven by human forces, that shift landscapes dominated by vegetation to agricultural then urban areas (Meyer et al., 1992). In addition, the effects of global warming are represented by the sea-level rise (IPCC, 2014), as well as the intensification of the cyclone/storm-induced wave and rain (Moussavi et al., 2011). These hazards are becoming ever more threatening as the population grows all over the world. Coastal areas constitute the flagships of the challenges tied to the Anthropocene era, because they host all the components of the erosion and flooding risks, exacerbated by an increasing exposure due to the population densification (Neumann et al., 2015).

Optical Data at Very High Resolution Over Regional Extent
The remote sensing and monitoring of the biosphere spatial patterns and processes require some specific technology, able to capture large-scale information but provided with high resolution. Satellite imagery has been efficient to track changes in air, terrestrial and ocean temperature, in sea level, or in land use / land cover (LULC) at a global scale (Lu et al., 2004). However, these data were constrained by a coarse-grained filter, as a logical trade-off between the spatial extent and resolution. Despite the meaningful products of the global normalized difference vegetation index, normalized difference water index, or rate of urbanization, the LULC planning mostly necessitates finer-scale information as claimed by the local managers and stakeholders, tasked with the spatial adaptation to anthropogenic changes. Solving global issues needs local management provided with very high resolution (VHR), that is to say equal to or less than 1 m pixel size.
Airborne imagery, such as the passive multispectral camera borne by manned or unmanned aerial vehicle (UAV), as well as the active light detection and ranging (LiDAR) system, suits well with the local expectations of the pixel grain (Yan et al., 2015;Collin et al., 2012). However, the spatial extent covered by airborne platforms remains either too small (less than a few km 2 for drones) or too costly (per km 2 for manned aerial campaign, given the aircraft charges and the time-consuming flight planning, Mury et al., 2019). The advent of the VHR satellite imagery has enabled to address the local demand by bridging the large-scale extent with the fine-scale resolution (Collin et al., 2013). IKONOS and QuickBird-2, launched in 1999 and 2001, respectively, were the first spaceborne sensors able to deliver 1 m pixel size (for the panchromatic band) within regional scenes. In addition to the spatial capabilities of these sensors, multispectral optical information has been collected: blue, green, red, (BGR) and near-infrared (NIR). LULC mapping at VHR has been improved thanks to the NIR, given the high absorbance by water and the high reflectance by vegetation .

Spaceborne 3D Land Cover and Land Use Mapping
More recently, the VHR satellite remote sensing dedicated to Earth thematic mapping was significantly enhanced by the stereo-acquisition of the following sensors WorldView-1, GeoEye-1, WorldView-2, Pleiades-1, Worldview-3, Superview-1, and WorldView-4, respectively launched in 2007, 2008, 2009, 2011, 2014(e.g. Collin et al., 2018. Beyond the spatial and spectral specificities of these state-ofthe-art sensors, their agility has allowed to sample regional extents from two (even three) incidence angles, thus offering the possibility to generate 3D point clouds and digital surface models (DSM) using the photogrammetry procedure (Collin et al., 2018). Leveraging submeter panchromatic bands, the accuracy of the vertical products has been quantified around 1 m, provided that calibration ground-truth measurements were implemented (Bagnardi et al., 2016).
The added value of the relief information to the spectral dataset for mapping VHR pixels or objects composing the LULC has been proven in data fusion, involving satellite optical imagery and airborne LiDAR information (Demarchi et al., 2016), or more recently, in photogrammetry-based unmanned aerial vehicle mounted with a multispectral sensor . However, such kind of data fusion has not been derived from a unique spaceborne sensor. Here, we propose to thematically map large-scale LULC at VHR using a satellite stereo-imagery, from which both the topographical and multispectral information will be extracted. The scene was captured by a stereo-imagery of the Pleiades-1 sensor over a coastal area along the French side of the Channel Sea ( Figure  1). The selected scene is a temperate complex area, displaying nine representative worldwide classes, such as tree, grass, soil, sand, mud, roof, road, freshwater, and seawater. Two supervised classifiers have been used: the maximum likelihood (ML) and the support vector machine (SVM). Four questions have been addressed: what is the contribution of the topography (i.e. DSM) into the standard BGR classification accuracy? What is the contribution of the NIR into the standard BGR classification accuracy? Is there a specific effect of the NIR and DSM contribution into the standard BGR classification accuracy? What is the influence of the supervised classifier into the previous contributions? taken with 18.04° incidence angle; (B) Image #2 taken with 18.95° incidence angle. Red rectangle represents the study site

Study Site
The study site is located along the French coast of the Channel Sea, called the Côte d'Emeraude (48.60°N, 2.00°W), between the bays of Mont-Saint-Michel and Saint-Brieuc, in Brittany (France; red polygon in Figure 1A and 1B). Stretching over a 64 km long coastline, this part of the Norman-Breton Gulf is hydrodynamically featured with a megatidal regime, reaching a tidal range of 14 m during the spring tide (Mahmoud, 2015). Composed of muddy-sandy estuarine, sandy shores and gneisses' cliffs, the studied part of the Côte d'Emeraude hosts the Frémur river and the bay of Beaussais ( Figure 2). This land/sea interface has attracted farmers, fishermen, shipowners and merchants for several centuries as attested by the presence of Saint-Jacut, Saint-Briac, and Saint-Lunaire medieval villages as well as the XIX-XX th seaside city of Dinard. Following the end of the World War II, the area has become increasingly attractive, which has mutated the territory from rural villages surrounded by crop fields bordered by hedgerows (still visible in hinterland) to seaside urban cities subject to high anthropogenic pressure on the coastline. Noteworthy is the increase in the population density, thus vulnerability, on exposed coastal areas, which exacerbates the effects of the ocean-climate hazards driven by the sea-level and temperature rise.

A B
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2020, 2020 XXIV ISPRS Congress (2020 edition)

Class Identification
Nine LULC classes were investigated given their representativeness and occurrence across the world (Table 1). A series of 1 000 pixels per class were manually selected in the form of spectrally-homogeneous polygons. This ground-truth dataset was splitted into 500 calibration and 500 validation pixels per class.

Airborne LiDAR Dataset
An airborne topographical LiDAR dataset was used as reference to calibrate/validate the digital elevation model (DEM) derived from the Pleiades-1A stereo-imagery. This dataset provides VHR information of the topographical characteristics of the study site with a vertical accuracy of 0.25 m and a horizontal point density of 10 points per m 2 . Collected in 2011 by the French national geographic institute, the point cloud has been interpolated in the form of a raster provided with 1 m pixel size, and constrained by the local RGF93 datum, Lambert 93 projection and IGN69 vertical reference. An array of 20 groundtruth points was built by focusing on ground features (neither trees nor built-up areas). This array was separated into 10 calibration and 10 validation XYZ points ( Figure 3).

A B
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2020, 2020 XXIV ISPRS Congress (2020 edition)

Satellite Stereo-Photogrammetry
Prior to the creation of the satellite-based DEM, some process steps needed to be realized. Panchromatic stereo-imageries were radiometrically corrected from digital number values to top-ofatmosphere radiance then reflectance values by accounting for the optical calibration factors and sun irradiance, respectively. The corrected panchromatic imageries were used to extract terrain features with a dense point matching algorithm, so as to compute a satellite-based DSM at 0.50 m pixel size (Xu et al., 2008). From the point cloud in LAS format, a projection change was carried out, transforming the X and Y coordinates in the WGS84 datum, UTM zone 30N to RGF93 datum, Lambert 93. The LAS tools utility enabled this transformation. Then, a conversion of the altitudinal parameters was also carried out in order to switch from the ellipsoidal elevation to an orthometric elevation from the new French vertical datum 2018 (RAF18). TcLAsConverter, a free tool developed by the Government of the Province of Navarra in Spain, has enabled this transformation to be implemented.
The root mean square error (RMSE) was calculated to validate the photogrammetry-driven Pleiades-1A DSM: RMSE= (1) where P = Pleiades-1 XYZ model points, O = LiDAR XYZ validation points, n = number of observations.

Pixel-Based Supervised Classification
The reflectance Pleiades-1 imagery (panchromatic and multispectral) provided with the smallest incidence angle (see Table 1) was subject to an orthorectification based on the rational polynomial coefficients. The radiometrically-/ geometrically-corrected dataset was then pansharpened so as to use the multispectral dataset at the panchromatic resolution (namely 0.53 m, since orthorectified). The accurate Gram-Schmidt procedure was performed (Collin et al., 2013). Based on the resulting BGRNIR at 0.53 m pixel size, 1 000 pixels per class were identified: 500 pixels for calibration, and 500 for validation. The classification validation was quantified using the overall accuracy (OA) and the class-level producer accuracy (PA), drawn from the confusion matrix.

Topographical and Spectral Contributions:
The DSM and NIR contributions to the nine-class LULC classification were computed from the BGR performance (OA). Following the BGR classification, were implemented the BGR+DSM, BGR+NIR, and finally BGR+NIR+DSM.

Algorithms Investigated:
Two supervised classification algorithms based on statistical learning theory were investigated: the common and fast ML, and the complexsolving but time-consuming SVM.
ML is a probabilistic method, assuming that the statistics for each class in each band are normally distributed, calculates the probability that a given pixel belongs to a specific class. Each pixel is assigned to the class that has the highest probability (Table 3).

Function Parameter Probability Threshold
Single value Data scale factor 1.00 Table 3. Description of the adjusted parameters with the probabilistic Maximum Likelihood (ML) classifier SVM is a non-probabilistic binary multiclassifier, that separates the classes with a decision surface (hyperplane), maximizing the margin between the classes (Table 4).

Function Parameter Kernel type
Radial basis function Gamma kernel 0.20 Penalty parameter 100.0 Pyramid levels 0 Classification probability threshold 0 Table 4. Description of the selected parameters with the Support Vector Machine (SVM) classifier

RESULTS AND DISCUSSION
Firstly, the DSM was reconstructed based on the Pleiades-1 stereo-imagery. Secondly, the contributions of the topographical DSM and spectral NIR to the classification accuracy of the nine LULC classes were evaluated.

Satellite-Based Digital Surface Model
The DSM, extracted from the stereoscopic information inherent to the panchromatic Pleiades-1 imagery, was in very good agreement with the LiDAR vertical validation ( Figure 4A, coefficient of determination of 0.99). Vertical (altimetric) and horizontal (planimetric) accuracies reached 0.65 m and 0.53 m, respectively. This result highlights the potential of Pleiades-1 stereo-imagery for the reconstruction of DSM with a high spatial resolution across a large scale ( Figure 4B). These figures corroborate the findings of a sub-meter coastal topography study, which obtained a vertical accuracy ranging from 0.35 to 0.48 m (highlighting some reconstruction problems due to shadows, Almeida et al., 2019). A DSM extracted from a tristereo Pleiades-1 could further improve the satellite-based topography, particularly for the shadowed tree and urban classes (Panagiotakis et al., 2018). The DSM product is thus very likely to provide additional information that spectral data cannot for mapping LULC (Sofia et al., 2015).

A
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2020, 2020 XXIV ISPRS Congress (2020 edition)

Overall Accuracy:
The classification performance of the BGR reference dataset provided a satisfactory OA of 84.64% and 76.13% ( Figure 5) with the ML ( Figure 6A) and SVM algorithm ( Figure 7A), respectively. On the one hand, the addition of the topographical predictor, DSM, increased the OA of 5.49% and 2.91% with ML ( Figures 5 and 6B) and SVM ( Figures 5 and 7B), respectively. On the other hand, the spectral predictor, NIR, augmented the OA of 6.78% with ML, Figure  6C and 3.89% with SVM, Figure 7C).  6D) and SVM (Figures 5 and 7D), respectively. The best classification combination consisted in the full dataset, followed by the BGR+NIR, then the BGR+DSM with the ML algorithm. Our findings reflect another study on the use of the complete Pleaides-1 dataset (panchromatic + multispectral) allowed to the best classification to be performed (Beguet et al., 2014).

Class-Level Accuracy:
The DSM and NIR contributions were assessed at the class-level for the nine LULC classes. Classification results of the four different combinations with both algorithms showed a class-level trend of a better performance when DSM and NIR were combined to BGR.
By analyzing the confusion matrix in detail using the PA, the DSM topographical predictor increasingly impacted the class detection with the ML algorithm: sand -3%, freshwater -1.8%, seawater 0%, road +0.8%, soil +5.2%, mud +6%, tree +10.6%, grass +15.2% and roof +16.4%. The SVM classification results were significantly correlated with those of the ML classification. These scores displayed the overall positive contribution of the DSM in addition to the basic BGR spectrum: road -53%, seawater +1.4%, soil +2.2%, mud +6%, sand +6%, tree +8.6%, freshwater +12.4%, grass +15% and roof +30.2%. For both algorithms, grass and roof classes obtained the largest DSM gains, which can be explained by the topographic component of these flat and high classes, respectively. The DSM can therefore help elucidate a topographical zoning .
Likewise, the contribution of the NIR was duly parsed. The gains with the ML were obvious: sand +0%, seawater +0%, road +1.6%, mud +2.4%, freshwater +4.4%, soil +6.8%, grass +11.2%, roof +11.6% and tree +23%. As for the SVM classification, the road class received a negative NIR impact of 18.4%. The discrimination of the other eight classes were increased by the NIR: seawater +1.4%, grass +2%, soil +2.2%, mud +3.4%, sand 4.2%, freshwater +9.4%, roof +14.4% and tree +19.4%. For both algorithms, roof and tree classes highly benefited from the NIR. The enhancement of the tree and roof classed by the NIR is greatly related to the chlorophyll pigments and silicate slates, which are strongly reflected in this part of the electromagnetic spectrum, beyond the visible .
Furthermore, our results can be corroborated with recent studies tasked with the LULC identification using multispectral UAVs . The UAV can indeed derive the DSM using the "structure-from-motion" photogrammetric approach (Casella et al., 2017), while leveraging a multispectral dataset (through dedicated Parrot, Sentera, MicaSense, etc. dedicated sensors), composed of BGR and NIR wavebands (Ahmed et al., 2017). Our on-going UAV findings have confirmed the Pleiades-1 results, namely that the LULC mapping is hierarchically increased by the DSM+NIR combination, followed by the NIR, and finally by the DSM (James et al., 2020).
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B2-2020, 2020 XXIV ISPRS Congress (2020 edition) Figure 8. Barplot of the Pleiades-1 DSM and NIR contributions to the BGR classification performance (Producer's Accuracy) built with the maximum Likelihood (ML) and Support Vector Machine (SVM) algorithms

Classifier Contribution
Given their specificities, the ML and SVM algorithms have been investigated to decipher their influence on the classification performance. The probabilistic ML systematically outperformed the non-probabilistic SVM overall findings for BGR, BGR+DSM, BGR+NIR, as well as the full combination. This examination was also true at the class-level (Figure 8). It can be hypothesized that the Pleiades-1 dataset over the coastal area is topographically and spectrally-simple enough to benefit the ML, rather than the SVM. These results suitably match with those derived from UAV studies on UAVs which also highlighted a better discrimination of coastal objects (beach, slikke, schorre) when the ML algorithm was built on the combination of spectral bands + topography (Le Poulain, 2019). It is also important to evaluate the algorithms based on their computation time requirement (Table 5) for the sake of transferability. Once again, the probabilistic ML surpassed the non-probabilistic SVM in terms of the computation processing for the four classifications. The full combination (BGD+DSM+NIR) needed for 6 minutes against 720 min for ML and SVM, respectively. The radial basis function used for the SVM kernel allowed the algorithm to solve for non-linear, complex and noisy data, but demanded more computation effort. In addition, the SVM computation processing increased with the pairwise classification strategy for multiclass classification, based on the combination of several binary SVM sub-classifiers.

CONCLUSIONS
This study shows the capabilities of an unique Pleiades-1 stereo-imagery to produce very satisfactory classifications of the coastal area composed of nine LULC classes, by (1) extracting an accurate DSM (0.53 and 0.65 m horizontal and vertical precision), by (2) combining the resulting DSM and NIR information to the BGR reference classification. The best classification performance was reached with the BGR+DSM+NIR dataset (93.55% of OA, Figure 9), followed by the BGR+NIR (91.42% of OA), then BGR+DSM (90.13% of OA), and finally the BGR (84.64% of OA). The best performances resulted from the use of the ML (compared to the SVM), which was, in addition, more efficient in computation time (6 min versus 720 min). This trend was also highlighted at the class-level: the DSM favored more the grass and roof classes; the NIR highly contributed to the roof and tree classes; the DMS+NIR significantly improved the roof and tree classes. Further research about the single use of the Pleiades-1 could evaluate the capabilities of a tri-stereo dataset to improve the identification of the LULC classes. Such examination could integrate the derived DSM and spectral bands (e.g. rugosity index, slope, normalized difference vegetation index), and object-oriented classification based on the convolutional neural network. Figure 9. Best classification score (BGR+DSM+NIR with the ML algorithm) draped on the DSM, both derived from a Pleiades-1 stereo-imagery