Potential of High-Resolution Pléiades Imagery to Monitor Salt Marsh Evolution After Spartina Invasion

Abstract

An early assessment of biological invasions is important for initiating conservation strategies. Instrumental progress in high spatial resolution (HSR) multispectral satellite sensors greatly facilitates ecosystems’ monitoring capability at an increasingly smaller scale. However, species detection is still challenging in environments characterized by a high variability of vegetation mixing along with other elements, such as water, sediment, and biofilm. In this study, we explore the potential of Pléiades HSR multispectral images to detect and monitor changes in the salt marshes of the Bay of Arcachon (SW France), after the invasion of Spartina anglica. Due to the small size of Spartina patches, the spatial and temporal monitoring of Spartina species focuses on the analysis of five multispectral images at a spatial resolution of 2 m, acquired at the study site between 2013 and 2017. To distinguish between the different types of vegetation, various techniques for land use classification were evaluated. A description and interpretation of the results are based on a set of ground truth data, including field reflectance, a drone flight, historical aerial photographs, GNSS and photographic surveys. A preliminary qualitative analysis of NDVI maps showed that a multi-temporal approach, taking into account a delayed development of species, could be successfully used to discriminate Spartina species (sp.). Then, supervised and unsupervised classifications, used for the identification of Spartina sp., were evaluated. The performance of the species identification was highly dependent on the degree of environmental noise present in the image, which is season-dependent. The accurate identification of the native Spartina was higher than 75%, a result strongly affected by intra-patch variability and, specifically, by the presence of areas with a low vegetation density. Further, for the invasive Spartina anglica, when using a supervised classifier, rather than an unsupervised one, the accuracy of the classification increases from 10% to 90%. However, both algorithms highly overestimate the areas assigned to this species. Finally, the results highlight that the identification of the invasive species is highly dependent both on the seasonal presence of itinerant biological features and the size of vegetation patches. Further, we believe that the results could be strongly improved by a coupled approach, which combines spectral and spatial information, i.e., pattern-recognition techniques.

1. Introduction

Biological invasions often represent an important source of change in a system. They are mostly perceived as a possible threat to local biodiversity [1,2,3,4]. Occasionally, the arrival of a new species can also carry significant benefits to ecosystem services [3,5,6,7]. The consequences of biological invasions are usually difficult to predict. Their impact on ecosystems may significantly change in space and time [8] because of a variable adaptive response of the host communities. This response can induce a slowdown of the invasive species (IS) proliferation or, conversely, a sudden rapid expansion. The detection and mapping of the earliest stage of IS development and expansion often turns out to be the best strategy for anticipating and controlling their impact [9,10].

Due to the frequency of physical disturbance, tidal flats and salt marshes present a high susceptibility to biological invasions [11]. This vulnerability particularly increases with human pressure [12]. These systems are characterized by complex interactions between geomorphological patterns and ecological structures [13] and present a highly dynamic and extensive character, which, being associated with access difficulties, make remote sensing an appealing monitoring tool [14,15,16].

Several studies have already demonstrated the strong potential of using optical satellite imagery in mapping intertidal vegetation [10,13,17,18,19], and particularly, for the monitoring of biological invasions [20,21,22]. Intertidal vegetation follows a natural spatial zonation, according to its tolerance to the different tidal levels [23,24]. This natural zonation represents an advantage for remote sensing assessment, because it provides a clear physical separation between different species. Vegetation tends to grow in patches of clearly demarcated borders, presenting a relative homogeneity from a radiometric point of view, and the meadows of each species tend to be constant in time [14]. However, despite the relative intra-patch homogeneity, different degrees of vegetation cover and mixing can occur when the pixel size increases [14,25]. The environmental background noise, due to the soil characteristics, water and biofilm, then generates a significant per-pixel variability of the upwelling radiative signal [26]. Further, there is also a significant time dependency on plant phenology that can largely influence species differentiation in spaceborne imagery [22,26,27,28].

For these reasons, it is essential to select the most suitable satellite on the basis of instrumental characteristics. Especially, remotely sensed multispectral satellite missions can provide data with significantly different temporal, spatial and spectral resolutions [25,29]. The selection of one resolution over the other needs to be weighted according to the objectives, as current orbiting sensors cannot yet provide simultaneous high temporal, spectral and spatial resolutions [30]. For instance, Belluco et al. [15] compared multispectral and hyperspectral remote sensing data efficiency in mapping salt marsh vegetation in the Venice lagoon (Italy). They showed that spatial resolution affected the pixel classification accuracy in a much more important way than spectral resolution and that hyperspectral data can contain largely redundant information. Consequently, they suggested that, for such studies, preference should be given to spatial, rather than spectral, resolution. However, a possible way to overcome these limitations is through the combination of images from sensors, providing different spatial and spectral products by the application of pan-sharpening techniques [31,32]. The application of such techniques enables one to simultaneously take advantage of panchromatic (higher spatial resolution) and multispectral (higher spectral resolution) images. The identification and differentiation of vegetation species is commonly based on a similarity in pixel clustering, using supervised or unsupervised methods [10,15,21,28,33]. Better results are most frequently obtained using supervised classification techniques [34]. However, in the case of muddy coastal zones, to achieve good results, these methods may require extensive and impracticable ground truth verifications.

Even with the current advances in remote sensing products and techniques, the differentiation of individual intertidal species, an important aspect in the assessment of biological invasions, is still a difficult task [10] and usually comprises large errors and confusion with other features. In this study, we aim to assess the potential of multispectral images, acquired by the Pléiades-1 satellite constellation, in order to discriminate the invasive cordgrass, Spartina anglica, from the remaining marsh vegetation in the Bay of Arcachon (SW Atlantic coast of France). The evaluation exercise is based on NDVI maps and the use of pixel classification methods—either supervised (Random Forests) or optimized unsupervised (simulated annealing) approaches—applied to Pléiades images, which were acquired in different seasons. Any advances in the accurate identification of the different salt marsh structures will enable a better understanding of the mechanisms and extent of colonization by invasive species and allow for the monitoring of its evolution, facilitating an early intervention if needed.

2. Study Site

The Bay of Arcachon is located in the South-western French Atlantic coast (Figure 1). This coastal lagoon is separated from the open ocean by the Cap Ferret sand spit. Its main flow entrance is divided into two channels, with sandy banks and shoals between them [35]. The presence of a mouth bar and inter-channel banks limits the internal propagation of external ocean waves. Therefore, inside the lagoon, waves are mainly wind generated and fetch limited, especially by the presence of the “Ile aux Oiseaux”. Wind speed and wind direction vary seasonally; they are variable and weak during summer and mostly westerly or northerly in winter, reaching speeds > 10 m/s. Tides are semidiurnal, ranging from 0.8 m during neap tides to 4.5 m during spring tides [36], and the high tide surface can reach an area of 174 km², mainly composed of a network of narrow channels, which is reduced to 40–50 km² at low tide [37]. The main source of freshwater is the Leyre river in the South-east, and the Porge Canal in the North. The total freshwater input represents less than 1% of the tidal prism [37].

For this study, we chose to consider a study area located on the intertidal zone at Andernos (Figure 1), on the north-eastern side of the Bay. The particularity of the marsh in this zone is that, even though it presents a stable niche of native vegetation, there is a significant presence of the exotic Spartina in small random patches. This is representative of an initial stage of the invasion and is an ideal scenario for our objective of identifying this species from its early stage of territory occupation.

Vegetation in the Bay of Arcachon

In the Bay of Arcachon, low elevations of the tidal flat (slikke) are mainly colonized by the seagrass, Zostera noltei, whereas the margins of the lagoon (schorre) are occupied by typical salt marsh vegetation (Figure 2a), principally Spartina sp. (Spartina maritima—Figure 2b–e and Spartina anglica—Figure 2f–j). The arrival of the invasive species, Spartina anglica, in the 1980s [38] provoked concern regarding threats to local species diversity and, especially, the eventual replacement of the native Spartina (Spartina maritima) by its exotic congener. Even though Spartinas constitute the dominant vegetation type in the marsh zone, there are many other species, such as Salicornia spp., Juncus gerardii, Limonium vulgare, Halimione portucaloides, and Puccinellia maritima [39,40] (Figure 2k–m). Usually, there is a well-defined separation between different vegetation species, but different degrees of mixture, as well as different degrees of density within the same type of vegetation, can also occur.

More than other marsh species, the invasive Spartina anglica has the particular ability to spread to offshore positions within the tidal flat, at the inner side of the lagoon. For this reason, the vegetated inner tidal flat zone was included, with the intention of identifying the isolated Spartina anglica patches that might occur in this zone.

Seasonal differences in the intertidal vegetation’s development and phenology can help to determine their identification and discrimination using remotely sensed images [22,27]. Notably, Spartina maritima is a species known for its continuous growth and absence of significant seasonal die-back [41]. In the Bay of Arcachon, this was verified by Proença et al. [42], who estimated a very low productivity for this species at different intertidal levels. On the contrary, the congener Spartina anglica presented a strong aerial productivity of aerial biomass, which means that it suffers from a strong die-back during winter, and it completely recovers during the growing season, particularly during summer, due to its exploitative growth strategy, which is better than that of the native. An important difference concerning the phenology of growth is also observable, with a maximum of growth during spring for the native and during summer for the exotics (Appendix A). Additionally, the cordgrasses seem to present some sensitivity both to continuously high and continuously low temperatures, expressed as a color change to yellow/brownish shades during summer and winter seasons, respectively. The seagrass, Zostera noltei, was also found to present a pronounced seasonality. Auby and Labourg [43] found a minimum aboveground biomass for this species between November and February and minimum shoot lengths in springtime. The maximum aboveground biomass was measured between June and August, but high values of shoot length were found during Fall. Concerning algae, several species can be found in the Bay at different times of the year [44]. According to these authors, in the inner parts of the Bay, there is a strong presence of the green alga, Monostroma obscurum, that tends to decline during the period between March and September. However, they found that a clear seasonal trend is difficult to predict.

3. Datasets and Methods

In this study, we used five high-resolution multispectral (MS) images, acquired by the Pléiades-1 satellite constellation, aerial and drone photographs, centered on the study area, Global Navigation Satellite Systems (GNSS) data, collected as ground truth validation, radiometric measurements of the features likely to appear in the classifications and georeferenced photographs.

3.1. Aerial Photographs

The aerial photograph dataset was obtained from French National Institute for Geographic and Forestry Information (IGN) database, and it comprised 41 orthorectified images, taken between 1934 and 2004, with a spatial resolution ranging from 0.10 to 0.85 m. From this dataset, 31 photographs were selected, based on visual evaluation, for their quality and potential for exploitation. They were then georeferenced. A subsection of the study zone was then extracted, and polygons of the identifiable Spartina meadows were obtained, distinguishing between the two species based on the GNSS data and georeferenced ground photographs. The number of identified Spartina patches (polygons) is variable, depending on the year and quality of the photograph. In order to consider the error associated with the visual determination of the Spartina meadows’ limits, three replicates of polygons were drawn for each patch. We then used the polygon surface of the three replicates to estimate the mean surface and standard deviation error of each patch.

The most recent aerial photograph was taken during a drone flight that was performed on 21 June 2016 over the study zone. We used a DJI Phantom 2 drone, equipped with a GoPro Hero 4 Black Edition, mounted on a DJI Zenmuse H3-3D gimbal to counterbalance unwanted camera rotations and absorb vibrations. The flight was performed at a very low altitude, 30 m, in order to obtain high-resolution images. Ground control points were taken over the study area and georeferenced using a GNSS system. The image dataset was aligned by a structure-from-motion algorithm through Agisoft Photoscan software, from which we obtained a 3D model that was optimized and georeferenced using the ground control points. Finally, a high-resolution (5 cm) orthophotography was extracted from the 3D-model. This high-resolution orthophotography was used both as a reference image and ground truth for the validation of feature identification in the satellite images.

3.2. GNSS Data

A field campaign was carried out on 26 October 2016 to determine the extent of the difference between plant species present in the study site using a differential Global Navigation Satellite Systems (DGNSS) technique. Two GNSS stations, composed of Leica GR25 receivers and AR10 antennae, were installed, one on the shore, as a base station, and a second one, mobile, carried in a backpack. The maximum distance between the two stations never exceeded 300 m during the field campaign. The GNSS station operated at multi-frequencies (L1, L2, L5) on multi-constellation signals. They were acquiring GNSS (GPS and GLONASS) data at a sampling rate of 1 Hz.

The RTK open-source program package (version 2.4.2) [45,46] was used to process the GNSS raw data. This freeware is available at: http://www.rtklib.com/. The satellite coordinate product provided the GNSS orbit and clock offset data used in this study, which were made available by the International GNSS Service (IGS) at ftp://cddis.gsfc.nasa.gov/gps/products/. We also used the Antenna Exchange Format (ANTEX) file (igs08_wwww.atx), made available by the International GNSS Service (IGS—https://igscb.jpl.nasa.gov/igscb/station/general/igs08.atx), which provides a precise and consistent set of phase center offset (PCO) and variation (PCV) values for both satellites and receivers’ antennae, taking into account these corrections. As the two stations were continuously acquiring data at both the L1 and L2 frequency bands, GNSS data were processed in ionosphere-free mode, with tropospheric corrections from the Saastamoinen model, and solid Earth and pole tides, followed by ocean tide loading, were applied [47].

3.3. Radiometric Measurements

On 21 July 2016, a field campaign was conducted in the study area to collect radiometric measurements, associated with the different types of substrate. Hyperspectral measurements were performed in the 400–900 nm spectral range (every 3 nanometres) from two radiometrically-calibrated TriOS-RAMSES sensors. The radiance sensor pointed downward to measure the upward signal of the substrate (L_u(λ), W·m⁻²·sr⁻¹·nm⁻¹), while the irradiance sensor pointed upward to measure the downward irradiance (E_d(λ), W·m⁻²·nm⁻¹). Radiance and irradiance measurements were performed simultaneously, as recommended in Reference [48]. The remote sensing reflectance (R_rs(λ), sr⁻¹) was then calculated as the ratio of L_u(λ) to E_d(λ). Note that, for brevity, we do not include the reference to the spectral dependency in the following text. Several features were chosen according to the different aspects of the vegetation of interest found within the study zone. For each feature, a run of 10 radiometric measurements was performed for three different targets in order to compute the average spectrum of R_rs.

Finally, in order to compare the field data with the Pléiades data, the in situ R_rs spectra were spectrally transformed using the spectral response function (SRF) of the Pléiades-1A and -1B bands in order to make them compatible with the spectral characteristics of the sensors. The SRF-corrected R_rs spectra (R_{rs_SRF}) were computed as follows:

$$R_{rs\_SRF}\left({\lambda }_0\right)=\frac{\int R_{rs}\left(\lambda \right)SRF\left(\lambda \right)d\lambda }{\int SRF\left(\lambda \right)d\lambda }$$

(1)

where λ₀ represents the central wavelength of the different MS bands of Pléiades.

3.4. High Resolution Pléiades Images

Standard orthorectified Pléiades images, acquired on five dates (25 April 2013, 3 August 2016, 6 October 2016, 24 Mai 2017 and 7 October 2017—

Table 1. Pléiades satellite data set information.

Satellite	Acquisition Date	Season	Acquisition Time (UTC)	Time of Low Tide (UTC)
Pléiades-1A	25/04/2013	Spring	11h15	12h04
Pléiades-1A	03/08/2016	Summer	11h15	12h45
Pléiades-1A	06/10/2016	Autumn	11h23	14h59
Pléiades-1B	24/05/2017	Spring	11h04	11h06
Pléiades-1A	07/10/2017	Autumn	11h08	13h27

), were used to determine the vegetation type in the study area. The Pléiades satellite constellation (Pléiades-1A was launched on 17/12/2011 and Pléiades-1B on 02/12/2012) acquired PAN (wavelength in the 470–830 nm range) images, with a ground sampling distance of 0.7 m in the nadir direction, which were resampled at 0.5 m and in four MS bands—blue (B: 430–550 nm), green (G: 500–620 nm), red (R: 590–710 nm), and Near-Infrared (NIR: 790–940 nm)—at a spatial resolution of 2 m and over a swath width of 20 km at nadir.

3.5. Pre-processing of the Pléiades Images

All images were not provided in the same coordinate system. Therefore, all images were first re-projected in the Universal Transverse Mercator (UTM), Zone 30 North—World Geodetic System 1984 (WGS84). Radiometric calibration and atmospheric correction were then applied to each image in order to convert the Level 1 (L-1) images, provided in digital counts (DC), to images in the surface reflectance values. Then, the surface reflectance was divided by π to obtain R_rs. Atmospheric corrections were performed using the optical calibration module from OTB—Monteverdi software [49], converting all images to the top of the canopy reflectance. Both sun and viewing angles (elevation and azimuth), as well as the calibration coefficients (gain/biases and solar illuminations), were obtained from the images’ metadata files. The applied aerosol model was chosen, depending on the corresponding season of each image, following the methods in Reference [50]. The atmospheric pressure, aerosol optical thickness and ozone and water content values were obtained from NASA’s Aerosol Robotic Network (AERONET) for the Arcachon station at the corresponding acquisition time of the satellite image. The atmospheric corrections in the blue band (B) present the highest residual errors. For this reason, and because the B-band does not provide useful additional information in the identification of vegetation, this band was removed from all images.

For each image, we computed the Normalized Vegetation Index (NDVI, [51]) as:

$$NDVI=\frac{R_{rs}\left({\lambda }_{NIR}\right)-R_{rs}\left({\lambda }_R\right)}{R_{rs}\left({\lambda }_{NIR}\right)+R_{rs}\left({\lambda }_R\right)}$$

(2)

Other indexes proposed in the literature for assessment of vegetation cover were tested, but they did not present significant differences relative to the NDVI. The NDVIs were also estimated for the field reflectance data, and the values obtained for each species were compared to the ones estimated from the satellite images. Based on the analysis of the NDVI maps, for each image, we performed a prior assessment of their capability of identifying the different types of habitats and species using classification methods.

3.6. Pixel-based Classification

To attempt to identify the two types of Spartina and differentiate them from each other, we tested pixel classification algorithms, using both unsupervised and supervised approaches.

Unsupervised methods provide the advantage of being more flexible [52], and this approach was considered with a view to extending the vegetation mapping to the entire Bay. Unsupervised approaches recognize patterns of spectral response using cluster analysis techniques. The resulting classification must be analyzed to regroup and label the classes using ground truth information. Even if there is a perspective of the long-term monitoring of salt marshes of the entire Bay, in this work, we will only focus on the results for the study zone.

Pixel-based unsupervised classification was performed using the Simulated Annealing Algorithm [53], a probabilistic method to find the class centers, which minimizes the cumulated distance from every pixel to their center of class in the spectral space. This is a non-linear algorithm that, with the appropriate parameters, provides a converged solution to the global optimized minimum. To improve the search for the best solution, a start temperature (a control parameter that determines the range of the initial search for the solution) is given, and its decrease, with the number of iterations, corresponds to a decrease in the probability of the acceptance of a worse solution than that previously found.

For the supervised classification approach, we considered different supervised classifiers, where preliminary tests showed that all techniques provided similar results. For this reason, we further proceeded with the Random Forests (RF) classifier [54], which is indicated in the literature as a suitable algorithm, when using a high number of variables and small training sets [4,10]. When applying the RF classifier, the classification map is accompanied by an accuracy map, which contains the indicator of the confidence degree of the classification. More precisely, the RF classifier is composed of a set of different individual base classifiers (decision trees). Each tree makes a decision, and the final decision corresponds to the majority vote. The degree of confidence is obtained as the percentage of the individual classifiers that provide a decision equivalent to the final one.

Pixel classification was performed on the three-band MS images. Supervised classification was performed using MATLAB^® on a standard personal computer, an Intel^® Core™ i5-4690 processor at 3.50 GHz, with 8 Go memory. For images with dimensions of 281 × 300 at a 2 m resolution, the algorithm’s computation time was under 5 minutes. As for the unsupervised classification, the simulated annealing algorithm was written in Fortran. It ran on a bi-socket server using AMD Opteron processors at 2.30 GHz, and it was parallelized on 24 cores using the OpenMP library. The longest computation was approximately 10 minutes.

All image processing steps are summarized in Figure 3.

4. Results and Discussion

4.1. Spectral Signature of Vegetated Structures

Figure 4 shows the mean R_rs spectra of different characteristic substrate types, collected at the study site on 21 July 2016. The spectral signature of the two Spartina species is very clearly distinguishable from other benthic habitat types. This difference is particularly obvious for dense and green Spartina meadows, without biofilm or sediment. Their R_rs spectra display a narrow peak center at around 550 nm and a sharp increase in the signal between the red and NIR bands, which is associated with high R_rs values over 780 nm. Further, Spartina anglica and Spartina maritima are also spectrally distinguishable. In green wavelengths, Spartina anglica shows a more pronounced peak than Spartina maritima, while in the red and NIR bands, Spartina anglica shows higher R_rs values than Spartina maritima. For instance, the R_rs value is 0.013 sr⁻¹ and 0.004 sr⁻¹ at 650 nm, and 0.071 sr⁻¹, and 0.042 sr⁻¹ at 840 nm, respectively (

Table 2. Mean (standard deviation) values of field hyperspectral, field SRF-corrected Pléiades, and Pléiades remote sensing reflectance (sr⁻¹) in the Green, Red, and NIR bands, and in the NDVI pseudo-band for the same dense and green Spartina anglica and Spartina maritima patches (see Figure 6). Field data were acquired on 21 July 2016, and the Pléiades image was acquired on 3 August 2016.

Spectral Band (nm)		Field Rrs (sr−1)		Field Rrs_SRF (sr−1)		Pléiades Rrs (sr−1)
		S. Maritima.	S. Anglica	S. Maritima	S. Anglica	S. Maritima.	S. Anglica
Green	560 500–620	0.007 ± 0.002 -	0.019 ± 0.003 -	- 0.006 ± 0.004	- 0.017 ± 0.004	- 0.015 ± 0.001	- 0.018 ± 0.001
Red	650 590–710	0.004 ± 0.002 -	0.013 ± 0.004 -	- 0.005 ± 0.006	- 0.014 ± 0.005	- 0.007 ± 0.001	- 0.009 ± 0.001
NIR	840 740–940	0.042 ± 0.008 -	0.071 ± 0.007 -	- 0.039 ± 0.009	- 0.064 ± 0.010	- 0.029 ± 0.003	- 0.034 ± 0.003
NDVI	$\frac{\mathit{N}\mathit{I}\mathit{R}\mathbf{-}\mathit{R}\mathit{e}\mathit{d}}{\mathit{N}\mathit{I}\mathit{R}\mathbf{+}\mathit{R}\mathit{e}\mathit{d}}$	0.83 ± 0.03	0.69 ± 0.02	0.77 ± 0.05	0.64 ± 0.04	0.60 ± 0.074	0.57 ± 0.061

). The NDVI, obtained from in situ Rrs, also indicates a significant difference between Spartina anglica and Spartina maritima, with a value of 0.69 ± 0.02 and 0.83 ± 0.03, respectively.

However, field measurements show a high variability of the spectral signatures of both Spartina species, depending on phenology, vegetation cover, biofilm and sediment deposition. For instance, the algae deposition on Spartina maritima causes a significant decrease in NDVI (0.54 instead of 0.83) due to an increased signal in the red and a decreased signal in the NIR. On the other hand, the sediment deposition on Spartina anglica shows a negligible impact on the radiometric signal in the red, but a strong impact in the NIR, which is responsible for a moderate decrease in NDVI (from 0.69 to 0.61).

Figure 5 provides the R_rs values in the green, red and NIR channels for the MS Pléiades images, acquired on 3 August 2016. The images zoom in on the Spartina anglica and maritima meadows, sampled during the field campaign conducted on 21 July 2016. MS images exhibit a high intra-patch variability of their spectral signatures, as previously observed using field observations. This variability is directly related to the spatial heterogeneity of the ground cover percentage in the meadows of Spartina species. In order to compare field and satellite remote sensing reflectances, we computed the field SRF-corrected R_rs from Eq.1 (Table 2). A general smoothing effect on R_{rs_SRF} due to the wide Pléiades spectral bands can be noted. However, small differences were recorded between R_{rs_SRF} and R_rs. The lowest and highest relative differences were found in the red band for Spartina anglica (7%) and Spartina maritima (22%), respectively.

A comparison of the satellite and field SRF-corrected R_rs shows that they exhibit a similar general spectral behavior, but with significant differences in the magnitude and spectral shape. For instance, the NDVI values associated with Spartina maritima were 0.77 and 0.60 for R_{rs_SRF} and Pleiades R_rs, respectively. Further, it is worth noting that the difference between the Spartina maritima and Spartina anglica NDVI values is significantly lower in the satellite data (0.03) than in the field data (0.13). This reduced capability of discriminating Spartina maritima from Spartina anglica using satellite observations may be related to the scaling factor [26]. Due to the small scale of the Spartina meadow and high spatial heterogeneity of the vegetation cover, the lower the spatial resolution, the lower the radiative contribution of Spartina to the per-pixel R_rs and the higher the contribution of the soil background signal.

Maps of the NDVI were generated for the five different dates of the Pléiades images (Figure 6). This index highlights strong seasonal changes, which appear to be primarily related to species phenology and the ephemeral presence of itinerant biological structures (e.g., the presence of biofilm and algae). The global high marsh zone can be clearly identified from the other types of habitat, particularly during spring (April and May), when the NDVIs of Spartina sp. vegetation show particularly higher values than sand and sediment. The results also highlight both seasonal and inter-annual variations of the NDVI signal for the Spartina vegetation. The NDVI does not, by itself, seem to be sufficient to differentiate the two Spartina species, as they present values of the same order of magnitude (Table 2). However, due to a small phase shift in their respective seasonal development, which generates considerable differences in aerial biomass density and color (Figure 2b–j), autumn appears as the best period to discriminate the two species, where average NDVI values for Spartina anglica and Spartina maritima are 0.8 and 0.9, respectively. Previous studies have already shown the advantages of a temporal assessment of marsh species identification using this index [55,56]. For instance, Sun et al. [55], considering the monthly NDVI of multiple marsh species in the middle coast of Jiangsu, found important temporal windows that made the species distinguishable.

Further, relevant information regarding the general intertidal vegetation presence during the growing season can be withdrawn from this parameter. In particular, remarkable differences appear between both spring images (April 2013 and May 2017). The NDVI values for both intertidal and supratidal vegetation are above 0.8 for the first date and much lower for the second one, ranging from 0.3 to 0.6. This indicates that potentially significant meteorological differences might have occurred within the four years, affecting the vegetation development, namely, in terms of temperatures and sunlight availability. However, such analysis goes beyond the scope of the present study, and further correlations with meteorological data would need to be conducted. For the image corresponding to the beginning of August 2016, at the peak of summer, NDVI values are no longer at their maximum, which could possibly be related to a particularly hot summer, during which high temperatures were reached. A relative senescence of the vegetation during this period was observed during field verifications, which is probably the reason for the marsh vegetation NDVI values being under 0.8 during this period of the year. Finally, for both fall images (2016 and 2017), not only the marsh vegetation structure shows particularly high NDVI values, but there is also the presence of a large vegetated cover over the tidal flat. This is probably related to the fact that, by this time of the year (October), temperatures have not yet significantly decreased, and vegetation growth is still occurring, for the seagrass in particular [43]. Indeed, it seems that this corresponds to an optimum period for Zostera noltei meadows to colonize the intertidal flat. Additionally, field inspections on this date allowed us to verify that this corresponds to a particularly favorable period for algae that considerably cover a large part of the tidal flat surface. Several deposits of such algae were also found to be hooked on to the stiff marsh vegetation, which might also affect its spectral signal (Figure 2d,h,i). It is then important to note that the presence of biological deposits can significantly affect the NDVI of the target marsh vegetation, and eventually, alter the typical values presented by the Spartina species.

4.2. Long-Term Evolution of the High Marsh Zone and Ground Truth Data Validation

The long-term temporal evolution of the study area, as viewed through aerial photographs (Figure 7), shows that its marsh surface has been relatively stable, especially since the 1980s, the period from which the high marsh zone presents a configuration similar to what can be observed on the present date. Indeed, the most significant changes occurred between 1949 and 1979, the period of remarkable expansion of the high marsh vegetation, dominated by the European native, Spartina maritima.

The main marsh structure (identified by the red square in the 2016 image in Figure 7) presented a growth of 5700 m² (±85 m²) between 1949 and 2016. The expansion of this patch seems to have progressively increased over time. For instance, its surface was enlarged by 563 ± 124 m² between 1985 and 1993, which is more than twice its growth between 1949 and 1958 (236 ± 67 m²). Between 2004 and 2016, the patch surface increased by 1613.4 (±80 m²).

Considering the exotic Spartina anglica, although its presence in the Bay of Arcachon dates from the 1980s, within our dataset, the most identifiable S. anglica patch in the study area (blue rectangle in 2016 image—Figure 7) has only been visible since 1993. This is the only patch that can be monitored, and it presented a growth of 100 ± 2 m² between 1993 and 2016.

Hence, in this area of the Bay of Arcachon, the most significant changes in marsh evolution are associated with vegetation patch growth (occasionally followed by the merging with neighbor patches). This patch growth and subsequent merging is most remarkable for the native Spartina. For the exotic Spartina anglica, we notice the appearance of new patches that seem to mostly occupy empty niches within the higher marsh zone or over the tidal flat (Figure A1). Therefore, apart from growth, there are no major inter-annual (or even inter-seasonal) changes in patch vegetation dominance. Our ground truth data (drone flight, GNSS data, georeferenced photographs and a regular monitoring of this site from 2015 up to 2018 within the scope of the complementary studies) for patch species identification showed that this is a relatively stable site, and no major changes have recently occurred. This provides a correct validation for the images from 2016 (Figure 8). Comparing these field data with the oldest image from our dataset (2013), a good match is observed with the knowledge of species composition and patch positions from the most recent data. For this reason, we consider that our ground truth data are also valid for the image of April 2013.

4.3. Pixel Classification Using Unsupervised and Supervised Methods

Several pixel classification algorithms were applied for each of the five dates of MS Pléiades images, with the aim of identifying the dominant marsh species in the Bay of Arcachon, Spartina maritima and Spartina anglica, and monitor the changes resulting from the invasion of the latter.

Concerning the supervised classification, as previously mentioned, after preliminary tests with several classifiers, we proceeded with Random Forests, both because it is a commonly used classifier, and it provided slightly better results than the other tested classifiers. For all classes, we considered the same number of training zones, maintaining, as much as possible, a balanced number of pixels per training zone, depending on the size of the feature representative of the class. For the image of April 2013, this balance was relatively difficult to achieve, as the presence of the invasive Spartina was particularly reduced in the study zone at this date. For the simulated annealing algorithm, since it is an unsupervised method, it is important to note that the algorithm class color attribution is done randomly, and class colors need to be re-assigned by the user in order to achieve an agreement between different season images and a match with the previously considered supervised classification results. Different tests were performed for each image, considering the influence of the algorithm’s parameters, number of iterations, start temperature (a parameter associated with the initial range of search for pixel clustering), number of classes and transformation coefficient in order to estimate the position of the class center. Very similar results were obtained for the classifications, with differences under 5% for each class from one run to the other. Field surveys have shown that the meadow vegetation dominance did not change at the considered timescale. Therefore, the alteration of the spectral signatures and consecutive class attribution by the algorithm was mostly influenced by the seasonal appearance of transient biological features or the vegetation’s phenology. It is then necessary to perform a post-classification analysis to provide the best color match between images of different seasons.

Both the unsupervised and supervised pixel classifications of the overall images allowed for an accurate identification of the global high marsh zone, including all species. However, single-feature identification (species differentiation) remained difficult, and the degree of accuracy was strongly dependent on the presence of noise induced by episodic events or short-term perturbations that can alter the species spectral signal (Figure 9;

Table 3. Accuracy metrics from the supervised Random Forests classification on the five dates of 3-band (G, R and NIR) multi-spectral Pléiades images. Class assignment: C1—Spartina anglica, C2—Spartina maritima, C3—slikke vegetation dominated by Zostera, C4—mudflat, and C5—biological deposits (algae/microphytobenthos).

	April 2013					August 2016					October 2016					May 2016					October 2017
Overall accuracy	0.87					0.66					0.79					0.85					0.56
Kappa index	0.81					0.56					0.72					0.80					0.46
Average accuracy	0.69					0.72					0.79					0.83					0.60
Class accuracy	C1: 0.64 C2: 0.93 C3: 0.97 C4: 0.88 C5: 0.037					C1: 0.9 C2: 0.96 C3: 0 C4: 1 C5: 0.76					C1: 0.72 C2: 0.88 C3: 0.96 C4: 0.39 C5: 0.98					C1: 0.85 C2: 0.74 C3: 0.98 C4: 0.96 C5: 0.64					C1: 0.53 C2: 0.78 C3: 0.39 C4: 0.35 C5: 0.92
Confusion matrix (in %)	64	0	36	0	0	90	0	0	3	7	72	8	0	0	20	85	2.5	10	2.5	0	53	0	47	0	0
	2	93	0	0	5	0	96	0.4	0	3.6	10	88	0	0	1	14	74	11	0	1	6	78	15	0	1
	1	0	97	0	2	0	26	0	0	74	0.2	3.8	96	0	0	0.4	0	98	0	1.6	60	0	39	0	1
	0	0	0	88	12	0	0	0	100	0	0	0	0	39	61	0.4	0	2.4	96	1.2	5	0	0	36	59
	0	0	18	78	4	0	22	2	0	76	1	0	0	1	98	0.5	7	28	0.5	64	0	0	0	8	92

). Such difficulties in the isolation of spectral features due to the presence of water and a high humidity at low marsh levels were also found, in [21], in the realization of supervised classifications of Landsat images. Given the variability of field reflectance for one matching feature (Figure 6), such a result was not completely unexpected, and this was shown by the significant differences between classifications for all dates. We observed that the larger the noise induced by seasonal perturbations, the larger the variability of the classification within the patches, particularly considering the results obtained using the unsupervised classification algorithm.

Globally, the Spartina maritima-dominated patches were identified well with the supervised classifier but comprised a high variability in class assignment, depending on the density of the meadow and the presence of biological debris with the unsupervised one. Small, scattered patches of the invasive Spartina anglica extending over the tidal flat were difficult to identify because of both the variety of interference signals and the small dimensions of the meadows towards the interior of the bay. A comparison of the ground truth of the dominant vegetation, identified through the GNSS contours, acquired on October 2016, and the corresponding classified vegetated features in the Pléiades images provided the rate of correctly classified pixels for the two Spartina species, obtained with an unsupervised classifier. Overall, patches of Spartina maritima dominance were identified well with the supervised classifier, with a percentage of correctly identified pixels of 74% for the image of May 2017, which increased to 96% for that of August 2016. With the supervised classifier, the strong presence of ephemeral biological perturbation on some dates (August and October) did not significantly affect the accuracy of the identification of the native species, and it was always over 74%. Using the unsupervised classification method, the correct classification of the native Spartina species attained 93% and 80% for the images of April 2013 and October 2017, respectively. The identification of the invasive Spartina anglica proved to be much more difficult, using both the supervised and unsupervised classifiers. Considering the results obtained using the supervised classifier, the identification of the ground truth patches was fairly successful, with accuracies ranging from 53% to 90% for the images acquired in October 2017 and August 2016, respectively. However, this class was frequently mistaken for the extensive zones of the vegetated mudflat (Zostera), particularly within the image from October 2017, and we obtained an overestimation of the areas assigned to the invasive Spartina. Belluco et al. [15] also found similar uncertainties in the identification of Spartinas that would be confused with the typical spectral signature of soil containing microphytobenthos, which contains significant quantities of chlorophyll. Indeed, the mudflat zone can either be vegetated in different density degrees by seagrasses or present deposits of itinerant features, like algae and microphytobentos, of several sorts. The combined presence of all these features provides complex spectral signatures that can be difficult to characterize. Thus, in this case, the only assurance we have that they do not correspond to the same species is field knowledge. We are confident that such a strong extrapolation could be reduced by increasing the number of ground truth samples of this species. Additionally, the accuracy maps of the classifications obtained for all dates (Figure 9–right panel) shows that regions classified by a higher uncertainty degree (accuracy uncertainty under 60 %) correspond to channels and zones in the tidal flat, where the mixing of features is high (vegetated and un-vegetated mud, water with microphytobenthos, sediment with different degrees of humidity, etc.). Considering the accuracy metrics provided in Table 3, only images corresponding to Springtime (April 2013 and May 2017) provide a Kappa index of 0.8, reinforcing that the lower level of biological noise within these images positively affects the performance of the classifier.

With the unsupervised algorithm (Figure 9), on the other hand, the correct classification of meadows dominated by the Spartina anglica species was between 10% and 25%. Such a low performance in the identification of the invasive species is probably primarily due to the fact that we consider several patches of small dimensions (2–4 m diameter). Secondarily, pixels within meadows of this species were often mistaken for pixels classified as vegetated mudflat and water with microphytobenthos, which is most likely due to the fact that Spartina anglica vegetation cover tends to be reduced outside of the growing season. Thus, there is a high ground exposure within zones colonized by this type of vegetation, which significantly affects its radiometric signal. However, when considering only patches of larger dimensions (5–10 m diameter), the identification increases to 69.5% for the unsupervised classification of the image acquired on May 2017.

The comparison of field species reflectance is measured in the NIR and Green bands (Figure 6; Table 2), and the classification results reveal an overall satisfactory class distinction, especially concerning the distinction of the native Spartina maritima from the remaining vegetation. Since, in the considered area, the presence of Spartina anglica was mostly reduced to isolated patches towards the inner bay or to zones of mixed vegetation in the higher marsh, the identification of the invasive species was not evident, particularly in images with a higher noise level. Considering the image of April 2013, where the vegetation did not yet present a strong seasonal development, several individual Spartina anglica patches presented a spectral signal under 0.01 sr⁻¹, which is much lower than the radiometric measurements. Indeed, the phenology of this species is rather late, and at this time of the year, the plants were still withered, and the meadows presented a very low density. Other studies [14,15] have also shown the interference of visible soil within vegetation in its reflectance, significantly lowering it. The classification of the native Spartina, on the contrary, provides center class values in the range of the radiometric field measurements, between 0.02 and 0.04 sr⁻¹ (NIR band), mostly depending on the meadow density. The better accuracy of the classification for the native Spartina was not so unexpected, considering that this species presents a much earlier development than that of the invasive congener.

The pixel-based classification allowed for the estimation of a surface of the main Spartina maritima patch of 5685.5 m², whereas with the manual identification of the patch, a surface of 5378.17 (±87) m² was obtained for the year of 2013. Considering the monitored patch of Spartina anglica, followed in the long-term aerial photographs, we obtained a surface of 67.25 m², using the supervised classification, and 74.31 (±1) m² using manual identification for the same date.

5. Conclusions

Technological evolutions of HSR multispectral satellite sensors greatly facilitates ecosystems’ monitoring capability at an increasingly smaller scale. In order to make progress in the knowledge of biological invasion dynamics in salt marsh environments, we evaluated the potential of Pléiades images to assess salt marsh evolution, after a biological invasion, and discriminate the native species Spartina maritima from the invasive species Spartina anglica for the study site of the Arcachon Bay (South-west France). Different techniques for land use classification were investigated, and their performance was evaluated using a large set of ground truth data. As expected, the supervised methods exhibit, on average, better results than the unsupervised algorithms. Moreover, in this study, we demonstrated that classification accuracy is highly dependent on two main factors, which are the size of vegetation patches and the species phenology.

Pléiades images provided a robust and consistent identification of the global marsh zone. We showed that, in the Bay of Arcachon, the arrival of the exotic species was relatively controlled and did not lead to significant changes in the marsh composition, which remained dominated by the native Spartina species. On the other hand, distinguishing between species remains challenging, particularly when considering small features, such as the invasive species at an early stage of the invasion, which represents relevant information concerning invasion control issues. Our results showed that the current pixel classification methods applied to Pléiades images provide a poor performance, when the order of magnitude of the Spartina sp. meadow size is lower than 10 m. However, lags in plant phenology provide temporal windows, making species distinction possible. The results highlight that a high-performance classification requires good knowledge of the field conditions and the vegetation’s development. The biggest challenge is revealed to be the dissociation of established vegetation from ephemeral biological features. In periods of low biological productivity, the level of biological noise can be significantly reduced, and hence the success in species differentiation increases. Thus, the seasonal interference of other biological features and stages of development of the target vegetation can alter its spectral signature or influence the exposure of the colonized environment (bare ground).

In the Bay of Arcachon, the two congener Spartina species present a phenology lag, where their spectral signature is considerably different between Autumn and early Spring. The strong presence of algae between late Spring and Autumn shortens this temporal window of remotely sensed assessment. We found that a suitable time frame for the remote sensing exploitation of the Spartina species distinction should be between late winter and early spring, when Spartina anglica is still withered enough to be differentiated from Spartina maritima, and algae blooms have not yet occurred.

The current pixel classification methods are highly sensitive to environmental noise. We highlight here that the consideration of particular temporal windows might be crucial for an improved assessment of the different types of marsh vegetation. Additionally, 2-meter spatial resolution data still generate significant uncertainties in the identification of vegetation patches characterized by a decametric scale. In order to better perform the identification of the invasive species, it will be necessary to evaluate the pan-sharpened products of HSR sensors, provided by Pléiades or WorldView 3 and 4, as well as coupled approaches based on pixel classification and pattern-recognition techniques.