2.1.

Study Area

The study area is Taihu Lake, the second largest freshwater lake in China, located between 30°56′ to 31°33′ N and 119°55.3′ to 120°53.6′ E and having an area of $2427.8{\mathrm{km}}^2$ and an average depth of 2.12 m. This lake is large and the water in the lake is highly turbid. The eutrophication of the lake is serious, with the hypereutrophic area mainly covering Meiliang Bay, Wuli Lake, and the Western Lake.²⁵ The Chla in the lake has obvious seasonal and annual variation, with the highest concentration in July to August when algae blooms occur.

Four datasets were used in this study, including July to August of 2004 and 2005, as well as September of 2010 and 2011, and the sampling distributions are shown in Fig. 1. In July to August of 2004 and 2005, water samples were collected at the Taihu Lake monitoring sites and the spectra were measured on the first 10 days of each month. On September 19, 2010, and September 3, 2011, the sampling positions mainly covered the hypereutrophic area of Meiliang Bay and the central lake, with spectral measurements taken near noon within 1 day.

Fig. 1

Sampling distribution in July to August of 2004 and 2005 (a), September of 2010 (b) and September of 2011 (c) in Taihu Lake, China.

2.2.

Data Acquisition

Using an Analytical Spectral Devices Field Spectroradiometer (Analytical Spectral Devices Inc., Boulder, Colorado) with 512 bands ranging from 350 to 1050 nm with increments of 1.5 nm, the hyperspectral reflectance in the study area was measured at 0.5 to 1 m high above the water surface, with a probe field angle of 10 deg. The instrument was positioned at a specific viewing geometry to avoid the effects of direct solar radiation and prevent the ship from interfering with the water surface.²⁶ Ten curves were acquired for each location, and the median value of the repeated measurements was used to calculate the remote sensing reflectance. The reflectance of a standard gray plate is 30%, and the spectra with wavelengths shorter than 400 nm or longer than 900 nm were discarded owing to noise. The spectra were resampled to 1-nm interval and smoothed using the kernel regression smoothing method with a window width of 5 nm.²⁷

Discarding samples under windy or cloudy conditions, a total of 85 samples were finally kept, of which 24 samples in July to August of 2004 were used for model building and 20 samples in July to August of 2005, 25 samples in September of 2010, and 16 samples in September of 2011 were used for model validation.

Chla was measured according to the Chinese national standard three-color spectrophotometry (SL88-1994). First, the water samples collected in the field were filtered through a Whatman GF/C membrane, after which the membrane was kept in darkness in a refrigerator overnight. After the removal from refrigeration, Chla was extracted using 90% acetone. The extracted liquid was centrifuged for 10 min, and the supernatant was spectrophotometrically analyzed using a Shimadzu UV-2550 UV-Vis spectrophotometer. The Chla ($\mathrm{mg}/{\mathrm{m}}^3$) was calculated using the absorbance at 750, 663, 645, and 630 nm, according to the standard formula.

The total suspended sediment (TSS) concentration ($\mathrm{mg}/\mathrm{L}$) was determined gravimetrically according to the Chinese national standard (GB11901-89, 1990).

2.3.

Spectral Index

The spectral index is the mathematical combination of reflectance at the visible and NIR bands,²⁸ including vegetation index, fluorescence index, and trilateral index. A total of 106 spectral indices were collected in this study (Tables 1 and 2).

Table 1

List of the spectral vegetation indices.a

Table 2

List of spectral fluorescence and triangular indices.a

2.4.

Model Building and Accuracy Assessment

Regression analysis was used to build the model between Chla and the spectral index. Both the three-band^9,10,15 and the four-band¹² algorithms were used for model comparison. The formulas of the three- and four-band models are as follows:

Based on the study by Zimba and Gitelson,⁷⁵ the three bands were searched in the range of 450 to 800 nm and the initial iterative positions of $\lambda 1$ and $\lambda 3$ were 675 and 750 nm, respectively. Iterative calculation will stop when the root mean square error (RMSE) comes to its lowest. The optimal four bands were searched using the same iteration method, and the initial positions of $\lambda 2$, $\lambda 3$, and $\lambda 4$ were 700, 720, and 750 nm, respectively.¹²

The RMSE and average relative error (ARE) were used to evaluate the model accuracy, and their respective formulas are as follows:

The determination coefficient $R^2$ was used to evaluate the goodness of fit in the regression model, and the $F$-test and $p$ value were used to evaluate the model significance. To guarantee the model assumption, the residual plots were used to diagnose the regression model before the model was used for prediction.⁷⁶

Two residual plots are important: Q-Q (Quantiles-quantiles) plots and scatter plots between the residuals and the variable, both of which were used in this study. The former was used to check the normality of the residuals, whose scatter points should form an approximately straight line if the distribution is close to the standard normal distribution. The latter was used to check the residual heteroscedasticity: the model must be modified if the residual is a certain function of the estimated Chla.⁷⁷

2.5.

Hyperspectral Data Simulation

To avoid the uncertainty of atmospheric effect, the band reflectance of hyperspectral sensors EO-1/Hyperion and PROBA/CHRIS that are currently in orbit was simulated using the in situ spectra and the spectral response function of each sensor so as to validate the application of the model in satellite images.

The spectral response function of Hyperion can be simulated by a Gaussian spectral response function since each band’s full width at half maximum (FWHM, nm) of this hyperspectral image is narrow.⁷⁸ Assuming the Gaussian peak value at the center wavelength is 1, the spectral response function of Hyperion can be calculated using Eq. (5), where the subscript $i$ represents the sensor band, ${\overline{\lambda }}_i$ is the center wavelength (nm), ${\sigma }_i$ is the band width (nm) and can be calculated from FWHM.

The spectral response function of CHRIS is a strip function determined by the center wavelength and band width [Eq. (6)], whereas each channel’s center wavelength (${\lambda }_i$, nm) and band width ($\mathrm{\Delta }{\lambda }_i$, nm) values can be found in the technical documentation of CHRIS.⁷⁹

Using the spectral response function and in situ spectra data with wavelength interval of 1 nm, the reflectance of channel $i$ in hyperspectral sensors can be calculated according to Eq. (7), where $r_i$ represents the spectral reflectance of channel $i$, ${\lambda }_{si}$ is the initial wavelength of channel $i$, ${\lambda }_{ei}$ is the end wavelength of channel $i$, $r(\lambda )$ is the in situ reflectance at wavelength $\lambda$, and ${\phi }_i(\lambda )$ is the spectral response factor at wavelength $\lambda$ within channel $i$, which can be calculated from the spectral response function.

3.1.

Spectral Reflectance and Constituent Concentration

The spectral magnitudes and shapes of the four datasets after smoothing (Fig. 2) showed similar characteristics to that of typically turbid water:¹⁵ the relative low reflectance in the blue range (400 to 500 nm) was a result of high absorption by water constituents; a slight reflectance trough at around 620 nm was formed by the absorption peak resulting from phycocyanin in water containing blue-green algae;⁸⁰ a second reflectance trough at around 675 nm corresponded to Chla absorption; and a distinct peak at around 700 nm mainly resulted from both chlorophyll fluorescence and minimum absorption by optically active constituents and water.

Fig. 2

Reflectance spectra in July to August of 2004 (a), July to August of 2005 (b), September of 2010 (c), and September of 2011 (d).

Table 3 shows the statistical characteristics of Chla and TSS across the four datasets. The TSS in July to August of 2004 and 2005 can refer to the TSS data in July to August of 1998 to 2003. The CDOM is not considered in this study.

Table 3

Statistical characteristics of Chla in July to August of 2004 and 2005, Chla and total suspended sediment (TSS) in September of 2010 and 2011, and TSS in July to August of 1998 to 2003.

The datasets used encompass widely varying optical conditions and constituent concentrations (Table 3), showing that the water in Taihu Lake is extremely turbid and productive. The Chla ranged from 4 to $192\mathrm{mg}/{\mathrm{m}}^3$ in the four datasets, with averages of 49.88, 25.43, 51.02, and $39.62\mathrm{mg}/{\mathrm{m}}^3$ in the datasets of 2004, 2005, 2010, and 2011, respectively. The high Chla indicates serious eutrophication in Taihu Lake. The TSS was also very high, averaging 49.20, 36.61, and $34.36\mathrm{mg}/\text{L}$ in the datasets of 1998 to 2003, 2010, and 2011, respectively.

3.2.

Spectral Index Sensitive to Chla

Correlation analysis was used to determine the sensitivity of spectral index to Chla. The workflow included determining the relationship between Chla and the spectral reflectance and calculating the correlation coefficient.⁸¹ Given that Pearson’s correlation coefficient only implies a linear relationship, the relation between spectral index and Chla has to be first verified using a scatter plot. A logarithmic relation was found in these datasets; therefore, Chla was transformed into its natural logarithm, denoted as lnChla.

The spectral indices and their correlation coefficients with lnChla in the datasets of 2004 and 2005 were calculated (Fig. 3). The same tendency indicated the consistent change of spectral indices with Chla in these two years.

Fig. 3

Correlation coefficients between spectral indices and lnChla calculated by datasets of 2004 and 2005, Taihu Lake.

Correlation coefficients in 2004 were sorted from high to low and the top 10 are listed in Table 4, where $p<0.001$. With regard to wavelengths, the sensitive reflectance to Chla mainly focused at 450, 550, 670 to 700, 800 nm or the nearby location. In terms of spectral indices, RARSa had the highest sensitivity to Chla. The top three highly correlated spectral indices were all composed by the band ratio of the reflectance peak near 700 nm and the reflectance trough near 670 nm, demonstrating the effectiveness of these two bands for Chla estimation in turbid water.

Table 4

Spectral indices highly correlated with lnChla in the dataset of 2004 (p<0.001).

The combinations of spectral indices, including band difference, band ratio and band normalization, were calculated based on the dataset of 2004.

The top five highly correlated difference and ratio combinations are shown in Table 5, where $p<0.001$. Results show that the difference and ratio of spectral indices generally had higher correlation with Chla than a single spectral index, and the performance of band ratio was slightly better than that of band difference. The ratio of the red/green pigment index³¹ and the ratio analysis of reflectance spectra⁵⁶ (RGI/RARSa) has the highest correlation (0.94) with lnChla, better than RARSa ($-0.91$), as shown in Table 4.

Table 5

Difference and ratio combinations highly correlated with lnChla in the dataset of 2004.

The correlation coefficients between lnChla and the normalized combination of the spectral index pairs were calculated and the top three are listed in Table 6, denoted as NR1, NR2, and NR3. They were highly correlated with lnChla ($\text{correlation coefficient}>0.93$, $p<0.001$).

Table 6

Three normalized combinations highly correlated with lnChla in the dataset of 2004.

3.3.

Estimation of Chla Based on Spectral Index

Based on the same dataset used in Sec. 3.2, the three normalized combinations, NR1, NR2, and NR3 (Table 6), as well as the RARSa (Table 4) and RGI/ RARSa (Table 5) were used for model building. Results are shown in Table 7.

Table 7

Regression models between the spectral indices (NR1, NR2, NR3, RARSa, RGI/RARSa) and lnChla in the dataset of 2004 (p<0.0001).

In the comparison between the $R^2$ value and RMSEs of NR1, NR2, and NR3, NR1 had the best fitting accuracy and NR3 had the lowest RMSE (Table 7). The fitting accuracy of NR1 was slightly better than NR3, but its RMSE and residual variation were larger. The residual diagnostics plots of NR1 and NR3 (Fig. 4) indicated that the residuals of NR1 showed heteroscedasticity as Chla changed, and more points in NR1 deviated from the straight line of its Q-Q plot [Fig. 4(a)]. Comparatively, NR3 fulfilled the requirement of regression analysis better than NR1 [Fig. 4(b)].

Fig. 4

Residual diagnostics of the regression models of NR1 (a) and NR3 (b). The residual diagnostics include the homogeneity test of variance and the normality test of residuals.

The results of RARSa and RGI/RARSa were then compared, and the results of the comparison are also shown in Table 7. Results showed that the RGI/RARSa had a higher fitting accuracy ($R^2=0.93$) than NR3 ($R^2=0.92$). However, the scatter plot between RGI/RARSa and NR3 with lnChla (Fig. 5) showed that the variation of RGI/ RARSa was larger than NR3, indicating that the value of the RGI/RARSa ratio was higher and had a greater variation in sample distribution and that the NR3 value was lower and had a smaller variation. The NR3 was more stable with the change in lnChla.

Fig. 5

Correlation between lnChla and the spectral indices of RGI/RARSa and NR3.

Therefore, NR3 was selected as the optimal index and called the normal chlorophyll index (NCI), $\mathrm{NCI}=(R690/R550-R675/R700)/(R690/R550+R675/R700)$. The Chla estimation model was $\mathrm{Chla}=\exp (7.6334\times \mathrm{NCI}+3.3325)$, with ARE of 18.27%, and RMSE of $14.36\mathrm{mg}/{\mathrm{m}}^3$.

3.4.

Model Validation

Datasets of 2005, 2010, and 2011 were used for model validation. First, Chla was directly estimated using the model of $\mathrm{Chla}=\exp (7.6334\times \mathrm{NCI}+3.3325)$, and their models were denoted as DV1, DV2, and DV3. Then, the model was reparameterized by rebuilding the regression models between the lnChla and NCI derived from the three datasets and then the Chla was estimated. The models were denoted as RV1, RV2, and RV3. Figure 6 shows the scatter plots between the estimated and measured Chla of the three datasets.

Fig. 6

Validation results of the NCI model by datasets of 2005 (a), 2010 (b), and 2011 (c). Validation process includes the direct Chla estimation from the 2004 model (DV1, DV2, and DV3), and the Chla estimation after reparameterization of the NCI model based on new datasets (RV1, RV2, and RV3).

By using the model of 2004 directly, the validation result by dataset of 2005 was satisfactory [DV1 in Fig. 6(a)], with RMSE of $11.17\mathrm{mg}/{\mathrm{m}}^3$ and ARE of 34.75%. This condition may be attributed to the consistency between these two datasets in July to August of 2004 and 2005. The standard deviations of the relative error in the NCI models for 2004 and 2005 were 0.216 and 0.242, respectively. However, although the distribution trend of the model was consistent with the measured Chla, Chla was generally underestimated when the model was directly validated using the datasets of 2010 and 2011 [DV2 in Fig. 6(b) and DV3 in Fig. 6(c)].

When the coefficients of the model were refitted using the calculated NCI and Chla in the new datasets, the estimated and measured Chla in the three datasets were consistent with each other (RV1, RV2, and RV3 in Fig. 6) and most samples were within the error line of $10\mathrm{mg}/{\mathrm{m}}^3$. The RMSEs between the estimated Chla and the measured Chla for 2005, 2010, and 2011 were 10.39, 11.87, and $12.65\mathrm{mg}/{\mathrm{m}}^3$, respectively. These results demonstrated the usability of NCI in the new datasets for the remote sensing mapping of Chla in a water body.

3.5.

Model Application

To verify the application of NCI index model in hyperspectral sensors, the Hyperion and CHRIS band reflectances were first simulated using the in situ spectra data in 2004 and 2005 according to the method mentioned in Sec. 2.5. Both Hyperion and CHRIS have four spectral channels consistent with the 550, 675, 690, and 700 nm of NCI, the basic parameters of which are shown in Table 8.

The NCI was first calculated using the Hyperion band reflectance simulated by the in situ data in 2004 and 2005, whereas data in 2004 were used for model building and data in 2005 were used for model validation. Figure 7(a) shows the regression model between NCI and lnChla in 2004 and also the model accuracy, including $R^2$, RMSE, and ARE. Figure 7(b) shows the validation results by directly using this model in 2005 data and also the validation accuracy of RMSE and ARE, calculated by the estimated results and the measured Chla. Similarly, Fig. 7(c) shows the NCI regression model using the CHRIS band reflectance simulated using the in situ data in 2004 and Fig. 7(d) shows the model validation result in 2005.

Fig. 7

The NCI model and validation results using simulated Hyperion/CHRIS data. Model result in 2004 (a) and validation result in 2005 (c) using simulated Hyperion data; Model result in 2004 (b) and validation result in 2005 (d) using simulated CHRIS data.

The accuracy of the NCI model using hyperspectral data is satisfactory, with $R^2>0.92$ and $\mathrm{RMSE}<17\mathrm{mg}/{\mathrm{m}}^3$. The accuracy of the model built by the simulated hyperspectral data is better than that of using the in situ spectra (Table 7). Meanwhile, the estimation results are also satisfactorily in line with the measured Chla in 2005 using simulation data from both sensors [Figs. 7(c) and 7(d)], with $\mathrm{RMSE}<10\mathrm{mg}/{\mathrm{m}}^3$.

These results indicated that the NCI model can be used to estimate the Chla in Taihu lake using the Hyperion/EO-1 and CHRIS/PROBA spectra data with an accurate atmospheric correction of water body. This index will also have good usability in other hyperspectral satellite images with similar channel settings.

4.1.

Model Performance of NCI

With regard to Taihu Lake, Le et al.¹² stated that the optimal positions for the three-band model are 660, 692, and 740 nm, and the optimal positions for the four-band model are 662, 693, 740, and 705 nm. When the three- and four-band models were directly used without reparameterization of the new datasets in this study, the performance was poor (not shown). This result is similar to the results from previous studies.^82,83

Based on the data obtained in July to August, 2004, the three- and four-band combinations were compared with NCI. The regression parameters between $(1/R660-1/R692)\times R740$ and $(1/R662-1/R693)/(1/R740-1/R705)$ were first refitted with Chla in 2004, and the optimal three- and four-band positions were then tuned to calibrate the model. Table 9 shows the results of the three- and four-band models after reparameterization and model calibration, and Fig. 8 shows their residual diagnostic plots.

Table 8

Hyperion’s and CHRIS’s band parameters near 550,675,690, and 700 nm.

Table 9

Model and test results of the three- and four-band models after reparameterization and calibration in the dataset of 2004 (p<0.0001).

Fig. 8

Test of variance homogeneity and residual normality of the three- and four-band algorithms after model reparameterization (a, c) and band tuning (b, d).

Table 9 shows that the models of the previous three- and four-band combinations after reparameterization ($R^2=0.88$; $\mathrm{RMSE}>15.9\mathrm{mg}/{\mathrm{m}}^3$) were less significant than the NCI model ($R^2=0.92$; $\mathrm{RMSE}=14.36\mathrm{mg}/{\mathrm{m}}^3$). Both the three- and four-band models achieved superior results after model calibration by band tuning, whereas the four-band model ($R^2=0.93$; $\mathrm{RMSE}=12.14\mathrm{mg}/{\mathrm{m}}^3$) had a better estimation than NCI.

However, the residual diagnostics plots (Fig. 8) showed that the residuals of the three- and four-band models were clearly non-normal and increased with Chla. The samples with lower and higher Chla had a larger deviation. Moreover, the error variance was correlated with Chla, indicating the residuals heteroscedasticity and the necessity for data transformation.⁷⁷ Comparatively, the residual diagnostics of the NCI model [Fig. 4(b)] was better, indicating its robust performance in Chla prediction.

The robustness of the NCI model can be explained from two aspects. The first is the logarithmic transformation of Chla before model building, which guarantees data normality. After data transformation, the model is better than directly using regression analysis. Previous studies^83,84 demonstrated the same result.

Regression models between the three- and four-band combinations and lnChla were also constructed to further verify the effect of data transformation. The reparameterized and calibrated three- and four-band models had poorer residual diagnostics results than the NCI, except for the calibrated three-band model (not shown). However, in this calibrated model $(1/R664-1/R702)\times R679$, the third band $\lambda 3$ found in the range of 450 to 800 nm did not confirm the assumption that it is the spectral region minimally affected by pigment absorption and can compensate for the variability in the backscattering between samples.⁹ If the searching ranges of $\lambda 1$, $\lambda 2$, and $\lambda 3$ were reset into 658 to 676 nm, 691 to 735 nm, and 723 to 780 nm, respectively, according to Gitelson et al.,^10,15 although the 664, 716, and 768 nm would be the last optimal bands, the residual diagnostic result of the calibrated model would still be poorer than that of NCI.

The second aspect that can explain the robustness of the NCI model is the stable wavelengths used in this integrated index. Using the characteristic bands obtained from the vegetation index tested by Chappelle et al.⁵⁶ and Zarco-Tejada et al.,³¹ the NCI proposed in this study was robust because it achieved satisfactory results when validated by multiple datasets. Comparatively, the three- or four-band combinations based on a small number of datasets did not perform well during validation, thus requiring band tuning in new datasets. Numerous results, including those of this study and Gitelson et al.,^10,15 showed the necessity of model calibration of the three-band algorithm to make it suitable to new datasets. The three- and four-band models constructed by specific dataset may not fulfill the model assumptions when applied to other water bodies with different optical properties, thus producing unfavorable model diagnostics and inapplicable validation in new datasets. The estimation model based on NCI had a better residual diagnostics result, indicating that it is robust enough to be used in multiple datasets in Taihu Lake.

4.2.

NCI and Chla

The spectral reflectance of green plants is primarily affected by leaf color, cell structure, and plant moisture. The typical spectral curves of green plants are as follows: a small reflectance peak near 550 nm with reflectance of 10% to 20%; two obvious reflectance troughs near 450 and 670 nm owing to the strong absorption of pigment concentration; and a sharp increase of reflectance from 700 to 800 nm, leading to an obvious slope called “red edge.” Given a wide range of leaf greenness, the maximum sensitivity to Chla was found at 550 to 560 nm and 700 to 710 nm;^85,86 its correlation with Chla at these two bands was larger than 675 nm, and the reciprocal reflectance near 550 and 700 nm was proportional to Chla.^18,85 The reflectance at 670 nm decreased sharply when Chla increased up to 3 to $5\mathrm{mol}/{\mathrm{cm}}^2$. Thereafter, R670 was almost pigment-concentration independent. R670 can thus be used as a reference. However, the reflectance peak in the NIR area (larger than 750 nm) was insensitive to Chla.⁸⁵ The ratio may amplify the differences between the spectra at specific bands due to absorption maxima and minima of the photosynthetic pigments.⁵⁶ The band ratio of 675 and 700 nm (RARSa) was found to have a strong linear relationship with Chla ($R^2=0.93$) in the leaves of soybean.⁵⁶ The red and green pigment index RGI ($R690/R550$) was used to estimate the pigment concentration of the leaf.³¹

The spectrum above water surface is affected by the absorption and scattering of water and particles in it, primarily including phytoplankton, suspended sediment, and CDOM. The four characteristic bands related to chlorophyll are the blue, green, red, and NIR bands in water color remote sensing.^7,87–89 The ratio of the reflectance peak near 700 nm and the reflectance trough near 675 nm are proven to be most sensitive to Chla in turbid waters of different trophic states,^6,17 which makes RARSa ($R700/R675$) reasonable. Moreover, the R690 in NCI can be regarded as the fluorescence peak because the fluorescence peak is around 690 nm, but can vary from 685 to 695 nm,⁹⁰ or be constant at $\sim 685\mathrm{nm}$.⁹¹ Studies showed that normalizing the fluorescence peak near 700 nm to the reflectance value of the global maximum of the spectrum at green peak can accurately predict Chla,^6,70 making RGI ($R690/R550$) reasonable.

Although RARSa and RGI are proven to be theoretically effective in Chla estimation, the normalized combination of these two indices-NCI is optimal because of its robust performance as Chla changes (Fig. 5). The first advantage of NCI is the information supplementary of low Chla in the hyperspectral reflectance that cannot be fully expressed by the red-NIR spectrum by using the green band. Accurate estimates of low Chla using fluorescence algorithms, such as FLH and maximum chlorophyll index,⁹² or the red-NIR algorithms,¹⁷ are unavailable in oligotrophic and some mesotrophic lakes. Thus, the green band was additionally used because the blue and green bands in the OC2-OC4 algorithms were successfully applied to retrieve 0 to $10\mathrm{mg}/{\mathrm{m}}^3$ of Chla in optically complex waters. The second advantage of the NCI is that the feature wavelengths of 550, 675, 690, and 700 nm have all been set in hyperspectral sensors, such as EO-1/Hyperion, HJ1/HSI, PHILLS/HICO, and so on. Compared with the NCI, the commonly used three- or four- band models always require a reflectance longer than 750 nm. However, obtaining reliable estimates of reflectance at these wavelengths is a difficult task based on existing atmospheric correction schemes for turbid waters.¹⁶ Moreover, unlike previous studies on Chla estimation based on the spectra data obtained from July to August, 2004, in Taihu Lake, including the $\mathrm{NIR}/R$ ratio,⁹³ the three-band combination⁹⁴or the band normalized combination,⁸⁴ the NCI model proposed in this study was tested by regression diagnostics and validated by multiple datasets, thus guaranteeing more reliability.