Yumeng Song and Jing Zhang

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

De Gruyter | Published online: April 22, 2021


We integrated hyperspectral and field-measured chlorophyll-a (Chl-a) data from the Kristalbad constructed wetland in the Netherlands. We developed a best-fit band ratio empirical algorithm to generate a distribution map of Chl-a concentration (Cchla) from SPOT 6 imagery. The Cchla retrieved from remote sensing was compared with a water quality model established for a wetland pond system. The retrieved satellite results were combined with a water quality model to simulate and predict the changes in phytoplankton levels. The regression model provides good retrievals for Chl-a. The imagery-derived Cchla performed well in calibrating the simulation results. For each pond, the modeled Cchla showed a range of values similar to the Chl-a data derived from SPOT 6 imagery (10–25 mg m−3). The imagery-derived and prediction model results could be used as the guiding analytical tools to provide information covering an entire study area and to inform policies.

1 Introduction

The increase in the industries, population, and human activities lead to increasing demands for using water. This will cause a destruction of water ecosystem and an augmentation of water pollutions. And these are the main sources of nutrient supplements in water environment. When nutrients exceed the limitation in a water ecosystem, it will lead to the eutrophication, heave metal contamination, and damage to the public health [1,2]. All kinds of heavy metals entering into the water are easily absorbed, complexed, or coprecipitated by suspended solids or sediments, resulting in sediment deposition or lake bed elevation in the water environment, posing a long-term threat to the aquatic ecosystem [3]. Many ways are involved in removing the nutrients in the aquatic system. One deliberated approach used for effluent (contains high nutrients) purification is utilizing the wetlands (both natural and artificial). Wetlands as biological treatment systems, which contain various aquatic plants and phytoplanktons, are used for treating the water and purifying polluted water from high nutrients [4]. At the same time, this method contains three complex processes, namely, biological, physical, and chemical processes. Therefore, the wetland aquatic vegetation plays an important role in purifying the urban effluents. The initial function of removing the nutrient loads is achieved through filtering and setting of organic and inorganic particles that are associated with the nutrients to make them slowly pass through the wetland. The aquatic vegetation also has a substrate for attaching the decomposed microorganisms, which behave like a filter to block the dissolved organic matters. Furthermore, aquatic vegetation influences the nitrification and denitrification by controlling the dissolved oxygen concentration of the wetland within the rhizosphere. And they can also provide bacteria to fix the N in root nodules. All the abovementioned give a good evidence to prove a great capability of aquatic vegetation in purifying wastewater [5,6].

The use of constructed or artificial wetlands is an efficient method for treating polluted water containing substantial nutrients, and this process has been well described [7,8,9,10,11]. In the Netherlands, the constructed wetlands are used more and more, often following the water harmonica concept [12]. An effective case of artificial wetland use is Kristalbad, which was designed to provide multiple water functions and ecosystem services. However, over the years, increasing human activities, domestic wastewater, and industrial policies have drastically increased urban sewage emissions and the resulting nutrient loading of Kristalbad, affecting its ecological health and function. Phytoplanktons are an important component of wetland microbes for purifying wastewater. The phytoplankton community plays an important and complex role in removing nutrients, organic matter, and other pollutants [13]. Dynamic changes in the ecofunction of Kristalbad can be indicated by phytoplankton population data. [14]. The distribution of phytoplankton can be mapped by characterizing the spatial distribution of chlorophyll-a (Chl-a), since chlorophyll is the main pigment in the phytoplankton cells. Therefore, these dynamic changes in the ecofunction of Kristalbad, as indicated by Chl-a concentration (Cchla), can be monitored using phytoplankton population data.

For continuous and intensive management of the water quality in wetland reserves, integrating in situ sampling, modeling, and remote sensing could help to lower the costs and time spent on measurement and provide complementary results. Moreover, using water quality modeling, we can forecast Cchla dynamics for a certain period [15,16].

Using remote sensing to derive Cchla is a rapid and low-cost method of phytoplankton monitoring and has thus been studied for large-scale applications such as oceans [17,18,19] and coastal areas [20,21,22]. However, because the constructed wetlands have a relatively small size, high-resolution satellite images are needed to extract Cchla for phytoplankton monitoring. Geospatial and high-resolution data from satellites such as FORMOSAT-2, SPOT 6, and LANDSAT-8 can be used for remote monitoring of water systems, in addition to traditional analyses using the simplified hydrologic model.

The main objective of this research was to develop a low-cost and efficient approach to monitoring and to model the spatial distribution of phytoplankton and predict changes to their concentration in constructed wetlands such as Kristalbad. Wetlands, as biological treatment systems that contain various phytoplanktons, are used for treating water and purifying water polluted by heavy nutrient loading. The method contains complex biological, physical, and chemical processes. Therefore, we used remote sensing and other spectral-measured data to retrieve the distribution of Cchla. We combined the satellite retrieval results with the DUFLOW water quality model to simulate and predict the dynamic changes in phytoplankton in a given period. In this article, we aim to (1) determine the relationship between the Cchla and spectral properties of artificial wetlands, using field spectrometry and water sampling, (2) map the spatial distribution of phytoplankton by integrating field spectrometry, SPOT 6 images, and water sampling, (3) build a prototype ecohydraulic model to simulate the hydrologic behavior of artificial wetlands, and (4) build a water quality prediction model to simulate and predict the changes in Cchla during a certain period.

2 Materials and methods

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.

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 Cchla, by combining the satellite data and spectral measurements obtained with a Trios/Ramses spectrometer.

Figure 2 Location of water sampling sites in Kristalbad.

Figure 2

Location of water sampling sites in Kristalbad.

2.2.2 Laboratory analysis of C chla

The water samples were analyzed in the University of Twente laboratory, on the day after field sampling. Water samples were kept in the dark in a refrigerator, at −20°C, filtered with a Whatman GF/C filter, and analyzed for Cchla according to method 446.0, i.e., in vitro determination of chlorophyll by visible spectrophotometry [23]. Chl-a pigment was extracted using 90% acetone at 4°C overnight in the dark. Filtered water sample of 700 mL and extraction solution of 100 mL were required, and 4 mL water sample was absorbed by the filter paper.

We measured the sample absorbance at each wavelength using a spectrophotometer (visible, multi-wavelength, resolution <2 nm). We warmed the spectrophotometer and used 90% acetone to zero the instrument at each selected wavelength. To determine Chl-a, the absorbance at 750, 664, 647, and 630 nm was used in the Jeffrey and Humphrey trichromatic equations. That is, we subtracted the absorbance value at 750 nm from those at 664, 647, and 630 nm, using

(1) C E , a = 11.85 ( Abs  664 ) 1.54 ( Abs  647 ) 0.08 ( Abs  630 )
where C E,a is the concentration of Chl-a (mg/L) in the extraction solution.

We calculated the total concentration of Chl-a in each water sample using the general equation:

(2) C s,a = C E,a × Extract volume (L) × DF Sample volume ( L ) × Cell length ( cm )
where C s,a is the concentration of Chl-a (mg/L) in the water sample; extract volume is the volume of the extract solution (L; in this case 0.104); DF is any dilution factor; sample volume is the volume of the entire filtered water sample (L; in this case 0.7); and cell length is the optical path length (cm) of the cuvette used in the analysis (typically 1 cm).

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

(3) L w ( λ ) = L u ( λ ) ρ × L sky ( λ )

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

(4) R rs ( λ ) = L w ( λ ) E d ( λ )

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.

Figure 3

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

2.4.2 Data processing

A SPOT 6 image of the study area was orthorectified using a digital elevation model, its rational polynomial coefficient, and ERDAS IMAGINE (an image processing and GIS software package). Orthorectification is the process of removing distortion from an image. SPOT 6 image was also used to provide detailed background information of the Kristalbad wetland, which helped to build the DUFLOW Modeling Studio. Next the Kristalbad water body was identified through unsupervised classification. Through ERDAS IMAGINE, the NIR band (band 4) of SPOT 6 was classified into two categories: water areas and nonwater areas. Using a recode module to renumber the classes, water areas were set to “1” while other classes were set to “0.” We then extracted the water areas using a mask module in the ERDAS software. This step was performed to select water areas on the image and minimize the effect of land, vegetation, or artificial structures on image analysis. The Kristalbad water area is shown in Figure 3b.

For quantitative remote sensing, an appropriate atmospheric correction must be applied to convert the digital number (DN) value of image pixels into a value representing remote sensing reflectance from radiance. For high-resolution imaging data, a number of atmospheric correction algorithms have been developed, such as atmosphere correction, imaging spectrometer data analysis system, and high-accuracy atmospheric correction for hyperspectral data. In our case, fast line-of-sight atmospheric analysis of spectral hypercubes (FLAASH) was applied to SPOT 6 atmospheric correction. FLAASH is based on an atmospheric radiative transfer model, namely, moderate resolution atmospheric transmission (MODTRAN 5). Atmospheric correction was executed using ENVI software. The input parameters used for atmospheric correction are listed in Table 1.

Table 1

Input parameter values used for atmospheric correction

Parameters Value
Scene center location 52.25
Sensor type SPOT 6
Sensor altitude 695 km
Ground elevation 0.011 km
Pixel size 6 m
Flight date October 1, 2015
Flight time GMT 10:00:00
Atmospheric model Subarctic summer
Aerosol model Rural
Aerosol retrieval 2-Band (K-T)
Initial visibility 25

In this study, we used an unsupervised classification method to separate water from land and evaluated the classification results. Unsupervised classification involves classification without prior (known) classification criteria. It is based on the reflectance feature vector of image features in the feature space. Therefore, when ERDAS carries out unsupervised classification, only the classification number, the maximum iteration number, and the transformation threshold need to be set. The classification is defined mainly through field investigation and visual interpretation. Through field investigation and image comparison, the extracted wetland water area was found to be consistent with the actual situation.

For retrieval of C chla from the satellite image, the first step was to match the sample locations with specific pixels on the image. Each pixel was matched with one corresponding sample site, without applying the window method to average the pixel (e.g., 3 × 3, 5 × 5). Because Kristalbad is very small and many things can cause interference with the sampling stations, such as weirs, bank vegetation, or station shadows, the pixel value can be affected when an automatically averaged n × n window is applied to the image. The total number of sample sites and their linked pixel reflectance values are shown in Table 2. Various parameters are needed to model the water of the Kristalbad wetland, including water flow parameters (such as discharge and temperature), water components (concentrations of NH4, PO4, and BOD), and other external variables (water surface intensity). All these parameters were acquired from the local water authority and other research.

Table 2

Sampled pixel reflectance values

Sample code Pixel reflectance values
Band 1 Band 2 Band 3 Band 4
SK1 0.068 0.096 0.045 0.101
SK2 0.065 0.093 0.038 0.133
SK3 0.063 0.095 0.04 0.044
SK4 0.067 0.094 0.053 0.087
SK5 0.074 0.1 0.056 0.08
SK6 0.068 0.092 0.044 0.115
SK7 0.064 0.133 0.056 0.08

2.5 Regression modeling

2.5.1 Resampling of hyperspectral data

This section mainly describes the algorithmic method used to retrieve the spatial distribution of C chla from SPOT 6 satellite image, by integrating R rs ( λ ) with in situ Chl-a. To develop a reliable regression model, the following two data sets are needed: one for regression modeling and the other for model validation. The first data set was used to derive the relationship between the in situ C chla and R rs ( λ ) obtained with the radiometer. Then we applied this mathematical relationship to SPOT 6 image to derive the spatial distribution of C chla in Kristalbad. This method effectively correlated the field spectra reflectance with remote sensing reflectance obtained from the atmospheric corrected satellite image.

The R rs ( λ ) hyperspectral data measured with the TriOS Ramses radiometer were recorded in 172 bands within the visible range from 380 to 950 nm, with variable spectral resolution from 2.6 to 3.5 nm. All these data were resampled to four SPOT 6 spectral bands R rs ( i ) . The resampling method was specific to the transmittance functions of the remote sensing imager (RSI) onboard SPOT 6 [27]. The following equation is the resampling function [28]:

(5) R rs ( i ) = λ = 400 λ = 1 , 000 R rs ( λ ) T i ( λ ) λ = 400 λ = 1 , 000 T i ( λ )
where R rs ( i ) is regarded as the remote sensing reflectance obtained by the RSI, i represents the blue, green, red, and NIR bands with values 1, 2, 3, and 4, respectively. T i ( λ ) was set from 400 to 1,000 nm at 1 nm intervals based on prelaunch data.

Table 3 shows the results of resampling the TriOS Ramses R rs ( λ ) into R rs ( i ) for each sample station. In general, the resampled hyperspectral data retained most of the spectral features linked to the optical properties of various water types. In this instance, the ratios of the green to blue bands and NIR to red bands are both proportional to C chla .

Table 3

Remote sensing reflectance R rs ( i ) for each sample station

Sample code R rs ( 1 ) R rs ( 2 ) R rs ( 3 ) R rs ( 4 )
SK1 0.0023 0.0027 0.0013 0.0005
SK2 0.0016 0.0028 0.0015 0.0009
SK3 0.0012 0.0021 0.0009 0.0007
SK4 0.0024 0.0038 0.0021 0.0010
SK5 0.0017 0.0018 0.0008 0.0009
SK6 0.0020 0.0026 0.0012 0.0008
SK7 0.0114 0.0233 0.0342 0.0259

2.5.2 Water leaving reflectance and Cchla

In nature, the R rs ( λ ) of a lake is affected by all components of the water body. Because Kristalbad is a very shallow wetland, the bottom may contribute to R rs ( λ ) . It is therefore necessary to determine which part of the spectral signal can be attributed to Chl-a alone. In retrieving C chla from spectral reflectance, the optical properties of Chl-a can provide spectral information. The absorption of blue and red light could indicate the presence of the Chl-a pigment. Further, the maximum reflectance in the green- (560–570 nm) and red-edge (near 710 nm) bands can also be considered. We developed an empirical band ratio regression model to retrieve C chla , because the ratios can highlight latent information when there is an inverse relationship between spectral responses to the same biophysical phenomenon [29].

Multiband ratio algorithms have been proposed and used for estimating C chla because they produce the largest satellite sensor signal-to-noise ratio (SNR) for different C chla values [17]. For example, the ocean chlorophyll 2-band (OC2) algorithm, which uses the ratio between blue and green reflectance, is commonly used to retrieve phytoplankton concentrations between 0 and 10 mg m−3 [30]. Spectral bands with longer wavelengths, that is, the NIR, red bands, and the associated band ratios, are often used to investigate C chla in a turbid environment. A three-band ratio between blue (482 nm), red (660 nm), and NIR (825 nm) can be used to retrieve Chl-a in very shallow and turbid lakes. Our specific objective was to test the feasibility of multiband ratio algorithms and then find a best-fit model for Chl-a retrieval. All models tested in this research were from Singh et al. [31]. We also tested the performance of the resampled R rs ( λ ) in a regression analysis of Chl-a. These models can be calculated using the following equations:

(6) Two-band model : C chla R rs 1 ( λ 1 ) × R rs ( λ 2 )
(7) Two-band model: C chla R rs 1 ( λ 3 ) × R rs ( λ 4 )
(8) Three-band model: C chla [ R rs 1 ( λ 3 ) [ R rs 1 ( λ 1 ) ] × R rs ( λ 4 )

Here λ 1 = 485 nm, λ 2 = 560 nm, λ 3 = 660 nm, and λ 4 = 825 nm. R rs ( λ 1 ) ,   R rs ( λ 2 ) , R rs ( λ 3 ) , and R rs ( λ 4 ) are the blue, green, red, and NIR ranges, respectively, for resampled R rs ( λ ) and SPOT 6 image bands 1, 2, 3, and 4, respectively. Therefore, the relationship between the measured C chla and each spectral ratio was calculated by all possible linear and nonlinear functional forms in the regression analysis, such as exponential, polynomial, logarithmic, and power. The performances of the equations used to retrieve Chl-a from SPOT 6 image were assessed using the determination coefficient (R2), root mean square error (RMSE), and relative percentage difference (PRD). The PRD was used to quantify the absolute difference between imagery-derived values and the corresponding observations.

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 Cchla 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 Results and discussion

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−3, with a mean of 10 mg m−3. 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−3 Chl-a. The temperature, suspended matter, turbidity, and measured Cchla 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, Cchla, 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 Cchla was 10 mg m−3 and the maximum Cchla was 50.25 mg m−3. Therefore, using the OECD classification, the Kristalbad wetland can be defined as eutrophic.

Table 4

Laboratory analysis results of water quality for each sample station

Sample code Geographic location Time Temp (°C) Turbidity (NTU) Laboratory-measured Cchla (mg m−3) Laboratory-measured SPM (mg/L)
(m) (m)
SK1 351681.75 5790262.84 11:21:03 11 1.94 14.80 6
SK2 351541.69 5790540.77 11:50:55 11.1 2.26 18.80 12
SK3 351524.65 5790553.74 12:04:29 11.3 2.272 25.70 22
SK4 351467.87 5790614.26 12:23:00 11.4 2.686 14.38 14
SK5 351447.68 5790639.58 12:36:31 11.4 1.7346 11.42 16
SK6 351225.6 5790600.06 13:00:54 12 1.632 12.05 22
SK7 351141.52 5790564.23 13:23:00 13.5 5.365 50.25 No data

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:

(9) Two-band model: R rs 1 ( 485 ) × R rs ( 560 )
(10) Chl-a = 29.609 ( x ) 24.898

Table 5

Results of two-band ratio regression modeling showing the functional relationship between Rrs and Cchla

Band ratio Model equation R2
[ R rs 1 ( 660 ) R rs 1 ( 485 ) ] × R rs 1 ( 825 ) y = 2 E 05 ( x ) + 29.631 0.304
y = 25.851 e 1 E 06 x 0.334
y = 5 E 11 ( x 2 ) 6 E 05 x + 32.372 0.353
R r s ( 560 ) / R r s ( 485 ) y = 29.609 ( x ) 24.898 0.79
y = 8.9138 X 1.7681 0.70
y = 43.604 ln ( x ) + 3.1994 0.70
y = 2.961 e 1.1766 x 0.77
y = 40.3192 ( x 2 ) 102.14 x + 76.546 0.78
R rs ( 825 ) × R rs 1 ( 660 ) y = 3.1099 ( x ) + 18.865 0.004
y = 19.609 X 0.1587 0.012
y = 6.2916 ln ( x ) + 23.606 0.03
y = 18.494 e 0.008 x 0.002
y = 125.26 ( x 2 ) 209.87 x 57.07 0.405
[ R rs ( 485 ) R rs ( 660 ) ] / R rs ( 825 ) y = 11.959 ( x ) + 27.986 0.606
y = 23.838 e 0.448 x 0.59
y = 7.9457 ( x 2 ) 20.974 x + 25.072 0.77
[ R rs ( 485 ) R rs ( 660 ) ] / [ R rs ( 485 ) + R rs ( 660 ) ] y = 44.656 ( x ) + 25.365 0.86
y = 21.543 e 1.641 x 0.815
y = 34.922 ( x 2 ) 38.531 x + 21.979 0.86
[ R rs 1 ( 660 ) R rs 1 ( 565 ) ] × R rs ( 485 ) y = 26.487 ( x ) + 39.023 0.72
y = 36.04 e 0.992 x 0.5
y = 23.601 ( x 2 ) 50.704 X + 40.614 0.85
Figure 4 Spectral reflectance of resampled Rrs(λ){R}_{\text{rs}}(\lambda ) for four SPOT 6 bands, with varying phytoplankton Chl-a concentrations, from sample station 3.

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 R2 = 0.79, RMSE = 5.80 mg m−3 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−3) than the center (<5 mg m−3). 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.

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 Cchla with in situ measured Cchla (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−3 while the RPD was 2.25. In Figure 6b, a strong correlation was observed between the imagery-derived and laboratory-measured data (R2 = 0.8). However, differences were still noticed because some uncertain factors were included.

Figure 6 Comparison between Kristalbad Cchla measured in laboratory and that derived from SPOT 6 imagery. (a) Comparison in situ, (b) the correlation of them, with R2.

Figure 6

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

Table 6

Estimation error between laboratory-measured and SPOT 6 imagery-derived Chl-a concentration

Sample code Laboratory-measured Chl-a (mg m−3) Imagery-derived Chl-a (mg m−3) Δ (difference) RMSE (mg m−3) RPD
SK1 14.80 16.90 2.105978811 6.13 2.25
SK2 18.80 17.475 −1.328168589
SK3 25.70 19.75 −5.952915388
SK4 14.38 16.64 2.259301772
SK5 11.42 15.16 3.74415318
SK6 12.05 16.13 4.080985714
SK7 50.25 36.63 −13.62072542

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 Cchla 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−3 and then remained steady at approximately 15 mg−3 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−3.

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

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−3 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.

Figure 8

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

3.6 Cchla 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−3 in Kristalbad, with a higher value predicted for pond 3 (3–40 mg m−3) 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−3), 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.

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−2 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 Cchla, 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.

4 Conclusion

The principal conclusions are as follows:

  1. (1)

    In the study area, Cchla values ranging from 11.41 to 50.25 mg m−3 were observed and used to develop the regression model. Using the OECD classification system, the entire wetland can be defined as eutrophic since the average Cchla is above 8 mg m−3, and the maximum value is above 25 mg m−3.

  2. (2)

    Results from regression modeling indicate that the OC2 model, adopting blue and green spectral channels, is the best-fit model for Cchla retrieval. The regression model provides good retrievals for Chl-a, with an RMSE of 5.8 mg m−3. After applying the best-fit regression model to SPOT 6 image, the Cchla can be derived from the image with an RMSE of 6.13 mg m−3. The estimation accuracy could be enhanced by improving the accuracy of the atmospheric correction and surface reflectance data. The FLAASH atmospheric correction method requires a large amount of calculation and more parameters. For example, the water vapor content in the atmosphere, ozone content, spatial distribution, aerosol optical characteristics, etc. In conventional atmospheric calibration, such measurements are difficult to implement. The difficulty of atmospheric correction lies in determining these parameters. How the parameters are determined directly affects the calculation accuracy. Quick atmospheric correction (QUAC) is a fast atmospheric correction algorithm that has been extensively verified and evaluated. It is based on the empirical finding that the average reflectance of a collection of diverse material spectra, such as the end-member spectra in a scene, is essentially scene independent. In the future, we will try to use this QUAC method to improve the estimation accuracy.

  3. (3)

    Each sampling site showed a good agreement between the measured data and the map of modeled Chl-a distribution. Near the bank of the ponds, the level of Cchla was higher than that in the center, which was caused by the adjacency effect.

  4. (4)

    The imagery-derived Chl-a data, combined with the hydraulic data, were applied to the DUFLOW water quality modeling (EUTROF 1). The imagery-derived Cchla data were used to calibrate the simulation results.

  5. (5)

    The modeled Cchla for each pond in Kristalbad showed a similar range of 10–25 mg m−3 as the Chl-a data derived from SPOT 6 imagery, and increased gradually from the inlet to the outlet of each pond. The prediction of Cchla dynamics for each pond was also modeled, with Cchla ranging from 5 to 35 mg m−3.

Overall, C chla estimation integrating high-resolution remote sensing images and water surface reflectance is feasible for rapid monitoring and modeling of C chla , but the estimation accuracy can still be improved. The impact of SPM, CDOM, or nutrients on Cchla estimation in wetland waters should be further quantified. The performance of the best-fit algorithmic model could be improved by adopting an advanced OC4 model and a variation in the multiband model that uses multispectral sensors, such as HYPERION. A simple eutrophication model was used to simulate Chl-a and algae concentrations. A more complex ecological water quality model would be used to present the real wetland eutrophication conditions. The accuracy of DUFLOW water quality modeling could be improved by considering more environmental parameters, such as water discharge, precipitation, and water storage. For future research, it should be noted that the monitoring and modeling were only validated in Kristalbad for one season (autumn). To reduce the uncertainty of the research, it is imperative to consider different environmental regimes and seasons.


This work was supported by the National Natural Science Foundation of China (41271004). We thank Elsevier Language Editing Services for editing the English text of a draft of this manuscript.

    Author contributions: Y. S. carried out the experiments. Y. S. and J. Z. performed the data analysis and wrote the draft manuscript. Y. S and J. Z. designed the experiments. All authors contributed to the discussions and reviewed the manuscript.

    Conflict of interest: The authors declare no conflict of interest.


[1] Ekoa Bessa AZ, Ngueutchoua G, Kwewouo Janpou A, El-Amier YA, Nguetnga OA, Kayou UR, et al. Heavy metal contamination and its ecological risks in the beach sediments along the Atlantic Ocean (Limbe coastal fringes, Cameroon). Earth Syst Environ. 2020. 10.1007/s41748-020-00167-5. Search in Google Scholar

[2] Bhardwaj LK, Jindal T. Persistent organic pollutants in lakes of Grovnes Peninsula at Larsemann Hill area, East Antarctica. Earth Syst Environ. 2020;4(2):349–58. 10.1007/s41748-020-00154-w. Search in Google Scholar

[3] Kumar R, Bahuguna IM, Ali SN, Singh R. Lake inventory and evolution of Glacial lakes in the Nubra-Shyok Basin of Karakoram Range. Earth Syst Environ. 2020;4(1):57–70. 10.1007/s41748-019-00129-6. Search in Google Scholar

[4] Sharip Z, Abd Razak SB, Noordin N, Yusoff FM. Application of an effective microorganism product as a cyanobacterial control and water quality improvement measure in Putrajaya Lake, Malaysia. Earth Syst Environ. 2020;4(1):213–23. 10.1007/s41748-019-00139-4. Search in Google Scholar

[5] Bhardwaj LK, Jindal T. Persistent organic pollutants in lakes of Grovnes Peninsula at Larsemann Hill area, East Antarctica. Earth Syst Environ. 2020;4(2):349–58. 10.1007/s41748-020-00154-w Search in Google Scholar

[6] Sharip Z, Razak SB, Noordin N, Yusoff FM. Application of an effective microorganism product as a cyanobacterial control and water quality improvement measure in Putrajaya Lake, Malaysia. Earth Syst Environ. 2020;4(1):213–23. 10.1007. Search in Google Scholar

[7] Chan E. The use of wetlands for water pollution control; 1982. Search in Google Scholar

[8] Gersberg RM, Elkins BV, Goldman CR. Nitrogen removal in artificial wetlands. Pergamon. 1983;17(9):1009–14. Search in Google Scholar

[9] Coleman J, Hench K, Garbutt K, Sexstone A, Bissonnette G, Skousen J. Treatment of domestic wastewater by three plant species in constructed wetlands. Water, Air, Soil Pollut. 2001;128(3–4):283–95. Search in Google Scholar

[10] Vymazal J. Removal of nutrients in various types of constructed wetlands. Sci Total Environ. 2007;380(1–3):48–65. Search in Google Scholar

[11] Zhang S, Song H-L, Yang X-L, Huang S, Dai Z-Q, Li H, et al. Dynamics of antibiotic resistance genes in microbial fuel cell-coupled constructed wetlands treating antibiotic-polluted water. Chemosphere. 2017;178:548–55. Search in Google Scholar

[12] Kampf R, van den Boomen RM. Waterharmonica’s in the Netherlands (1996–2012), natural constructed wetlands between well treated waste water and usable surface water. Foundation for Applied Water Research, The Netherlands; 2013. 10.13140/2.1.3930.1128. Search in Google Scholar

[13] Iamchaturapatr J, Yi SW, Rhee JS. Nutrient removals by 21 aquatic plants for vertical free surface-flow (VFS) constructed wetland. Ecol Eng. 2007;29(3):287–93. Search in Google Scholar

[14] Ahmed A, Wanganeo A. Phytoplankton succession in a tropical freshwater lake, Bhoj Wetland (Bhopal, India): spatial and temporal perspective. 2015;187(4):192. Search in Google Scholar

[15] Huntley M, Brooks ER. Effects of age and food availability on diel vertical migration of Calanus pacificus. Mar Biol. 1982;71(1):23–31. Search in Google Scholar

[16] Dekker AG, Hoogenboom HJ, Goddijn LM, Malthus TJM. The relation between inherent optical properties and reflectance spectra in turbid inland waters. Remote Sens Rev. 1997;15(1–4):59–74. Search in Google Scholar

[17] O’Reilly J. Ocean color chlorophyll a algorithms for SeaWiFS, OC2, and OC4: Version 4. SeaWiFS Postlaunch Calibration Valid Analyses. 2000;11:9–23. Search in Google Scholar

[18] Subramaniam A, Brown CW, Hood RR, Carpenter EJ, Capone DG. Detecting Trichodesmium blooms in SeaWiFS imagery. Deep-Sea Res Part II. 2001;49:1. Search in Google Scholar

[19] Melin F. Global distribution of the random uncertainty associated with satellite-derived Chla. IEEE Geosci Remote Sens Lett. 2010;7(1):220–4. Search in Google Scholar

[20] Gitelson AA, Dall’Olmo G, Moses W, Rundquist DC, Barrow T, Fisher TR, et al. A simple semi-analytical model for remote estimation of chlorophyll- a in turbid waters: validation. Remote Sens Environ. 2008;112(9):3582–93. Search in Google Scholar

[21] Moses WJ, Gitelson AA, Berdnikov S, Povazhnyy V. Estimation of chlorophyll-a concentration in case II waters using MODIS and MERIS data-successes and challenges. Environ Res Lett. 2009;4(4):045005. Search in Google Scholar

[22] Al-Naimi N, Raitsos DE, Ben-Hamadou R, Soliman Y. Evaluation of satellite retrievals of chlorophyll-a in the Arabian Gulf. Remote Sens. 2017;9(3):301. Search in Google Scholar

[23] Arar EJ. In vitro determination of chlorophylls a,b, C1 + C2 and pheopigments in marine and freshwater algae by visible spectrophotometry. EPA Method. 1997;446:1–22 (Print). Search in Google Scholar

[24] Lee ZP, Carder KL, Steward RG, Peacock TG, Davis CO, Patch JS. An empirical algorithm for light absorption by ocean water based on color. J Geophys Research-Oceans. 1998;103(C12):27967–78. Search in Google Scholar

[25] Mobley CD. Estimation of the remote-sensing reflectance from above-surface measurements. Appl Opt. 1999;38(36):7442–55. Search in Google Scholar

[26] Dominique C, Aline B, Emmanuel K, Rachid H, Olivier H, Olivier M, et al. Assessing the potentialities of FORMOSAT-2 data for water and crop monitoring at small regional scale in South-Eastern France. Sens (Basel, Switz). 2008;8(5):3460–81. Search in Google Scholar

[27] Ma SB, Yang WF, Zhang K. Study of key technology of SPOT6 satellite image processing. Remote Sens LResour. 2015;27(3):30–5. Search in Google Scholar

[28] Chih-Hua C, Cheng-Chien L, Ching-Gung W, I-Fan C, Chi-Kin T, Ching-Shiang H. Monitoring reservoir water quality with Formosat-2 high spatiotemporal imagery. J Environ Monitor: JEM. 2009;11(11):1982–92. Search in Google Scholar

[29] Campbell JB. Introduction to Remote Sensing. Boca Raton: CRC Press; 2002. Search in Google Scholar

[30] Odermatt D, Gitelson A, Brando VE, Schaepman M. Review of constituent retrieval in optically deep and complex waters from satellite imagery. Remote Sens Environ. 2012;118:116–26. Search in Google Scholar

[31] Singh K, Ghosh M, Sharma SR, Kumar P. Blue-red-NIR model for chlorophyll-a retrieval in hypersaline-alkaline water using landsat ETM plus Sensor. IEEE J Sel Top Appl Earth Obser Remote Sens. 2014;7(8):3553–9. Search in Google Scholar

[32] Vieira JMP, Pinho JLS, Duarte AALS. Eutrophication vulnerability analysis: a case study. Water Sci Technol. 1998;37(3):121–8. Search in Google Scholar

[33] Carder KL, Cannizzaro JP, Lee Z. Ocean color algorithms in optically shallow waters: limitations and improvements. SPIE Opt + Photo. 2005. Search in Google Scholar

[34] Hu L, Liu Z, Liu Z, Hu C, He MX. Mapping bottom depth and albedo in coastal waters of the South China Sea islands and reefs using Landsat TM and ETM + data. Int J Remote Sens. 2014;35(11–12):4156–72. Search in Google Scholar

Received: 2020-06-24
Revised: 2021-02-14
Accepted: 2021-03-11
Published Online: 2021-04-22

© 2021 Yumeng Song and Jing Zhang, published by De Gruyter

This work is licensed under the Creative Commons Attribution 4.0 International License.