Multi-Algorithm Indices and Look-Up Table for Chlorophyll-a Retrieval in Highly Turbid Water Bodies Using Multispectral Data

Many approaches have been proposed for monitoring the eutrophication of Case 2 waters using remote sensing data. Semi-analytical algorithms and spectrum matching are two major approaches for chlorophyll-a (Chla) retrieval. Semi-analytical algorithms provide indices correlated with phytoplankton characteristics, (e.g., maximum and minimum absorption peaks). Algorithms’ indices are correlated with measured Chla through the regression process. The main drawback of the semi-analytical algorithms is that the derived relation is location and data limited. Spectrum matching and the look-up table approach rely on matching the measured reflectance with a large library of simulated references corresponding to wide ranges of water properties. The spectral matching approach taking hyperspectral measured reflectance as an input, leading to difficulties in incorporating data from multispectral satellites. Consequently, multi-algorithm indices and the look-up table (MAIN-LUT) technique is proposed to combine the merits of semi-analytical algorithms and look-up table, which can be applied to multispectral data. Eight combinations of four algorithms (i.e., 2-band, 3-band, maximum chlorophyll index, and normalized difference chlorophyll index) are investigated for the MAIN-LUT technique. In situ measurements and Medium Resolution Imaging Spectrometer (MERIS) sensor data are used to validate MAIN-LUT. In general, the MAIN-LUT provide a comparable retrieval accuracy with locally tuned algorithms. The most accurate of the locally tuned algorithms varied among datasets, revealing the limitation of these algorithms to be applied universally. In contrast, the MAIN-LUT provided relatively high retrieval accuracy for Tokyo Bay (R2 = 0.692, root mean square error (RMSE) = 21.4 mg m−3), Lake Kasumigaura (R2 = 0.866, RMSE = 11.3 mg m−3), and MERIS data over Lake Kasumigaura (R2 = 0.57, RMSE = 36.5 mg m−3). The simulated reflectance library of MAIN-LUT was generated based on inherent optical properties of Tokyo Bay; however, the MAIN-LUT also provided high retrieval accuracy for Lake Kasumigaura. MAIN-LUT could capture the spatial and temporal distribution of Chla concentration for Lake Kasumigaura.


Introduction
Ocean color satellites have provided synoptic observation for oceans, coastal areas and inland water bodies [1,2].Coastal Zone Color Scanner (CZCS) was the first ocean color satellite, which enables the ocean community to map the spatial distribution of ocean properties [3], such as the diffuse attenuation coefficient [4] and chlorophyll-a (Chla) concentration [5].After a gap of 10 years from the end of CZCS mission, the Sea-viewing Wide Field-of-view Sensor (SeaWiFS) was launched in 1997.Since then, continuous monitoring of water bodies is achieved through many satellites [6,7].Most of these satellites provided multispectral data with reasonable spatial and temporal resolution [8], and it is freely available [9,10].Monitoring of Chla concentration in Case 1 waters (e.g., open oceans) using satellite sensors has been established due to the fact that the phytoplankton is governing the optical properties [11,12].Several approaches were developed for Case 1 waters such as quasi-analytical algorithm (QAA) [13], the Garver-Siegel-Maritorena (GSM) algorithm [14], and the band ratio in the blue-green wavelengths as ocean color V4 (OC4E) algorithm [15].Both QAA and GMS derive inherent optical properties of open oceans, while the difference between then is that GMS requires the absorption spectral shape of water constituents, while QAA does not require these data [16].
Lee et al. [48] developed a semi-analytical inversion model for shallow water to derive water depth, and water column properties (i.e., absorption and backscattering coefficients) of the water column by optimizing the difference between measured reflectance (R rs ) and simulated R rs .Many researchers upgraded the model of Lee et al. [48] to involve many bottom types instead of sand-bottom type [28,49] and to retrieve Chla, NAP, and CDOM concentrations [30].These inversion models were mainly applied to shallow, and relatively clear water with high CDOM concentrations in the Bahamas and Moreton Bay in eastern Australia using airborne hyperspectral data with a high spatial resolution [50][51][52][53].The airborne hyperspectral data is much more expensive with small area coverage compared with multispectral satellite dates that are freely available for many satellites with routine global coverage [54].Regarding spectrum matching approach, a matchup process is performed between the measured R rs and a library of simulated R rs covering different water conditions to find the closest simulated R rs with its known properties [55].Similarly, a spectrum matching approach requires hyperspectral data to perform the matching process.
Accordingly, a new technique is required to combine the advantages of semi-analytical algorithms (i.e., providing indices correlated with Chla concentration) with the merits of the look-up table (i.e., providing dataset corresponding to various water properties) and applicability to multispectral data.The current study proposed multi-algorithm indices and the look-up table technique.The aim of this study is to (1) develop a new technique for Chla retrieval based on algorithms' indices and look-up tables; (2) select the proper combination of algorithms for the proposed technique; and (3) evaluate the performance of proposed technique using in situ measurements and Medium Resolution Imaging Spectrometer (MERIS) sensor data over turbid and highly turbid water bodies.

In Situ and Laboratory Measurements
Data collected from Tokyo Bay (35 • 25 N; 139 • 47 E), and Lake Kasumigaura (36 • 03 N; 140 • 25 E) were used to validate the proposed technique (Figure 1).In situ remote sensing reflectance spectra along with water quality parameters (e.g., Chla concentrations) were collected for Tokyo Bay from June 2010 to August 2013 (sixteen campaigns) and for Lake Kasumigaura in 2016 (eight campaigns).Figure 2 illustrates the remote sensing reflectance spectra for both sites.The red reflectance spectra shown in Figure 2b represent 12 stations of Tokyo Bay where additional inherent optical measurements (i.e., absorption and backscattering coefficients) were executed.Furthermore, water quality parameters were also collected monthly since 1977 by the National Institute for Environmental Studies (NIES) [56].Monthly campaigns of NIES were compared with acquired MERIS data to extracted matched dates, resulting in 77 stations.Table 1 summarized the basic characteristics of collected data, revealing that Lake Kasumigaura is a highly turbid water body comparing with Tokyo Bay.The Chla concentrations were fluorometrically measured using a Turner Designs 10-AU fluorometer [57].A 20 mL of water samples was filtered through 25-mm Whatman GF/F filters of a 0.7-µm pore size, then the pigments were extracted by soaking the filter in 6 mL of N,N-dimethylformamide (DMF) and storing it in dark at 4 • C for 4 h [58].The total suspended solids (TSS), organic suspended solids (OSS), and inorganic suspended solids (ISS) were determined gravimetrically [59].

In Situ and Laboratory Measurements
Data collected from Tokyo Bay (35°25′N; 139°47′E), and Lake Kasumigaura (36°03′N; 140°25′E) were used to validate the proposed technique (Figure 1).In situ remote sensing reflectance spectra along with water quality parameters (e.g., Chla concentrations) were collected for Tokyo Bay from June 2010 to August 2013 (sixteen campaigns) and for Lake Kasumigaura in 2016 (eight campaigns).Figure 2 illustrates the remote sensing reflectance spectra for both sites.The red reflectance spectra shown in Figure 2b represent 12 stations of Tokyo Bay where additional inherent optical measurements (i.e., absorption and backscattering coefficients) were executed.Furthermore, water quality parameters were also collected monthly since 1977 by the National Institute for Environmental Studies (NIES) [56].Monthly campaigns of NIES were compared with acquired MERIS data to extracted matched dates, resulting in 77 stations.Table 1 summarized the basic characteristics of collected data, revealing that Lake Kasumigaura is a highly turbid water body comparing with Tokyo Bay.The Chla concentrations were fluorometrically measured using a Turner Designs 10-AU fluorometer [57].A 20 mL of water samples was filtered through 25-mm Whatman GF/F filters of a 0.7-µm pore size, then the pigments were extracted by soaking the filter in 6 mL of N,N-dimethylformamide (DMF) and storing it in dark at 4 °C for 4 h [58].The total suspended solids (TSS), organic suspended solids (OSS), and inorganic suspended solids (ISS) were determined gravimetrically [59].In situ remote sensing reflectance spectra were collected using two different instruments: (1) three TriOS-RAMSES hyperspectral radiometers (Oldenburg, Germany) in the spectral range between 350 nm and 900 nm at a 2 nm spectral interval and a field-of-view of 7° used for Tokyo Bay field observation, and (2) Fieldspec HandHeld 2 Spectroradiometer (Boulder, CO, USA) in the spectrum range from 350 to 1075 nm with spectral resolution of 1.0 nm and 25° field-of-view used for Lake Kasumigaura.The protocol proposed by Mobley [60] was adopted for HandHeld 2. The absorption coefficients of phytoplankton, non-algal particles, and colored dissolved organic matter were measured with the quantitative filter technique (QFT) approach [61] using a UV/VIS spectrophotometer (JASCO, V550, Hachioji, Japan).The vertical profiles of the total backscattering coefficients ( ( )) were measured at six wavelengths (420, 442, 488, 510, 550 and 676 nm) with a HydroScat-6 (HOBI Labs Inc., Tucson, AZ, USA).The backscattering coefficients for particles were estimated by subtracting the backscattering of pure water [62,63] from the total backscattering In situ remote sensing reflectance spectra were collected using two different instruments: (1) three TriOS-RAMSES hyperspectral radiometers (Oldenburg, Germany) in the spectral range between 350 nm and 900 nm at a 2 nm spectral interval and a field-of-view of 7 • used for Tokyo Bay field observation, and (2) Fieldspec HandHeld 2 Spectroradiometer (Boulder, CO, USA) in the spectrum range from 350 to 1075 nm with spectral resolution of 1.0 nm and 25 • field-of-view used for Lake Kasumigaura.The protocol proposed by Mobley [60] was adopted for HandHeld 2. The absorption coefficients of phytoplankton, non-algal particles, and colored dissolved organic matter were measured with the quantitative filter technique (QFT) approach [61] using a UV/VIS spectrophotometer (JASCO, V550, Hachioji, Japan).The vertical profiles of the total backscattering coefficients (b b (λ)) were measured at six wavelengths (420, 442, 488, 510, 550 and 676 nm) with a HydroScat-6 (HOBI Labs Inc., Tucson, AZ, USA).The backscattering coefficients for particles were estimated by subtracting the backscattering of pure water [62,63] from the total backscattering coefficients (b b (λ)) [64].

Bio-Optical Model
Subsurface remote sensing reflectance ( ) ( r rs λ , all symbols are defined in Table 2) was fitted with inherent optical properties (IOPs) using polynomial function by Gordon et al. [65], which can be Table 1.Descriptive statistics of the water quality parameters for Tokyo Bay, and Lake Kasumigaura.Tokyo Bay with IOPs (inherent optical properties) dataset refers to stations where additional measurements of absorption and backscattering properties were performed along with remote sensing reflectance measurements.Lake Kasumigaura Datasets 1 and 2 represent measurements collected in 2016 and during 2002-2012, respectively.

Bio-Optical Model
Subsurface remote sensing reflectance (r rs (λ), all symbols are defined in Table 2) was fitted with inherent optical properties (IOPs) using polynomial function by Gordon et al. [65], which can be simplified as where f is the light field factor; and Q denotes the light distribution factor.Kirk [66] found that f is a function of the solar altitude, which is expressed as a linear function of the cosine of the zenith angle of the refracted photons The average value of f /Q (0.09) was determined by using the bio-optical model in stations that have both measured R rs and IOPs in Tokyo Bay, which is consistent with that in Mishra et al. [11].Both a(λ) and b b (λ) are the sum of the contributions of pure water and three optically active constituents in water: phytoplankton, non-algal particles and colored dissolved organic matter.The a(λ), and b b (λ) can be expressed as where a w (λ) were taken from Buiteveld et al. [67] and b b,w (λ) were taken from Dall'Olmo et al. [63] and Morel [62].The radiance just above water surface was related to the radiance just below water surface by a factor of 0.544 [68].Hence, the remote sensing reflectance just above the water surface, R rs (λ), can be calculated as Table 2. Symbols and definitions.

Simulating the Remote Sensing Reflectance
Researchers have proposed many models to generate simulated R rs based on bio-optical models [47,69,70].The key difference among these models is the parameters of specific absorption and backscattering [9].The chlorophyll concentration and a ph (λ) have a proportional relationship [71] that can be expressed as The NAP and CDOM have similar trends for absorption spectrum, which can be described as exponential decay function with increasing wavelength.The a N AP (λ), and a CDOM (λ) are estimated from their reference values at 440-nm as follows: where the values of a * N AP (440) and a * CDOM (440) were 0.03483 and 1.0 m −1 , respectively.The slopes S N AP and S CDOM equaled 0.00899 and 0.01547, respectively.The backscattering coefficient of particles, b b,p (λ), at different wavelengths was estimated by fitting the measured backscattering at the six wavelengths with an exponential function: where n stands for the spectral slope of b b,p (λ) and b b,p (550) refers to the backscattering at the 550 nm band.The values for b b,p (550) and n were 0.02141 m −1 and 1.25848, respectively.The segregation of b b,p (λ) into the backscattering of phytoplankton (b b,ph (λ)) and non-algal particles (b b,N AP (λ)) was based on the assumption that their contributions to b b,p (λ) were proportionally correlated to their concentrations [47].Separating NAP from the TSS with laboratory analysis is impossible [72], so the approach of Gons et al. [73] was employed to divide the TSS into NAP and phytoplankton suspended solids (PSS).As a result, it was assumed that the 1.0 mg m −3 Chla concentration approximately equaled 0.0687 g m −3 TSS.The NAP concentration was calculated by subtracting the PSS (=0.0687 × Chla) from the TSS.The reflectance caused by chlorophyll fluorescence emission for Case 2 waters was also considered as proposed by Gilerson et al. [71].Fluorescence reflectance follows a Gaussian shape with a peak at 685 nm, a full-width half maximum of 25 nm and a standard deviation, σ, of 10.6 nm [74].Thus, the chlorophyll fluorescence emission reflectance R rs, f l (λ) was modeled as where E d (685) was 1.1 W m −2 nm −1 under clear sky conditions Based on the Hydrolight simulation [71], and f l(685) is the chlorophyll fluorescence at 685 nm, which was calculated as The total remote sensing reflectance was generated as the sum of the reflectance in Equations ( 5) and (10) [11].

MERIS Images Processing
MERIS L1b products of 300 m full resolution from European Space Agency (ESA) Earthnet Online [75] were processed using BEAM V5.0 software (Envisat/Brockman Consult, Hamburg, Germany).All clear and partially clear MERIS images that contain Lake Kasumigaura and matched the monthly field observations from 2002 to 2012 were extracted.The matchup process included all the images acquired at the same day or one-day difference from field observations.It was assumed that the changes in Chla concentrations within one day would not be significant due to the stability of weather conditions (e.g., wind speed) during that day (Table S1 in Supplementary Material).The matchup process resulted in a total of 13 images (i.e., five images at the same day and eight images within a one-day difference).
All images were pre-processed using a Smile Correction processor to adjust the variations in the spectral wavelengths among MERIS's five cameras and within the same camera [76].In order to convert MERIS L1b top-of-atmosphere radiances to water-leaving reflectance, the Case-2 Regional (C2R), the Eutrophic Lake (EUL) and Boreal Lake (BL) processors are neural network modules for performing atmospheric correction in Case 2 waters [35,77].The three processors share the same neural network atmospheric correction module, while the difference among them is the training data that used to trained inverse and forward coupled neural network to retrieve the IOPs and constituent concentrations from the atmospherically corrected reflectance [35,77].During this study, the water-leaving reflectance is only needed from MERIS data, while the Chla concentrations will be taken from the matched field observations, not from the neural network modules.Thus, the EUL processor was employed to perform atmospheric correction [77].The final number of matched stations decreased from 130 to 77 stations due to quality control flags in L1b product (e.g., invalid, bright, suspect, land, or glint risk) or during the atmospheric correction process.Figure 3 shows the water leaving reflectance of the 77 stations resulting from the matching process over Lake Kasumigaura.
Remote Sens. 2017, 9, 556 7 of 21 the images acquired at the same day or one-day difference from field observations.It was assumed that the changes in Chla concentrations within one day would not be significant due to the stability of weather conditions (e.g., wind speed) during that day (Table S1 in Supplementary Material).The matchup process resulted in a total of 13 images (i.e., five images at the same day and eight images within a one-day difference).All images were pre-processed using a Smile Correction processor to adjust the variations in the spectral wavelengths among MERIS's five cameras and within the same camera [76].In order to convert MERIS L1b top-of-atmosphere radiances to water-leaving reflectance, the Case-2 Regional (C2R), the Eutrophic Lake (EUL) and Boreal Lake (BL) processors are neural network modules for performing atmospheric correction in Case 2 waters [35,77].The three processors share the same neural network atmospheric correction module, while the difference among them is the training data that used to trained inverse and forward coupled neural network to retrieve the IOPs and constituent concentrations from the atmospherically corrected reflectance [35,77].During this study, the waterleaving reflectance is only needed from MERIS data, while the Chla concentrations will be taken from the matched field observations, not from the neural network modules.Thus, the EUL processor was employed to perform atmospheric correction [77].The final number of matched stations decreased from 130 to 77 stations due to quality control flags in L1b product (e.g., invalid, bright, suspect, land, or glint risk) or during the atmospheric correction process.Figure 3 shows the water leaving reflectance of the 77 stations resulting from the matching process over Lake Kasumigaura.

Accuracy Assessment of the Algorithms
The retrieval accuracy was assessed using three indices: coefficient of determination (R 2 ) [78] the root mean square error (RMSE) [79], and the mean absolute relative error (MARE) [79], as follows:

Accuracy Assessment of the Algorithms
The retrieval accuracy was assessed using three indices: coefficient of determination (R 2 ) [78] the root mean square error (RMSE) [79], and the mean absolute relative error (MARE) [79], as follows: where Chla i,meas and Chla i,retr refer to the measured and retrieved Chla concentrations, respectively; and N is the total number of samples.

Construction of the Multi-Algorithm Indices and Look-Up Table (MAIN-LUT)
The proposed technique relies on creating a library of simulated reflectance dataset covering wide ranges of Chla, NAP and CDOM concentrations.The simulated reflectance dataset was generated based on the bio-optical model and specific inherent optical properties (SIOPs) from Tokyo Bay (Sections 2.2 and 2.3).The bio-optical modeling is used to link remote sensing reflectance with optically active constituents in the water column (i.e., Chla, NAP, and CDOM).By changing the concentration of any of these three constituents, remote sensing reflectance will change.Consequently, a simulated dataset of 500,000 R rs spectra was generated covering wide ranges of Chla (1-200 mg m −3 ) (hereafter, reference Chla), NAP (1-200 g m −3 ), and CDOM (1-10 m −1 ) concentrations.Increments of 2.0 mg m −3 , 2.0 g m −3 , and 2.0 m −1 were used for Chla, NAP and CDOM concentrations, respectively.Each of the generated R rs spectra were tagged with the Chla, NAP, and CDOM concentrations, which were used during the generation of simulated R rs spectra.
The major limitation to perform spectral matching with multi-spectral data is associated with the limited number of available bands, particularly when the bands of the blue spectrum were excluded due to difficulties in atmospheric correction for that part of the spectrum [80].On the other hand, there was no relation between the reference Chla and individual bands of the multi-spectral data (e.g., 665-nm band) as shown in Figure S1 (Supplementary Material).Many algorithms have been proposed to retrieve Chla concentration in Case 2 waters [11, 21,24,25].These algorithms provide indices correlated with Chla concentrations [37].Four algorithms were proposed for the proposed technique (i.e., 2-band ratio (Equation ( 15)), 3-band algorithm (Equation ( 16)), maximum chlorophyll index (Equation ( 17)) and normalized difference chlorophyll index (Equation ( 18)): The wavelengths incorporated with each algorithm were summarized in Table 3.The eight algorithms mentioned in Table 3 represent the most accurate algorithms among 43 algorithms' combinations that were evaluated by Salem et al. [81] using measured and simulated datasets.In addition, these eight algorithms do not require hyperspectral R rs .Scatter plots between the reference Chla and the above-mentioned algorithms with different wavelengths (Table S3) were performed (Figure 4).Investigated algorithms provided correlation with reference Chla: high correlation (e.g., 2b-[665, 709]), moderated correlation (e.g., 3b-[680, 709, 754]), and low correlation (e.g., MCI [680, 709, 754]) (Figure 4).Thus, the index of these algorithms could be used to retrieve Chla concentration of measured R rs by comparing the index of measured R rs with the indices of simulated R rs to extract the simulated R rs with the same index.The problem is that there are many simulated R rs that have same measured index during the matchup process.Thus, multi-algorithm indices should be used to obtain the simulated R rs , which is close to the measured R rs .Different combinations of these four algorithms were examined with the MAIN-LUT technique to find the best combination.Table 4 summarized eight combinations that were investigated.The flowchart in Figure 5 illustrates the concept of MAIN-LUT using the 2-indices-[665] combination as an example, whereas the major steps of the MAIN-LUT as follows: Step (0) two things should be prepared and selected before using MAIN-LUT: (1) a library of simulated R rs should be provided through a look-up table; and (2) one combination of algorithms should be selected from Table 4.The selected combination is hereafter called Sample Combination through the coming steps.
Step (2) for each measured index estimated in Step (1), simulated reflectance spectra from a look-up table (LUT) that have the same index should be extracted.Of course, many simulated reflectance spectra with different tagged Chla will be extracted from this matching process.For example, the simulated reflectance (Rrs simu rs,200 , Rrs simu rs,1520 ) were extracted because they matched the 3b meas 665 .The 200 and 1520 represent the tagged numbers of simulated reflectance for a certain combination of Chla, NAP, and CDOM and these tagged numbers range from 1 to 500,000.
Step (3) grouping of simulated reflectance extracted from Step (2), which was generated from the matching process of each measured index.
Step (4) estimating algorithm indices for each of simulated reflectance grouped in Step (3).The number of required indices is based on the Sample Combination selected in Step (0).
Step (5) The RMSE is used to compare the indices of measured reflectance with simulated indices of each extracted reflectance in Step (3) as follows: where Index i,meas denote measured indices (i.e., 2b meas 665 and 3b meas 665 for 2b-indices-[665] combination); Index * Wavelengths incorporated with each algorithm were summarized in Table 3.  3.  3. represents the tagged number for each simulated Rrs that range from 1 to 500,000; and meas and simu denote measured and simulated, respectively.The blue and red colors highlight the measured and simulated parameters, respectively.

Validation Using In Situ Measurements
The in situ data collected from Tokyo Bay between 2010 and 2013 (71 stations) were used to investigate the retrieval accuracy of the MAIN-LUT technique.The MAIN-LUT technique does not require a calibration process because Chla concentrations are directly retrieved from the look-up table.Therefore, all of the measured Rrs passed through the MAIN-LUT steps (from step 1 to step 5 as explained in previous section) to retrieve Chla concentrations.For the eight algorithms'

Validation Using In Situ Measurements
The in situ data collected from Tokyo Bay between 2010 and 2013 (71 stations) were used to investigate the retrieval accuracy of the MAIN-LUT technique.The MAIN-LUT technique does not require a calibration process because Chla concentrations are directly retrieved from the look-up table.Therefore, all of the measured R rs passed through the MAIN-LUT steps (from step 1 to step 5 as explained in previous section) to retrieve Chla concentrations.For the eight algorithms' combinations proposed in Table 4, the Chla concentrations for each combination were retrieved from the MAIN-LUT and compared with the measured Chla (Figure 6).In general, the results of MAIN-LUT revealed that most of algorithms' combinations could accurately retrieve Chla concentrations (Figure 6).6g,h).This can be attributed to the fact that the 665 nm is located near the red Chla absorption maximum, which allow algorithms incorporated a 665 nm band to capture the variation of Chla better than that incorporated 680 nm band.
combinations proposed in Table 4, the Chla concentrations for each combination were retrieved from the MAIN-LUT and compared with the measured Chla (Figure 6).In general, the results of MAIN-LUT revealed that most of algorithms' combinations could accurately retrieve Chla concentrations (Figure 6).6g,h).This can be attributed to the fact that the 665 nm is located near the red Chla absorption maximum, which allow algorithms incorporated a 665 nm band to capture the variation of Chla better than that incorporated 680 nm band.4.
The performance of the MAIN-LUT was also examined with in situ measurements from Lake Kasumigaura collected in 2016 with total sampling points of 68 stations.Although the simulated R rs data of the MAIN-LUT technique was generated based on IOP values from Tokyo Bay, the MAIN-LUT produced higher Chla retrieval accuracy for Lake Kasumigaura than Tokyo Bay for different algorithms' combinations (Figures 6 and 7).This result can be attributed to the fact that the four algorithms incorporated with MAIN-LUT (Equations ( 15)-( 18)) are mainly developed for Case 2 waters with high Chla and TSS concentrations.Comparing the characteristics of these two water bodies revealed that water conditions in Lake Kasumigaura (Chla = 55.86 mg m −3 , TSS = 21.71g m −3 , Table 1) are more correlated with the algorithms' assumption than Tokyo Bay (Chla = 25.14 mg m −3 , TSS = 8.29 g m −3 ).In addition, the MARE values for all combinations were very high for Tokyo Bay comparing with the same combinations at Lake Kasumigaura (Figures 6 and 7).The reason is mainly due to the fact that the MARE is a relative error that assesses the Chla retrieval over a relatively low range of Chla concentration in Tokyo Bay.As a result, any small difference between the measured and retrieved Chla concentrations would cause a high relative error.
Remote Sens. 2017, 9, 556 13 of 21 The performance of the MAIN-LUT was also examined with in situ measurements from Lake Kasumigaura collected in 2016 with total sampling points of 68 stations.Although the simulated Rrs data of the MAIN-LUT technique was generated based on IOP values from Tokyo Bay, the MAIN-LUT produced higher Chla retrieval accuracy for Lake Kasumigaura than Tokyo Bay for different algorithms' combinations (Figures 6 and 7).This result can be attributed to the fact that the four algorithms incorporated with MAIN-LUT (Equations ( 15)-( 18)) are mainly developed for Case 2 waters with high Chla and TSS concentrations.Comparing the characteristics of these two water bodies revealed that water conditions in Lake Kasumigaura (Chla = 55.86 mg m −3 , TSS = 21.71g m −3 , Table 1) are more correlated with the algorithms' assumption than Tokyo Bay (Chla = 25.14 mg m −3 , TSS = 8.29 g m −3 ).In addition, the MARE values for all combinations were very high for Tokyo Bay comparing with the same combinations at Lake Kasumigaura (Figures 6 and 7).The reason is mainly due to the fact that the MARE is a relative error that assesses the Chla retrieval over a relatively low range of Chla concentration in Tokyo Bay.As a result, any small difference between the measured and retrieved Chla concentrations would cause a high relative error.The highest retrieval accuracy for Lake Kasumigaura was achieved by 4-indices-[2b & 3b] combination with R 2 of 0.866, RMSE of 11.3 mg m −3 , and MARE of 16.3% (Figure 7c).The top three combinations in terms of Chla retrieval were 4-indices-[2b & 3b], 8-indices, and 6-indices with average R 2 , RMSE, and MARE of 0.866, 11.23 mg m −3 and 16.47%, respectively.The second tier of algorithms' combinations were 3-indices-[665], 4-indices-[665], and 2-indices-[665] which introduced slightly high retrieval accuracy, with R 2 , RMSE and MARE of 0.840, 12.73 mg m −3 and 20.3%, respectively.The retrieval accuracy of a 4-indices-[680] combination significantly improved in terms of R 2 (from 0.287 to 0.703) and MARE (from 231.7 to 29.3%) (Figures 6g and 7g).The 3-indices-[680] combination was also improved compared with Tokyo Bay (Figures 6h and 7h).

Validation Using MERIS Data
MAIN-LUT technique was also validated using MERIS sensor data.A matchup process executed between field observation campaigns in Lake Kasumigaura and MERIS images acquired at same day or one-day differences during 2002-2012, resulting in 130 matched stations.However, the final number of stations that was used during the validation process reduced to 77 stations (Table 1-Dataset 2), due to flags occurring through image processing.Therefore, MERIS reflectance spectra were used as an input for MAIN-LUT to retrieve Chla concentrations, which were compared with the measured Chla during the field campaigns.In general, all algorithms' combinations provided relatively acceptable Chla retrieval, even for 4-indices-[680] and 3-indices-[680], which introduced the lowest retrieval accuracy for Tokyo Bay and Lake Kasumigaura (Figure 8).This can be attributed to the fact that the neural network atmospheric correction module does not include the inelastic reflectance (e.g., fluorescence) for the simulated dataset that was used to train the module, which will affect the accuracy of the water-leaving reflectance of MERIS data in the red and near infrared (Red-NIR) regions of the spectrum [76].
The 4-indices-[680] provided the highest accuracy in terms of RMSE of 28.9 mg m −3 .Furthermore, the disruption of scatter plot of the 4-indices-[680] combination was close to a 1:1 line (Figure 8g) compared to most of the other combinations that underestimated Chla >50 mg m −3 (Figure 8a-f).This is mainly due to the fact that the simulated dataset of the MAIN-LUT considered the fluorescence reflectance in contrast to the EUL atmospheric correcting module that neglected the fluorescence reflectance during the training process.Thus, the magnitude of atmospherically corrected reflectance, particularly in the Red-NIR, is lower than the simulated reflectance of MAIN-LUT for the same Chla concentration.The 2-indices-[665] outperformed other combinations with R 2 of 0.57.The retrieval accuracy of 2-indices-[665], 3-indices-[665], and 4-indices-[665] were the highest combinations in terms of R 2 and MARE (R 2 = 0.564, RMSE = 32.47 mg m −3 , MARE = 43.4%)(Figure 8d-f).However, the scatter plot distribution of 3-indices-[680] was also close to the 1:1 line (Figure 8h), and it produced the lowest retrieval accuracy (R 2 = 0.385, RMSE = 37.8 mg m −3 , MARE = 74.3%).
The MAIN-LUT technique was also applied to MERIS images (pixel by pixel) to retrieve and illustrate the spatial distribution of Chla concentration along Lake Kasumigaura.Figure 9a,b 4.

Comparing MAIN-LUT with Locally Tuned Algorithms
All algorithms incorporated with algorithms' combinations of MAIN-LUT (Table 4) were locally tuned and compared with MAIN-LUT.Each dataset (Tokyo Bay, Lake Kasumigaura datasets 1 and 2) shown in Table 1 was randomly divided to 70% and 30% for calibration and validation processes, respectively.Linear regression was adopted during the calibration process.Table 5 summarizes R 2 , slope and intercept derived during the calibration process for each dataset.The results of the MAIN-LUT technique are not shown in Table 5because the MAIN-LUT does not require the calibration stage.The result from the validation process showed that the locally tuned algorithms provided slightly higher retrieval accuracy compared with MAIN-LUT technique (Table 6).However, the MAIN-LUT has two major advantages over the locally tuned algorithms: (1) the MAIN-LUT does not require the calibration process, and (2) the MAIN-LUT introduced comparatively high retrieval accuracy for in situ and satellite datasets, even with a look-up table generated from Tokyo Bay.Moreover, locally tuned algorithms produced different relationships for different datasets at the same algorithms, which emphasize the limitation of derived relationships to be applied for different locations or conditions.

Comparing MAIN-LUT with Locally Tuned Algorithms
All algorithms incorporated with algorithms' combinations of MAIN-LUT (Table 4) were locally tuned and compared with MAIN-LUT.Each dataset (Tokyo Bay, Lake Kasumigaura datasets 1 and 2) shown in Table 1 was randomly divided to 70% and 30% for calibration and validation processes, respectively.Linear regression was adopted during the calibration process.Table 5 summarizes R 2 , slope and intercept derived during the calibration process for each dataset.The results of the MAIN-LUT technique are not shown in Table 5 because the MAIN-LUT does not require the calibration stage.In general, the algorithm that provided the highest correlation varied among datasets (Table 5).The 3b [680, 709, 754], 3b [665, 709, 754], and 2b [680, 709] algorithms outperformed other algorithms for Tokyo Bay, Lake Kasumigaura dataset 1 and Lake Kasumigaura dataset 2, respectively.All algorithms of Lake Kasumigaura dataset 1 comparatively introduced higher correlation than that for Tokyo Bay.These results were mainly due to the fact that these algorithms are much more suitable for high-turbidity water in Lake Kasumigaura than Tokyo Bay.
The result from the validation process showed that the locally tuned algorithms provided slightly higher retrieval accuracy compared with MAIN-LUT technique (Table 6).However, the MAIN-LUT has two major advantages over the locally tuned algorithms: (1) the MAIN-LUT does not require the calibration process, and (2) the MAIN-LUT introduced comparatively high retrieval accuracy for in situ and satellite datasets, even with a look-up table generated from Tokyo Bay.Moreover, locally tuned algorithms produced different relationships for different datasets at the same algorithms, which emphasize the limitation of derived relationships to be applied for different locations or conditions.

Conclusions
A multi-algorithm indices and look-up table (MAIN-LUT) technique was proposed to retrieve Chla concentration for optically complex Case 2 waters.The MAIN-LUT technique combines the advantage of a look-up table that covers wide ranges of trophic status with the ability of the existing Chla algorithms that provide indices correlated with Chla concentration.Multi-algorithms should be used to match algorithms' indices estimated from measured reflectance with the indices of simulated reflectance library.In situ measurements and MERIS data were used to validate MAIN-LUT over turbid and highly turbid water bodies.The results presented in this study reveal the ability of MAIN-LUT to accurately retrieve Chla concentration.Although the simulated dataset was generated based on IOPs from Tokyo Bay, the MAIN-LUT provided high retrieval accuracy for Lake Kasumigaura with in situ measurements and MERIS data.Eight algorithms' combinations were investigated for the MAIN-LUT.Six out of the eight combinations provided comparably high retrieval accuracy for Tokyo Bay (R 2 = 0.65, RMSE = 22.2 mg m −3 , MARE = 148.27%),Lake Kasumigaura (R 2 = 0.85, RMSE = 11.98 mg m −3 , MARE = 18.38%), and MERIS data (R 2 = 0.552, RMSE = 34.38 mg m −3 , MARE = 43.35%).Accordingly, the 2-indices-[665] combination is proposed for MAIN-LUT.In addition, the 4-indices-[680] combination is recommended for MERIS data to overcome the limitation of other combinations, which underestimated the Chla concentration >50 mg m −3 .The spatial and temporal distribution of Chla concentration was also captured using MAIN-LUT for MERIS images.

Figure 1 .
Figure 1.Location of the sampling sites.(a) Tokyo Bay; and (b) Lake Kasumigaura.

Figure 1 .
Figure 1.Location of the sampling sites.(a) Tokyo Bay; and (b) Lake Kasumigaura.

Figure 2 .
Figure 2. Remote sensing reflectance collected from (a) Lake Kasumigaura and (b) Tokyo Bay.The red lines highlighted the 12 stations of Tokyo Bay that had additional measurements of absorption and backscattering properties along with remote sensing reflectance measurements.

Figure 2 .
Figure 2. Remote sensing reflectance collected from (a) Lake Kasumigaura and (b) Tokyo Bay.The red lines highlighted the 12 stations of Tokyo Bay that had additional measurements of absorption and backscattering properties along with remote sensing reflectance measurements.

Figure 3 .
Figure 3. Atmospheric corrected remote sensing reflectance spectra for 77 stations matched the field observation in Lake Kasumigaura from 2002 to 2012.

Figure 3 .
Figure 3. Atmospheric corrected remote sensing reflectance spectra for 77 stations matched the field observation in Lake Kasumigaura from 2002 to 2012.
j i,simu refer to simulated indices calculated for each simulated R rs spectrum (i.e., 2b simu 665 j and 3b simu 665 j for 2b-indices-[665] combination), j is an index for simulated reflectance grouped in Step (3); and n represents the number of indices in the Sample Combination (i.e., n = 2 for 2b-indices-[665] combination).The minimum RMSE value representsthe closest simulated spectrum from the measured spectrum.As a result, the tagged Chla for that simulated R rs will be the retrieved Chla concentration of the measured R rs .

Figure 4 .
Figure 4. Scatter plots of reference Chla versus algorithms indices for the 500,000 simulated reflectance spectra.(a-h) represent indices of eight algorithms mentioned in Table3.

Figure 4 .
Figure 4. Scatter plots of reference Chla versus algorithms indices for the 500,000 simulated reflectance spectra.(a-h) represent indices of eight algorithms mentioned in Table3.

Figure 5 .; 35 ,
Figure 5. Concept and schematic flowchart of the multi-algorithm indices and the look-up table (MAIN-LUT) technique for the 2-indices-[665] combination.For other combinations in Table 4, a similar concept was employed.6652b

Figure 6 .
Figure 6.Validation plots between Chla measured in Tokyo Bay (71 stations collected between 2010 and 2013) and retrieved Chla from eight different algorithms' combinations assessed for MAIN-LUT.Algorithms used in each combination from (a) to (h) were mentioned in Table 4.

Figure 6 .
Figure 6.Validation plots between Chla measured in Tokyo Bay (71 stations collected between 2010 and 2013) and retrieved Chla from eight different algorithms' combinations assessed for MAIN-LUT.Algorithms used in each combination from (a) to (h) were mentioned in Table 4.

Figure 7 .
Figure 7. Validation plots between Chla measured in Lake Kasumigaura (68 stations collected in 2016) and retrieved Chla from different algorithms' combinations assessed for MAIN-LUT.Algorithms used in each combination from (a) to (h) were mentioned in Table 4.

Figure 7 .
Figure 7. Validation plots between Chla measured in Lake Kasumigaura (68 stations collected in 2016) and retrieved Chla from different algorithms' combinations assessed for MAIN-LUT.Algorithms used in each combination from (a) to (h) were mentioned in Table 4.
Figure 7. Validation plots between Chla measured in Lake Kasumigaura (68 stations collected in 2016) and retrieved Chla from different algorithms' combinations assessed for MAIN-LUT.Algorithms used in each combination from (a) to (h) were mentioned in Table 4.
show the spatial distribution of Chla concentration for the Lake Kasumigaura using two MERIS images on July 2004 and May 2009.The Chla concentrations were relatively high (up to 90 mg m −3 ) on May 2009, compared with July 2004 (up to 50 mg m −3 ) as shown in Figure 9a,b, respectively.In order to check the accuracy of these results, the monthly measured Chla concentrations at two stations (i.e., Station 1 (St. 1) and Station 2 (St.2)) during 2002-2011 were illustrated in Figure 9c.The measured Chla concentrations were low on July 2004 (St. 1 = 32 mg m −3 , St. 2 = 13 mg m −3 ) and relatively high on May 2009 (St. 1 = 79 mg m −3 , St. 2 = 58 mg m −3 ), which emphasized the ability of the MAIN-LUT technique to capture the spatial and temporal distribution of Chla concentration.In general, the Chla concentration of St. 1 was higher than that at St. 2 due to the nutrient-rich inflow from the Sakura River and the Hanamura River into the lake near St. 1.

Figure 8 .
Figure 8. Validation plots between Chla measured in Lake Kasumigaura (77 stations matched Medium Resolution Imaging Spectrometer (MERIS) images and collected during 2002-2012) and retrieved Chla from different algorithms' combinations assessed for MAIN-LUT.The red and blue dots highlighted stations that matched MERIS images at same day and one-day differences, respectively.Algorithms used in each combination from (a) to (h) were mentioned in Table4.

Figure 8 .
Figure 8. Validation plots between Chla measured in Lake Kasumigaura (77 stations matched Medium Resolution Imaging Spectrometer (MERIS) images and collected during 2002-2012) and retrieved Chla from different algorithms' combinations assessed for MAIN-LUT.The red and blue dots highlighted stations that matched MERIS images at same day and one-day differences, respectively.Algorithms used in each combination from (a) to (h) were mentioned in Table4.

Figure 9 .
Figure 9. Spatial distribution of Chla concentration for Lake Kasumigaura was generated using MAIN-LUT for MERIS images on (a) 7 July 2004 and (b) 14 May 2009; and (c) the time series of monthly measured Chla concentration in Lake Kasumigaura at station 1 and station 2 from 2002 to 2011.The two black lines on the time series chart represent the dates of the two images.

Figure 9 .
Figure 9. Spatial distribution of Chla concentration for Lake Kasumigaura was generated using MAIN-LUT for MERIS images on (a) 7 July 2004 and (b) 14 May 2009; and (c) the time series of monthly measured Chla concentration in Lake Kasumigaura at station 1 and station 2 from 2002 to 2011.The two black lines on the time series chart represent the dates of the two images.

Table 1 .
Descriptive statistics of the water quality parameters for Tokyo Bay, and Lake Kasumigaura.Tokyo Bay with IOPs (inherent optical properties) dataset refers to stations where additional measurements of absorption and backscattering properties were performed along with remote sensing reflectance measurements.Lake Kasumigaura Datasets 1 and 2 represent measurements collected in 2016 and during 2002-2012, respectively.

Max Mean Median SD RSD (%)
* Data was collected by the National Institute for Environmental Studies (NIES).Chla refers to chlorophyll-a concentration; TSS, ISS and OSS stand for total, inorganic and organic suspended solids, respectively; a ph (440), a NAP (440) and a CDOM (440) represent the measured absorption coefficients for phytoplankton, non-algal particles (NAP), and colored dissolved organic matter (CDOM) at 440 nm bands, respectively; b b,p (442) denotes the measured backscattering at 442 nm bands; n represents the number of samples; SD stands for standard deviation; and RSD denotes relative standard deviation = (SD/Mean × 100).

Table 3 .
The incorporated band for each algorithm.

Table 4 .
Algorithms' combination investigated with multi-algorithm indices and the look-up table (MAIN-LUT) technique.

Table 4 .
Algorithms' combination investigated with multi-algorithm indices and the look-up table (MAIN-LUT) technique.

Table 5 .
Algorithms' evaluation during the calibrations stage using Medium Resolution Imaging Spectrometer (MERIS) data.
Data highlighted in bold represents the best algorithm performance with the highest R 2 during the calibration stage.Lake Kasumigaura Datasets 1 and 2 represent measurements collected in 2016 and during 2002-2012, respectively.a0, and a1 represent slope and intercept for linear regression approach.The variable "N" denotes the number of samples used during the calibration stage.

Table 5 .
Algorithms' evaluation during the calibrations stage using Medium Resolution Imaging Spectrometer (MERIS) data.

Tokyo Bay (In Situ Data) Lake Kasumigaura (Dataset 1) Lake Kasumigaura (Dataset 2)
Data highlighted in bold represents the best algorithm performance with the highest R 2 during the calibration stage.Lake Kasumigaura Datasets 1 and 2 represent collected in 2016 and during 2002-2012, respectively.a 0 , and a 1 represent slope and intercept for linear regression approach.The variable "N" denotes the number of samples used during the calibration stage.

Table 6 .
Algorithms' retrieval accuracies during the validation stage using MERIS data.

In Situ Data) Lake Kasumigaura (Dataset 1) Lake Kasumigaura (Dataset 2)
Data highlighted in bold represents the best algorithm performance with the highest R 2 .