Monitoring and simulating the distribution of phytoplankton in constructed wetlands based on SPOT 6 images

2.1 Study area

“Kristalbad,” a recently constructed artificial wetland which came into use in 2012–2013, is located between Hengelo and Enschede in the Netherlands (Figure 1). It is a complex and challenging water management project because multiple water functions and ecosystem services are combined in a very limited area. These functions are storm water retention, water quality improvement, ecological connectivity, recreation (walking, bird-watching, etc.), and landscape management. The input water in Kristalbad is largely effluent from the urban sewage treatment plant of Enschede–West. This effluent contains nutrients, organic matter, and other pollutants and flows into the Elsbeek.

Figure 1

Geographical location of the Kristalbad in Enschede, the Netherlands.

Enschede and Hengelo have an oceanic climate. However, due to its inland location, the winters are less mild in this area than that in the rest of the Netherlands. The annual mean temperature is approximately 3°C during winter and approximately 16°C during summer. In the Netherlands, rainfall is evenly distributed over four seasons. During the dry season (January–June), the rainfall is about 40–60 mm/month; while in the summer season (July–December), it is about 60–80 mm/month. During 1981–2010, the annual average precipitation and evaporation were about 780 mm (SD = 8.82) and 460 mm (SD = 34.2), respectively.

2.2 In situ data collection

2.2.1 Water sampling

A total of 21 water samples (3 × 7) were collected at three locations in the Kristalbad wetland during fall (September 1, October 1, and November 1, 2015, respectively). All samples were taken along the banks of the ponds, at a depth of less than 0.3 m, immediately after spectral measurements were made. The field campaign was planned to be started during the satellite overpass time (around 10:00 am) to improve the accuracy of the satellite data. Locations of the sample sites were recorded using a handheld Trimble GPS receiver. Samples SK1 and S6 were taken at the input and output of the Kristalbad wetland, respectively. Sample SK7 had a color very different to the other samples (Figure 2). All samples were used for retrieving C _chla, by combining the satellite data and spectral measurements obtained with a Trios/Ramses spectrometer.

Figure 2

Location of water sampling sites in Kristalbad.

2.3 Surface reflectance measurements

R rs ( λ )   was measured on the same day and at the same location as the water sampling (Figure 2), using TriOS Ramses radiometers. The greatest discrepancies between the R rs ( λ ) spectra and the TriOS system were mainly caused by environmental conditions during data collection. However, on October 1, 2015, the weather was ideal for spectral measurement because wind speeds were <5 m/s, cloud cover was 1.7% (nearly clear) in the study area, and the sun angle was high. At each sampling site, R rs ( λ ) was measured using a radiometric measurement system consisting of three TriOS Ramses hyperspectral radiometers, which measure detailed spectral properties in the wavelength range of 350–950 nm. The instrument has a variable spectral resolution from 3.26 to 3.35 nm.

To obtain R rs ( λ ) , we computed the ratio between the total upwelling radiance above water and the total downwelling radiance. The procedure of measuring R rs ( λ ) is described in previous research [24,25]. Three parameters were measured, (1) total downwelling irradiance E d ( λ ) , (2) downwelling radiance from the sky L sky ( λ ) measured by a radiometer looking upward at 42° from the zenith, and (3) total upwelling radiance L u ( λ ) at 42° from the nadir with an azimuth angle of 138°. All these angles were taken from a previous study [19].

The water leaving radiance L w ( λ ) measured at an angle of 42° from the nadir and an azimuth angle of 138° can be calculated by

In equation (3), ρ is the air–sea interface reflection coefficient, which is dependent on wind speed (W) in m s⁻¹. Mobley et al. recommended that when the wind speed is not measured or is <5 m s⁻¹, ρ = 0.028 is acceptable. Subsequently, R rs ( λ ) at θ = 42° relative to the zenith is derived as

2.4 Remote sensing data

2.4.1 Satellite for earth observation (SPOT 6)

SPOT 6 is a commercial, high-resolution optical imaging earth-observation satellite operating from space. SPOT 6 images have a 2-m resolution in the panchromatic band and a 6-m resolution in the multispectral bands. SPOT 6 contains four spectral bands with wavelengths from 450 to 890 nm, including blue, green, red, and near-infrared (NIR). The scene coverage was 60 × 60 km. SPOT 6 has both high-spatial and -temporal resolution and can revisit the same area every day with a constant viewing angle. This sensor, with its repetitive acquisition of high-resolution images, is very useful for monitoring land surface dynamics [26]. The example image shown in Figure 3a was acquired on October 1, 2015, the first day of the field campaign (images of September 1 and November 1 were also used but are not shown here). The water leaving radiance L w ( λ ) measured at an angle of 42° from the nadir and azimuth angle of 138° was calculated.

Figure 3

Example images of Kristalbad acquired on October 1, 2015. (a) SPOT 6 imagery, (b) water area of Kristalbad after data processing.

2.5 Regression modeling

2.6 Water quality modeling

We used the Dutch modeling system DUFLOW to model water quality. To model the constituents of Chl-a, the eutrophication model EUTRO1 was used. The growth of phytoplankton is affected by environmental conditions in Kristalbad, including the nutrient load, surface light intensity, and water temperature. Since the main species of phytoplankton in Kristalbad is algae, we used nitrate and ammonia, as suggested by DUFLOW, to indicate the growth of phytoplankton. This assumption is also justified as a compromise between accuracy and model complexity. It was calibrated for Chl-a using the data sets derived from SPOT 6 imagery. These datapoints were selected at seven nodes (except the start and end nodes) in Kristalbad. At these points, the C _chla data used the average value for each section in each pond.

The model was calibrated by comparing the measured data with the simulated model results. The calibrated model parameters established in the default ranges from similar studies [32], and a Chl-a-to-carbon ratio of 30 µg Chl-a/mg C was then adopted.

After model calibration, we had to define the water conditions for the simulation. Our sampling was conducted during the dry season in the Netherlands, and this was taken into account when simulating the scenario. In selecting the scenarios, we considered (1) the average flow discharge and water level during the dry season and (2) variations in the limiting factors of phytoplankton (light energy and water temperature). The light energy and water temperature were defined using time series data for the period September 1 to November 1, 2015, which contained the daily mean of hourly values. This data set was acquired from the Royal Netherlands Meteorological Institute (KNMI). In this short-term analysis, different parts of Kristalbad (ponds 1, 2, and 3) were examined to show phytoplankton growth dynamics.

3.1 Laboratory analysis of Chl-a

On October 1, 2015, the concentration of Chl-a at seven sampling stations ranged from 11.42 to 50.25 mg m⁻³, with a mean of 10 mg m⁻³. SK1 and SK6 were located at the inlet and outlet of the Kristalbad. SK7 is an obvious colored and turbid water sample, so the concentration of pigments is very high, with 50.25 mg m⁻³ Chl-a. The temperature, suspended matter, turbidity, and measured C _chla are described in Table 4. The Organization for Economic Cooperation and Development (OECD) offers a classification system for qualifying eutrophication levels, which provides the specific boundary value of the total phosphor, C _chla, and Secchi depth. Some water bodies can be classified into a trophic condition based on one parameter. In the absence of absolute trophic standards, the overlap of ranges still makes trophic classification schemes subjective. Therefore, the terms oligotrophic and eutrophic can have different meanings in different limnological situations. In this research, the average C _chla was 10 mg m⁻³ and the maximum C _chla was 50.25 mg m⁻³. Therefore, using the OECD classification, the Kristalbad wetland can be defined as eutrophic.

3.2 Remote sensing reflectance properties

After resampling the in situ R rs ( λ ) to R rs ( i ) , the reflectance values were mostly smaller than 0.004 except for sample SK7. Spectral features of the resampled reflectance were determined in response to C chla of the water samples. In Figure 3, the four-band spectral reflectance with varying phytoplankton C chla at each sample site shows a reflection peak in the green band, an absorption peak in the red band, minor reflectance in the NIR band, and minor absorption in the blue band. These features are roughly similar to those reported by Singh et al. [31]. In the blue region, there was no obvious absorption; this indicates reflection from suspended particles and the bottom. The maximum reflection occurred in the green band (∼560 nm), which was related to absorption by the phytoplankton chlorophyll and the effects of dissolved organic matter. There was substantial variation in the green band, due to the variable concentration of Chl-a; a high concentration leads to a stronger reflection. The low-reflection regions in the red and NIR bands were mainly due to the maximum absorption by phytoplankton chlorophyll. Therefore, the OC2 two-band model could be used to link spectral features to varying C chla .

3.3 Chl-a estimation algorithm

All possible relationships between the measured C chla and multiband ratios were calculated by regression. The performances of each equation fit with all possible linear and nonlinear forms, including linear, exponential, logarithmic, power, and polynomial, were comparable in terms of R. Results are shown in Table 5. We found that the linear function of the two-band model using blue (485 nm) and green (560 nm) bands showed significant sensitivity to the measured C chla values (Figure 4). The equation is as follows:

Figure 4

Spectral reflectance of resampled R rs ( λ ) for four SPOT 6 bands, with varying phytoplankton Chl-a concentrations, from sample station 3.

The best-fit linear function of the three-band model produced the largest determination coefficient and smallest estimation error, with approximate prediction results of R ² = 0.79, RMSE = 5.80 mg m⁻³ Chl-a, and a relative percentage deviation (RPD) of 2.13. Based on the previous research that used the OC2 model, we also observed a strong sensitivity to retrieve Chl-a in inland water. The model used remote sensing reflectance values in the blue (485 nm) and green (558 nm) bands, which were obtained from high-resolution FORMOSAT-2 image data. We used blue (485 nm) and green (660 nm) reflectance values, owing to the location of SPOT 6 image bands. Values were acquired from the resampled spectral field data. Although the remote sensing reflectance was obtained from a different source, the results were very similar and could be applied to SPOT 6 to estimate the Chl-a spatial distribution.

3.4 Accuracy assessment of SPOT 6 Chl-a mapping

The best-fit OC2 model was applied to SPOT 6 image to map the spatial distribution of Chl-a, as shown in Figure 5. This calculation was performed using the ERDAS model maker module. On the Chl-a image, the water area displays different C chla values, with the bank of the ponds having a higher level (>20 mg m⁻³) than the center (<5 mg m⁻³). Kristalbad is an optically shallow wetland, and most of the water samples were taken close to the bank of the ponds. Therefore, the total remote sensing reflectance obtained from SPOT 6 could be contaminated by the vegetation on the banks, shadows, and the wetland bottom, yielding inaccurate inversion results. In most research using an empirical algorithm to estimate C chla , the bottom albedo is removed from the satellite image, but the cost is excessive [33,34]. To balance cost and accuracy, the proposed method does not involve removing the bottom albedo from the remote sensing images, but the resulting C chla is still similar to the measured results.

Figure 5

Example map of Chl-a of Kristalbad derived from SPOT 6 image taken on October 1, 2015.

We performed field measurements during SPOT 6 overpass time. An example accuracy assessment was performed by comparing the image-derived C _chla with in situ measured C _chla (Figure 6a), which were matched using GPS coordinates. The accuracy assessment was applied to seven samples of the imagery-derived and laboratory-measured Chl-a. The difference (Δ), RMSE, and RPD were calculated for three dates and averaged and are shown in Table 6. In this table, the samples show a large difference, possibly indicating the influence of other samples, or other accessory pigments, which could affect the absorption measurement in the area. The overall RMSE for the inversion was 6.13 mg m⁻³ while the RPD was 2.25. In Figure 6b, a strong correlation was observed between the imagery-derived and laboratory-measured data (R ² = 0.8). However, differences were still noticed because some uncertain factors were included.

Figure 6

Comparison between Kristalbad C _chla measured in laboratory and that derived from SPOT 6 imagery. (a) Comparison in situ, (b) the correlation of them, with R ².

The differences between the measured and imagery-derived C chla were due to the following: (1) strong absorption in the red region, which greatly reduces the signal to the sensor and thereby affects the SNR; (2) the R rs ( λ ) magnitude was affected by the total scattered radiance, which is directly related to the type of the suspended sediment; and (3) differences between SPOT 6 imagery-derived R rs and field spectral measurements caused by radiance (atmospheric) correction and time interval. The effects of the first two points are related to the specific inherent optical properties (IOPs) and environmental conditions of the Kristalbad water, so they must be scrutinized and reinitialized. For the influence of the third point, the effect of the time interval is negligible because the field measurement was conducted during SPOT 6 overpass time. Therefore, the deviation between SPOT 6 imagery-derived R rs and field spectral measurements was mainly caused by the radiance correction, in this case mainly the atmospheric correction. Hence, the efficiency of the inversion could be improved with better atmospheric correction. Other uncertain factors (specific IOPs and environmental conditions) could also be investigated.

3.5 C _chla modeling integrated with remote sensing

The C chla of Kristalbad, from the inlet flow from the Elsbeek stream to three wetland ponds and then the outlet, was simulated using a 1-D DUFLOW water quality model. Remote sensing, in the form of estimated Chl-a data, was used to calibrate the Chl-a model. Calibrated data were obtained from SPOT 6 image using nine estimated data sets (four for pond 1, three for pond 2, and two for pond 3). The Chl-a model results show similar trends of concentration and were derived from SPOT 6 image (an example for pond 1 is shown in Figure 7). Each pond (pond 1) had three sections, and each of which exhibited different Chl-a characteristics. From the section 1 to section 2 (50 to 250 m) of pond 1, C chla increased from 10 to 15 mg⁻³ and then remained steady at approximately 15 mg⁻³ from sections 2 to 3 (250 to 400 m). In the last section (400 to 800 m), the value increased dramatically from 15 to 25 mg⁻³.

Figure 7

Comparison between modeled and imagery-derived Chl-a concentration in pond 1.

Figure 8 presents the spatial distribution of Chl-a simulated in a steady state, for the three ponds on October 1, 2015. The simulated C chla ranged from 10 to 25 mg m⁻³ and gradually increased from the inlet to the outlet of each pond. The curve for each pond has a downward trend at the weir locations. Because a weir can alter the water level and discharge to reduce nutrient loading, C chla was affected by these structures. In each pond, Chl-a should increase with an increase in the distance to the water inlet, because sewage impact reduces. The results in Figure 8 show that the general trends of Chl-a in all three ponds confirm this trend. However, two downward trends (abnormal concaves) can be observed at distances of approximately 500 m in ponds 1 and 2. These locations coincided with the two weirs in ponds 1 and 2. Therefore, we speculate that these two abnormal concaves are caused by the weirs. At the same time, since pond 3 does not have a weir, no downward trend can be found in the pond 3 data.

Figure 8

Variation in DUFLOW-modeled Chl-a concentration with distance from inlet in each pond of Kristalbad.

3.6 C _chla prediction

Figure 8a presents the 1-month predictions of C chla dynamics in the three ponds, from October 1 to November 1, 2015. Overall, C chla was predicted to be 3–30 mg m⁻³ in Kristalbad, with a higher value predicted for pond 3 (3–40 mg m⁻³) than for the other two ponds. Figure 9a shows that the predicted trends of C chla were very similar in each pond and displayed the following features: (1) a positive trend from the October 1–8, (2) the peak of C chla on October 15 (>30 mg m⁻³), and (3) an obvious decrease from October 16 to the end of the period.

Figure 9

DUFLOW-modeled Chl-a concentration in Kristalbad from October 1 to November 1, 2015. (a) Comparison of the three ponds, (b) relationship with light energy for pond 1, (c) relationship with water temperature for pond 1.

Phytoplankton activities are strongly influenced by the intensity of photosynthesis which is sensitive to weather. Therefore, their changes may be attributed to the influences of light energy and temperature. Figure 9b and c shows the relationship between predicted C chla and light energy and temperature, respectively, from October 1 to November 1, 2015. During the prediction period, the light energy was 0–650 W m⁻² and the temperature was 1–18°C. It is clear that the light energy and water temperature directly affected C chla , but around October 14, 2015, temperature had the opposite effect. The results in Figure 9 show that light energy is a positive factor while temperature plays a negative role in phytoplankton growth in Kristalbad.

To explore the behavior of phytoplankton in Kristalbad, we used water quality mapping and modeling, integrating SPOT 6 high-spatial resolution image, field measurements of surface reflectance, laboratory measurements of C _chla, and the DUFLOW water quality model. For retrieval of Chl-a in the Kristalbad wetland, an empirical algorithm based on OC2 of SPOT 6 and field data was developed. In general, in the optically shallow inland water, we considered the two-band ratio model based on the red and NIR bands at 670 and 710 nm, respectively. This is because these two bands are primarily sensitive to chlorophyll pigments. For such small and shallow inland water, only a high spatial resolution image (±5 m) was appropriate for Kristalbad. Therefore, SPOT 6 satellite data set with only four bands was used to map the spatial distribution of C chla Once the time series data of environmental conditions were established in DUFLOW, the predicted C chla dynamics were modeled for each pond in Kristalbad from October 1 to November 1, 2015. To justify continued monitoring of plankton changes, a connection must be made between the variation in plankton and removal of nutrients or other chemicals. We aim to make this connection in the future studies.