A blended inherent optical property algorithm for global satellite ocean color observations

Water inherent optical properties (IOPs) can be derived from satellite‐measured normalized water‐leaving radiance (nLw(λ)) spectra. In this study, we evaluate the performance of the quasi‐analytical algorithm (QAA) and the near‐infrared (NIR)‐based IOP algorithm using a Hydrolight simulation data set covering a wide range of water types that span from clear open ocean to turbid coastal/inland waters. The NIR‐based algorithm produces significantly improved IOP retrievals over turbid coastal and inland waters, while the QAA algorithm performs well in the open ocean and less turbid coastal waters. Based on the advantages of the NIR‐based and QAA‐based algorithms, a combination of the NIR‐ and QAA‐based algorithm has been proposed using satellite‐measured nLw(745) as the threshold in order to produce accurate IOP products for both the open ocean and turbid coastal/inland waters. The new combined IOP algorithm can produce reasonably accurate IOP data for all water types, and can be easily implemented into the satellite ocean color data processing. The La Plata River Estuary region is used as an example to show the difference in performance of IOP retrievals from the Visible Infrared Imaging Radiometer Suite (VIIRS) measurements between 2012 and 2017 with the NIR‐based, QAA‐based, and NIR‐QAA combined IOP algorithms. We also demonstrate that the NIR‐QAA combined algorithm can be applied to VIIRS global ocean color observations to derive good quality IOP products in China's east coastal region, the US east coastal region, and the region of Mississippi River Estuary and tributaries.

Water inherent optical properties (IOPs) are the absorption and scattering of pure water, phytoplankton, color dissolved organic matter (CDOM), minerals, and so on. The normalized water-leaving radiance spectra (nL w (λ)) from the satellite ocean color observations are determined by IOPs (Gordon et al. 1988). Correspondingly, IOPs can also be inferred from the nL w (λ) spectra from satellite remote sensing. Retrievals of the water IOPs such as absorption coefficients of the phytoplankton, CDOM, and backscattering coefficients of the phytoplankton and nonalgal particles from the satellite ocean color remote sensing for both open oceans and the coastal waters are important to monitoring and assessing marine environments and natural hazards such as the biological impact of hurricanes and blooms of harmful algae. Satellite-derived water IOPs allow the study of the ocean's physical, optical, biological, and biogeochemical processes and their interactions, and permit an evaluation of the ocean's biological and biogeochemical responses to long-term global climate change.
Most of the current global satellite IOP algorithms for application in ocean color satellites such as the Sea-viewing Wide Field-of-view Sensor (SeaWiFS), the Moderate Resolution Imaging Spectroradiometer (MODIS) on the Terra and Aqua, and the Visible Infrared Imaging Radiometer Suite (VIIRS) on the Suomi national polar-orbiting partnership (SNPP) and NOAA-20 are based on the study of Gordon et al. (1988), which shows that a semi-analytical radiance model can be used to predict the upwelling spectral radiance in visible wavelengths with various materials in the water such as phytoplankton pigment concentration, dissolved organic and detrital matter, and so on. Specifically, the remote-sensing reflectance just beneath the ocean surface r rs (λ) can be written as where g 1 and g 2 are 0.0949 and 0.0794 sr −1 (Gordon et al. 1988). Coefficients g 1 and g 2 are also tuned and set to be 0.089 and 0.125 sr −1 for IOP retrievals with the quasi-analytical algorithm (QAA) (Lee et al. 2002). a t (λ) and b b (λ) are total absorption and backscattering coefficients, respectively. r rs (λ) can be computed from the remote-sensing reflectance above the *Correspondence: wei.1.shi@noaa.gov This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.
surface R rs (λ) or the normalized water-leaving radiance nL w (λ) (Morel and Gentili 1996;Gordon 2005;Wang 2006). The total absorption coefficient a t (λ) is the sum of the absorption coefficients from pure seawater a w (λ), phytoplankton a ph (λ), dissolved matter a g (λ), and nonalgal particles a d (λ). Coefficient b b (λ) is the sum of backscattering coefficients for pure seawater b bw (λ) and particulate matter b bp (λ). Specifically, the Garver-Siegel-Maritorena (GSM) IOP algorithm (Garver and Siegel 1997;Maritorena et al. 2002) fixes the absorption for dissolved and detrital matter (a dg (λ)) slope S to be 0.0206 nm −1 , and b bp (λ) power law exponential slope η to be 1.0037 in order to apply a nonlinear least-square technique to fit R rs (λ) at the four or five wavelengths of the satellite ocean color sensors. On the other hand, the QAA algorithm (Lee et al. 2002) starts with the computation of a t (λ 0 ) at a green/red wavelength using the empirical formula in Morel and Maritorena (2001). Coefficient b b (λ 0 ) is then calculated algebraically from r rs (λ) (Gordon et al. 1988;Hoge and Lyon 1996;Lee et al. 2002). The b bp (λ) power law exponential slope η is estimated from r rs (λ) at the blue and green bands with an empirical formula (Lee et al. 2002). Decomposition of a t (λ) into a ph (λ) and a dg (λ) is carried out with the empirical formula using a t (410), a t (443), and r rs (λ) at the blue and green bands (Lee et al. 2002). The generalized IOP (GIOP) algorithm (Werdell et al. 2013) is focused on the development of the GIOP model for the construction of different semi-analytical algorithms; thus, users can test different combinations of IOP models. It permits isolation and evaluation of specific modeling assumptions, construction of semi-analytical algorithms, and development of a regionally tuned IOP algorithm. The default configuration for b bp (λ) in GIOP is the same as in the QAA algorithm (Lee et al. 2002). Following the configuration of the GIOP algorithm, the retrievals of b bp (λ), a ph (λ), and a dg (λ) can be calculated via the linear or nonlinear least-square inversions of Eq. 1 using the r rs (λ) spectra.
It is noted that there are some restrictions for all three IOP algorithms. GSM is generally only valid in the open ocean with a b bp (443) range that is less than~1.0 m −1 , and an a dg (λ) that is less than~2.0 m −1 (Garver and Siegel 1997;IOCCG 2006). QAA and GIOP algorithms are only applicable to the open ocean and lessturbid coastal waters because both of these algorithms only use the r rs (λ) in the blue and green/red bands to estimate the b bp (λ 0 ) and b bp (λ) power law exponential slope η with the empirical formulas (Lee et al. 2002;Werdell et al. 2013). In turbid and highly turbid waters, studies have shown that r rs (λ) loses its sensitivity with the change of suspended particles in the visible wavelengths Wang 2009a, 2014;Shen et al. 2010a). Thus, it is necessary to retrieve ocean biological and biogeochemical properties such as chlorophyll a (Chl a) concentration, water turbidity, total suspended matter (TSM), and so on, with the satellite measurements in the near-infrared (NIR) wavelengths (Gitelson et al. 2007;Dogliotti et al. 2015;Shi et al. 2018).
Most of satellite ocean color products are derived from nL w (λ) spectra in the visible wavelengths. Unlike open oceans, coastal and inland water regions, particularly highly turbid regions like river estuaries, are dominated with inorganic materials such as clays, silts, and fine sands from the bottom or from the river flows (Bowers and Binding 2006). Suspended particles in coastal and estuarine waters are commonly in the form of aggregates or flocs (Eisma et al. 1990). For the open ocean, satellite-measured nL w (λ) spectra can be generally retrieved with uncertainty of~5% at the blue band; therefore, Chl a uncertainty is within~30% (Mcclain et al. 2004;Mcclain 2009). However, there are still some challenges in coastal regions, especially for highly turbid coastal and inland water regions. Such challenges include deriving accurate nL w (λ) spectra, i.e., atmospheric correction (Gordon and Wang 1994;Wang 2007;IOCCG 2010), Chl a concentration (O'Reilly et al. 1998;Hu et al. 2012;Wang and Son 2016), K d (490) (Wang et al. 2009), and the ocean's IOPs (Garver and Siegel 1997;Lee et al. 2002;Werdell et al. 2013). These challenges make it difficult to characterize and assess the physical, biological, and biogeochemical changes with satellite ocean color observations in coastal and inland water regions.
The normalized water-leaving radiance nL w (λ) spectra in the red, NIR, and shortwave infrared (SWIR) wavelengths can only be rarely used in the open ocean to derive ocean's biological and biogeochemical products. However, coastal and inland waters are generally featured with enhanced nL w (λ) at the NIR and SWIR wavelengths. As an example, the NIR ocean reflectance spectral shape represented by the ratio of nL w (λ) at the two NIR bands is highly dynamic and region-dependent (Doron et al. 2011;Shi and Wang 2014). The NIR spectral reflectance feature in the estuary of the Yellow River and Ancient Yellow River is found to be notably different from that of the Yangtze River (Shi and Wang 2014). This makes the satellite ocean color measurements at the red, NIR, and SWIR bands extremely valuable in conducting atmospheric correction (Gordon and Wang 1994;Wang 2007;Wang and Shi 2007), which characterizes and quantifies water properties in coastal and inland waters such as Chl a concentration (Gitelson et al. 2007), floating green algae blooms (Hu 2009;Shi and Wang 2009b), river plumes (Shi and Wang 2009c), TSM (Miller and Mckee 2004;Shen et al. 2010b;Son and Wang 2012;Shi et al. 2018), and the light diffuse attenuation coefficient at 490 nm (Wang et al. 2009;Zhang et al. 2012). It is also found that the spectral features of nL w (λ) at the red and NIR wavelengths are dynamic and regional dependent (Ruddick et al. 2006;Doron et al. 2011;Shi and Wang 2014). This can be attributed to the TSM concentrations, suspended particle shapes, size distributions, and compositions of the suspended particles (Stramski et al. 2007;Kostadinov et al. 2009;Shi and Wang 2017).
Based on the fact that, in coastal and inland waters, the absorption coefficient of seawater a w (λ) is the dominant constituent of a t (λ) at the NIR wavelengths; the semi-analytical radiance model as shown in Eq. 1 can be significantly simplified, thereby b b (λ) at the NIR wavelengths can be analytically calculated. Consequently, b bp (λ) in the visible wavelengths can be derived from b bp (λ) retrievals at the NIR wavelengths (Shi and Wang 2017). Hydrolight simulations show that b bp (λ) retrievals are accurate and robust for a wide range of both moderately and highly turbid waters in comparison to the Hydrolight input b bp (λ). In addition, Shi and Wang (2017) demonstrated that the NIR-based b bp (λ) algorithm can be applied to VIIRS-SNPP observations to produce b bp (λ) products in the global highly turbid waters such as the Bohai Sea and Yellow Sea in China's east coastal region, Amazon River Estuary, and Mississippi River Estuary and tributaries (Shi and Wang 2017).
In this study, we further extend our research to exploit VIIRSderived nL w (λ) at the visible and NIR bands using the NIR-SWIR combined atmospheric correction algorithm (Wang and Shi 2007) with VIIRS-SNPP observations between 2012 and 2017 to compute the total absorption coefficient a t (λ), color dissolved and detrital absorption coeffificent a dg (λ), and phytoplankton absorption coefficent a ph (λ) for the global ocean. A scheme is proposed to combine these NIR-based and QAA-based IOP algorithms in order to produce global IOP products from VIIRS-SNPP observations. The advantage of this scheme is first evaluated with Hydrolight simulation data. We demonstrate that this combined IOP algorithm can significantly improve IOP products in coastal and inland waters, and produce high-quality IOP products for both clear open oceans and turbid coastal/inland waters.

Methods
The data set of IOP and nL w (λ) spectra from Hydrolight simulations Due to the insufficient amount of in situ IOP data and nL w (λ) measurements covering a wide range of IOPs for both the open ocean and turbid coastal/inland water regions, Hydrolight simulations (Mobley et al. 1993;Mobley and Sundman 2013) were conducted to create the data set of IOPs and nL w (λ) spectra for a variety of water types. The Hydrolight radiative transfer model computes radiance distributions and related quantities (e.g., irradiances, reflectance, diffuse attenuation functions, etc.) in various water types using the absorption and backscattering of various constituents in the water, for given bottom depth, IOP property profiles, and atmosphere condition (Mobley et al. 1993).
To generate a representative IOPs and nL w (λ) spectra data set for the global ocean, IOPs from the International Ocean Colour Coordinating Group (IOCCG) synthesized data set for IOP algorithm development in IOCCG Report 5 (IOCCG 2006) were used as the baseline IOPs for the Hydrolight simulations. The IOPs and corresponding apparent optical properties (AOPs) can generally be used to develop, test, evaluate, and compare the IOP and AOP algorithms for satellite ocean color remote sensing. This data set contains 500 records of respective phytoplankton absorptions a ph (λ), dissolved matter absorptions a g (λ), nonalgal detrital absorptions a d (λ), and particle backscattering coefficients b bp (λ) at the wavelength interval of 10 nm between 400 and 800 nm in addition to the corresponding R rs (λ) spectra. However, the data set only covers the water types in the open ocean and less turbid coastal waters.
In order to develop an IOP algorithm for all water types, i.e., clear open ocean waters, less turbid and turbid coastal waters, estuarine, and inland waters, we further expanded the IOCCG simulated data set. Since turbid coastal waters feature enhanced particle backscattering coefficients b bp (λ), phytoplankton absorption coefficients a ph (λ), and dissolved and detrital matter absorptions a dg (λ), the last 25 records of IOPs in the IOCCG data set are chosen first. These 25 records have the highest values in b bp (λ), a d (λ), a g (λ), and a ph (λ) in the data set. With the assumption that nonalgal detrital backscattering b bd (λ) and a d (λ) are proportional to the detritus/mineral concentrations, and a g (λ) and a ph (λ) are proportional to the concentrations of the corresponding constituents, we boost detritus/mineral concentrations, a g (λ) and a ph (λ), by factors of 1, 2, 5, and 10 times. Thus, the total number of the possible IOP combinations can be as high 25 × 4 × 4 × 4 = 1600. Of the possible records of IOP data set, we randomly select a certain number of the records as the Hydrolight inputs to generate nL w (λ) spectra for the coastal and inland waters. These enhanced b bp (λ), a dg (λ), and a ph (λ) as well as the corresponding nL w (λ) spectra are randomly selected as the new IOP and nL w (λ) records for turbid coastal and inland waters.
In addition to the new IOP data records, we also ran Hydrolight simulations for the 500 records of IOPs in the synthesized data set from IOCCG Report 5 in order to produce nL w (λ) spectra at the VIIRS-SNPP bands of 410, 443, 486, 551, 671, 745, and 862 nm because the NIR nL w (862) is not covered in the IOCCG synthesized data set (IOCCG 2006). It has been verified that nL w (λ) spectra from the Hydrolight simulations are consistent with those in the IOCCG synthesized data set.
Hydrolight simulations can be used to build realistic records of enhanced IOPs and the corresponding nL w (λ) spectra in turbid waters with this scheme. Thus, the 500 records of IOPs and Hydrolight-simulated nL w (λ) spectra at the VIIRS wavelengths are combined with those of the enhanced b bp (λ), a dg (λ), and a ph (λ) along with the corresponding Hydrolight-simulated nL w (λ) spectra at the VIIRS bands, building a comprehensive data record of IOPs and nL w (λ) spectra. The new data set covers a variety of water types such as clear open ocean, less-turbid, turbid, and highly turbid coastal and inland waters. Therefore, this new data set can be used to evaluate the performance of various IOP algorithms, and develop a new IOP algorithm for producing good quality IOP products from VIIRS observations. VIIRS-derived nL w (λ) spectra from the NIR-SWIR based ocean color data processing VIIRS is one of the major sensors onboard the SNPP satellite, which was launched on 28 October 2011. It has 22 spectral bands, including 14 reflective solar bands, 7 thermal emissive bands, and a panchromatic day/night band. The spectral band specifications are similar to MODIS in order to observe Earth's atmosphere, land, and ocean properties (Goldberg et al. 2013). Ocean color environmental data records are one of the key product suites derived from VIIRS (Wang et al. 2013. VIIRS-SNPP has five visible bands (M1-M5) with nominal central wavelengths at 410, 443, 486, 551, and 671 nm, two NIR bands (M6 and M7) at wavelengths of 745 and 862 nm, and three SWIR bands (M8, M10, and M11) at wavelengths of 1238, 1601, and 2257 nm for satellite ocean color data processing. The spatial resolution of these spectral bands is 750 m. In addition, there are imaging bands (I-bands) with spatial resolution of 375 m. In particular, VIIRS I1 band at 638 nm is useful for various coastal and inland water applications (Wang and Jiang 2018). For satellite ocean color data processing, we also carried out the on-orbit vicarious calibration for VIIRS using the in situ nL w (λ) spectra from the Marine Optical Buoy (MOBY) (Clark et al. 1997) located in the waters off Hawaii .
To derive nL w (λ) spectra from satellite observations, one fundamental assumption for atmospheric correction is the black ocean assumption at the NIR and SWIR wavelengths for the specific water types. For the open ocean, atmospheric correction can be carried out using the two VIIRS NIR bands (745 and 862 nm) to derive nL w (λ) spectra (Gordon and Wang 1994), and consequently, ocean biological and biogeochemical property data can be obtained using nL w (λ) spectra. For turbid coastal and inland waters, however, nL w (λ) spectra from the visible to NIR wavelengths can be derived using the SWIR-based atmospheric correction algorithm (Wang 2007;Wang and Shi 2007) because the black ocean assumption is generally true for the SWIR bands (Shi and Wang 2009a), e.g., at the VIIRS SWIR bands of 1238, 1601, and 2257 nm.
VIIRS ocean color products are produced with the official NOAA VIIRS ocean color data processing system, i.e., the multisensor level-1 to level-2 (MSL12) (Wang et al. 2013. A method of ocean color data processing using the combined NIR and SWIR bands for atmospheric correction for satellite ocean color data processing (Wang and Shi 2007) was implemented in the MSL12 ocean color data processing system. With the NIR-SWIR combined atmospheric correction algorithm, ocean color data can be processed using the standard (NIR) atmospheric correction algorithm for the open ocean (Gordon and Wang 1994), whereas for the turbid coastal and inland waters, the SWIR atmospheric correction algorithm (Wang 2007) can be executed. The combined NIR-SWIR algorithm has been used to routinely produce global nL w (λ) spectra for VIIRS satellite ocean color data processing. In this study, nL w (λ) spectra derived from the NIR-SWIR atmospheric correction algorithm are used to evaluate and demonstrate IOP products with different IOP algorithms for a variety of ocean regions such as China's east coastal region, the La Plata River Estuary, and the Mississippi River Estuary region.

The QAA-based and NIR-based IOP algorithms
The QAA IOP algorithm was developed to derive the IOPs for optically deep waters (Lee et al. 2002). Briefly, the QAA algorithm includes two steps. The first step is to derive the total absorption coefficients a t (λ) and particle backscattering coefficients b bp (λ). An empirical formula is used to compute a t (λ) at the green/red wavelengths from R rs (λ) in the visible bands (Morel and Maritorena 2001). It is then applied to the water-leaving reflectance model to calculate b bp (λ 0 ) at the green/red wavelengths. The power law exponential slope η is also empirically estimated from R rs (λ) in the visible bands. With b bp (λ 0 ) at the green/red wavelengths and the power law exponential slopes η, b bp (λ) spectra are consequently derived, and they are applied to the remotesensing reflectance model to derive a t (λ) spectra.
The purpose of the second step in the QAA algorithm is to decompose a t (λ) into the dissolved and detrital absorption coefficient, a dg (λ), and the phytoplankton absorption coefficient, a ph (λ). In this step, the spectral ratio of a ph (410)/ a ph (443) is first estimated, and next the spectral slope S of a dg (λ) is calculated. These two parameters are both based on the empirical formulas from the in situ observations. With the known a ph (410)/a ph (443), S, a t (410), and a t (443) from step 1, a dg (443) and a ph (443) can be computed accordingly. Using a dg (443) and spectral slope S of a dg (λ), a dg (λ) spectra can be derived. Thus, a ph (λ) is consequently computed after deducting a dg (λ) and a w (λ) from a t (λ).
In comparison to the QAA IOP algorithm, Wang (2014, 2017) show that b bp (λ) can be computed analytically in turbid coastal and inland waters with the NIR-based b bp (λ) algorithm from VIIRS-SNPP observations. Briefly, the pure seawater absorption coefficients a w (λ) at the NIR bands are typically significantly higher than those from other constituents, i.e., a w (λ) >> a ph (λ), a g (λ), and a d (λ). For example, a w (862) is 5 m −1 , while a ph (λ), a g (λ), and a d (λ) at the NIR 862 nm are normally negligible compared to a w (862). Therefore, with the spectral features of a w (λ), a ph (λ), a g (λ), and a d (λ) at the VIIRS NIR 745 and 862 nm bands, b b (λ)/(a(λ) + b b (λ)) can be approximated as The NIR-based b bp (λ) algorithm is robust, and can be used to derive b bp (λ) for turbid waters when nL w (745) and nL w (862) are less than~6 mW cm −2 μm −1 sr −1 and~4 mW cm −2 μm −1 sr −1 (Shi and Wang 2017), respectively. After deriving b bp (λ) in the turbid waters, a t (λ) spectra can be computed with the IOPreflectance model as shown in Eq. 1. Following the same approach in the QAA algorithm, a t (λ) can be further decomposed into the dissolved and detrital absorption coefficient a dg (λ) and phytoplankton absorption coefficient a ph (λ) in turbid waters. Table 1 summarizes differences between the NIR-based and QAAbased IOP algorithms.
In this study, we evaluate and compare the performance of the NIR-based and QAA-based IOP retrievals using the new Hydrolight-generated comprehensive data set of IOPs and the corresponding nL w (λ) spectra. Since this data set covers various water types from clear open ocean waters, coastal less-turbid waters, and coastal, estuarine, and inland turbid and highly turbid waters, this evaluation can help to understand the performance of these two IOP algorithms when they are used to produce IOP products from VIIRS global ocean color observations.
In addition, a previous study revealed the advantages and disadvantages of the QAA-based and NIR-based b bp (λ) algorithms (Shi and Wang 2017). These two b bp (λ) algorithms actually can be complementary with each other. Thus, a combined NIR-and QAA-based b bp (λ) algorithm was proposed and demonstrated an improvement of global b bp (λ) for both open oceans and coastal turbid waters (Shi and Wang 2017). In this study, we further explore the possibility for a combined algorithm for all IOP products including a ph (λ), a dg (λ), and a t (λ) in order to produce high-quality IOP products from VIIRS global ocean observations. Specifically, the IOP products derived from the combined NIR-and QAA-based algorithm will be demonstrated in various regions.

Development of the combined IOP algorithm with the Hydrolight simulation data set
The new comprehensive Hydrolight simulation data set is used to evaluate the performance of the NIR-based, QAA-based, and NIR-QAA combined IOP algorithms. The total number of records for IOPs and nL w (λ) is 664. Figure 1 shows the comparisons of the NIR-derived b bp (λ) vs. known b bp (λ) (Fig. 1a), the NIRderived a ph (λ) vs. known a ph (λ) (Fig. 1b), the NIR-derived a dg (λ) vs. known a dg (λ) (Fig. 1c), and the NIR-derived a t (λ) vs. known a t (λ) (Fig. 1d). High uncertainties of the NIR-based b bp (λ) (Fig. 1a) and a t (λ) (Fig. 1d) occur for the clear open ocean waters with known b bp (λ) < 0.02-0.03 m −1 and a t (λ) < 0.1 m −1 , respectively. This is attributed to the extremely low nL w (λ) at the VIIRS NIR bands for the open ocean waters. On the other hand, the NIRbased IOP algorithm for b bp (λ) and a t (λ) shows high level of accuracy for high b bp (λ) and a t (λ) values as shown in Fig. 1a,d. Figure 1b,c shows that the uncertainties of a dg (λ) and a ph (λ) are high for low a dg (λ) and a ph (λ), while the performance of the NIR-based a dg (λ) and a ph (λ) in coastal turbid waters with high a dg (λ) and a ph (λ) values is reasonably good. This demonstrates that the NIR-based IOP algorithm for b bp (λ), a t (λ), a dg (λ), and a ph (λ) can be used to derive IOPs for coastal turbid waters with good accuracy, but it may lead to significant biases and uncertainties in clear open ocean waters. Figure 2 shows the performance of the b bp (λ) (Fig. 2a), a ph (λ) (Fig. 2b), a dg (λ) (Fig. 2c), and a t (λ) (Fig. 2d) derived from nL w (λ) spectra with the QAA IOP algorithm in comparison to the Hydrolight IOP inputs. Figure 2a shows that the QAAbased b bp (λ) performs well for b bp (λ) < 0.1 m −1 . However, the QAA-based b bp (λ) is significantly underestimated for b bp (λ) approximately greater than 0.2-0.3 m −1 . The uncertainty results in Figs. 1, 2 are consistent with the findings in Lee et al. (2010). Figure 1 indeed shows that uncertainties of a ph (λ) are larger than those of a dg (λ) when a ph (λ) and a dg (λ) are similar, even though a t (λ) retrievals agree well with a t (λ) input in the range of a t (λ) approximately greater than 0.5 m −1 . Figure 2 also shows that a t (λ) retrievals in the blue-green wavelengths are generally within~10% for 500 records of the IOCCG data set. It shows that uncertainties of a ph (λ) are larger than those of a dg (λ) for 500 records of the IOCCG data set.
Similarly, a t (λ) retrievals are generally consistent with the input a t (λ) < 1.0 m −1 (Fig. 2d), but significant underestimations can be found with increasing a t (λ). It performs better than a ph (λ) and a dg (λ), but notable uncertainties and biases still exist especially in the range of a t (λ) approximately greater than 1.0 m −1 . The QAA a ph (λ) retrievals (Fig. 2b) show that the QAA algorithm performed better with low a ph (λ) than those with high a ph (λ) (Fig. 2b). In fact, Fig. 2 shows that the QAA-derived a ph (λ) tends to be underestimated when a ph (λ) is greater than 0.5 m −1 . Similarly, a dg (λ) retrievals are also underestimated when a dg (λ) greater than 3 m −1 (Fig. 2c). In fact, the simulation performance evaluation for the QAA IOP algorithm shows that the QAAbased IOP algorithm is not suitable for turbid coastal waters, which normally feature high values of b bp (λ), a t (λ), a dg (λ), and a ph (λ).
Some of the different and complementary performances of the QAA-based and NIR-based IOP algorithms for various water types as shown in Figs. 1, 2 suggest that a combined IOP algorithm is necessary in order to produce high-quality IOP products for both open oceans and turbid coastal/inland waters. In fact, it has been demonstrated that the approach of the NIR-SWIR atmospheric correction algorithm can combine the NIR-and SWIR-based atmospheric correction algorithms in order to derive high-quality nL w (λ) spectra for the open   (Wang and Shi 2007). Similarly, a method was proposed for the combined NIR-and QAA-based b bp (λ) algorithm using nL w (745) as the criteria to determine which algorithm should be used for producing VIIRS b bp (λ) data (Shi and Wang 2017). Following the same strategy, a combined scheme is developed and proposed to derive all IOPs, i.e., b bp (λ), a ph (λ), a dg (λ), and a t (λ) in order to produce accurate IOP products for both open ocean and turbid coastal/inland waters. This scheme can be easily implemented into the MSL12 ocean color data processing system. With extensive tests and evaluations, the specific procedure for the NIR-QAA IOP scheme is outlined below.
With the proposed method, optimal IOP retrievals for b bp (λ), a ph (λ), a dg (λ), and a t (λ) can be achieved. Note that IOP retrievals from satellite observations are prone to various uncertainties and noise. For example, it is possible that IOP retrievals from one IOP algorithm can fail when deriving negative values. In this case, retrievals from another IOP algorithm are used for deriving the corresponding IOP values.
The performance of the NIR-QAA combined IOP algorithm for b bp (λ) (Fig. 3a), a ph (λ) (Fig. 3b), a dg (λ) (Fig. 3c), and a t (λ) (Fig. 3d) with the comprehensive IOP data set is significantly improved. The NIR-QAA combined IOP algorithm can derive reasonably high-quality IOPs for all water types. In comparison to the QAA-based algorithm (Fig. 2), all four IOP parameters show Furthermore, the statistical evaluation shows that the good retrievals for all four IOP parameters can be achieved with the NIR-QAA combined IOP algorithm (Table 2). High-quality retrievals can be achieved for the IOP products of b bp (λ) and a t (λ) with the corresponding coefficients of determination R 2 at 0.952 and 0.936, respectively. Reasonably accurate a ph (λ) and a dg (λ) can be achieved with the NIR-QAA combined IOP algorithm. The R 2 values for the a ph (λ) and a dg (λ) are 0.723 and 0.738, respectively. This shows that retrievals of a ph (λ) and a dg (λ) with the NIR-QAA IOP algorithm are less accurate than those of b bp (λ) and a t (λ). Since the Hydrolight IOP data set covers a wide range of water types, the NIR-QAA combined IOP algorithm can be applied to the global satellite ocean color observations for deriving the IOP parameters for both open ocean and turbid coastal/inland waters.

VIIRS-derived IOP products over the La Plata River Estuary
The La Plata River Estuary between Uruguay and Argentina is one of the most turbid regions in the world (Shi and Wang 2010). To further evaluate the performance of the IOP products from the QAA-based, NIR-based, and NIR-QAA combined IOP algorithms, the La Plata River Estuary region is selected to produce the IOP products from the three algorithms with VIIRS observations. The different b bp (λ) products from the NIR-based IOP algorithm (Fig. 4a-e), the QAA-based IOP algorithm (Fig. 4f-h), and the NIR-QAA combined IOP algorithm ( Fig. 4i-m) show that high-quality IOP products can be obtained with the NIR-QAA combined IOP algorithm. In general, QAA-derived b bp (λ) products are smooth and less noisy for the entire region. Note that there are no b bp (745) and b bp (862) data from the QAA-based algorithm. Figure 4a-e shows the NIR-based climatology b bp (λ) in this region. Indeed, b bp (λ) images with the NIR-based approach in the open ocean show high data noise especially for b bp (443) (Fig. 4a). Actually, b bp (443) (Fig. 4a), b bp (551) (Fig. 4b), and b bp (671) (Fig. 4c) are all derived from b bp (745) (Fig. 4d) and b bp (862) (Fig. 4e) with the NIR-based algorithm. Thus, small error and noise in nL w (745) and nL w (862) due to atmospheric correction can lead to significant uncertainties of b bp (λ) in visible bands in the offshore waters where nL w (745) and nL w (862) are close to 0. It is also noted that the NIR-based b bp (λ) in  highly turbid waters (Fig. 4a-f) in the visible bands is~4-5 times higher than those from the corresponding QAA-based b bp (λ) products ( Fig. 4f-h).
The QAA-and NIR-based b bp (λ) results in Fig. 4 clearly demonstrate the advantages and disadvantages of each IOP algorithm. The QAA-based algorithm usually underestimates b bp (λ) and even loses its sensitivity to b bp (λ) change in moderately to highly turbid waters (Fig. 4f-h), while the NIR-based b bp (λ) retrievals show high data noise in the offshore waters in comparison with those from the QAA-based algorithm. This is consistent with the performance assessment using the comprehensive Hydrolight simulation data set as shown in Figs. 1, 2. Figure 4i-m shows results of the NIR-QAA combined VIIRS climatology b bp (λ) at the bands 443 nm b bp (443) (Fig. 4i), 551 nm b bp (551) (Fig. 4j), and 671 nm b bp (671) (Fig. 4k), respectively. Note that the QAA IOP algorithm is specifically developed for IOP retrievals in the visible bands and not at the NIR bands. Since there are no QAA b bp (745) and b bp (862) data, b bp (745) (Fig. 4l) and b bp (862) (Fig. 4m) only show the areas with nL w (745) ≥ 0.2 mW cm −2 μm −1 sr −1 from the NIR-based IOP algorithm. In comparison with the NIR-based climatology b bp (λ), the NIR-QAA combined b bp (λ) products show significant improvement in the offshore region. In the estuarine regions, however, values of the NIR-QAA combined b bp (λ) become much higher than those from the QAA approach. Furthermore, no obvious data discontinuities can be found in the transition zones from the QAA b bp (λ) in the offshore open ocean to the NIR-based b bp (λ) in turbid coastal waters. This demonstrates that the proposed approach can provide the best b bp (λ) data quality for both open ocean and turbid coastal/inland waters.
( Fig. 6a-c), QAA-based IOP algorithm ( Fig. 6d-f), and NIR-QAA combined IOP algorithm ( Fig. 6g-i). Similar to Figs. 4a-c, 5a-c, the QAA-derived a dg (λ) (Fig. 6d-f) has much less noise than those from the NIR-based algorithm (Fig. 6a-c) in the open ocean. However, the QAA-based a dg (λ) values in the La Plata River Estuary are much less than those of the NIR-based a dg (λ). The QAAbased and NIR-based a dg (λ) in the La Plata River Estuary are consistent with the Hydrolight simulations in Figs. 1, 2, which show that the QAA-based a dg (λ) significantly underestimates a dg (λ), while the NIR-based a dg (λ) is close to the sea truth in turbid waters. On the other hand, high noise in the NIR-based a dg (λ) shown in Fig. 6a-c for the offshore waters is the retrieval artifact due to the insignificant nL w (λ) at the VIIRS NIR bands as suggested by the Hydrolight simulations (Fig. 1). In contrast, the QAA-based a dg (λ) in the offshore waters is reliable and accurate. The NIR-QAA combined a dg (λ) images in this region (Fig. 6g-i) show the high quality of a dg (λ) for both the turbid estuarine waters and the offshore waters. In fact, there are no a dg (λ) data gaps from the NIR-QAA combined algorithm. Figure 7 shows the VIIRS-SNPP climatology a t (443), a t (551), and a t (671) products derived from the NIR-based IOP algorithm (Fig. 7a-c), QAA-based IOP algorithm (Fig. 7d-f), and NIR-QAA combined IOP algorithm (Fig. 7g-i), respectively. Note that a t (λ) = a ph (λ) + a dg (λ) + a w (λ), and a w (λ) can play a significant role in the total absorption coefficient a t (λ) at the VIIRS 551 and 671 nm bands especially for the open ocean. The NIR-based a t (λ) indeed shows high noise a t (λ) retrievals in the offshore clear waters in Fig. 7a-c. The Hydrolight simulations for the NIR-based a t (λ) retrievals in Fig. 1d also suggest that a t (λ) derived from the NIR-based IOP algorithm in the offshore open ocean can be erroneous. In comparison, the QAA-based a t (λ) in the open ocean is good data quality in terms of both the noise level and data accuracy (Fig. 2d).
The QAA-based a t (λ) shows low noise for the coastal and estuarine turbid regions as well as for the open ocean. However, the QAA-based a t (λ) (Fig. 7d-f) in the La Plata River Estuary is significantly lower than the corresponding NIR-based a t (λ) (Fig. 7a-c). As an example, the QAA-based a t (551) (Fig. 7e) in the estuary is generally~0.5 m −1 or less, while the NIR-based a t (551) (Fig. 7b) in the estuary is~2.5 m −1 . This is also consistent with the Hydrolight simulations for turbid waters shown in Figs. 1d, 2d. Figure 7g-i shows IOP results from the NIR-QAA combined a t (λ) in the La Plata River Estuary. With the combination of the QAA-based a t (λ) for offshore oceans and NIR-based a t (λ) for coastal and estuarine waters, good quality a t (λ) data for the entire a ph (443) a ph (551) a ph (671) a ph (λ) (m -1 ) Fig. 5. VIIRS-measured climatology a ph (λ) at various VIIRS bands in the La Plata River Estuary derived using the approaches of (a-c) the NIR-based IOP algorithm, (d-f) the QAA-based IOP algorithm, and (g-i) the NIR-QAA combined IOP algorithm. (671) a dg (λ) (m -1 ) Fig. 6. VIIRS-measured climatology a dg (λ) at various VIIRS bands in the La Plata River Estuary derived using the approaches of (a-c) the NIR-based IOP algorithm, (d-f) the QAA-based IOP algorithm, and (g-i) the NIR-QAA combined IOP algorithm.
region are obtained. Specifically, combined a t (λ) keeps highquality a t (λ) from the NIR-based IOP algorithm in the La Plata River Estuary and high-quality a t (λ) from the QAA-based IOP algorithm in the offshore region. In the La Plata River Estuary, a t (443), a t (551), and a t (671) can reach~5 m −1 , 2.5 m −1 , and 1.5 m −1 , respectively. Indeed, with the NIR-QAA combined IOP algorithm, the less noisy and more accurate a t (λ) in offshore open ocean waters from the QAA IOP algorithm is chosen to produce high-quality a t (λ) for the entire region. A few studies have been conducted in the La Plata River Estuary region. The TSM concentration and the turbidity from the satellite ocean color observations and the in situ measurements (Moreira et al. 2013;Dogliotti et al. 2016) indeed show significant enhancements in the inner and middle La Plata River Estuary due to the river discharge and wind-driven sediment resuspension. The spatial pattern of the TSM concentration and turbidity in these studies are consistent with the b bp (λ) pattern in Fig. 4i-m and a dg (λ) pattern in Fig. 6g-i. Since both b bp (λ) and a dg (λ) are proportional to the particle concentrations in the water column, the consistence between the IOP products in this study and the measurements of TSM and turbidity qualitatively suggests that the IOP products in the La Plata River Estuary from the NIR-QAA combined method as shown in Figs. 4-7 are reasonably accurate.
The NIR-QAA combined IOP products in the other coastal regions VIIRS-SNPP IOP products b bp (λ), a ph (λ), a dg (λ), and a t (λ) derived from the QAA-based algorithm, NIR-based algorithm, and NIR-QAA combined algorithm in the La Plata Estuary region have different performances as shown in  Results show that the NIR-QAA blended IOP algorithm can produce high-quality b bp (λ), a ph (λ), a dg (λ), and a t (λ) for both offshore open ocean waters and turbid coastal waters. To further evaluate the applicability of the NIR-QAA combined IOP algorithm for the global ocean, China's east coastal region, the US east coastal region, and the Mississippi River Estuary region are selected to evaluate IOPs of b bp (λ), a t (λ), a dg (λ), and a ph (λ) using the NIR-QAA combined IOP algorithm.
The climatology b bp (λ), a ph (λ), a dg (λ), and a t (λ) from the combined IOP algorithm in China's east coastal region at the VIIRS-SNPP bands of 443, 551, and 671 nm are also derived (Fig. 8).
China's east coastal regions feature a high loading of sediment concentration Wang 2010, 2012) such as in the regions of the Hangzhou Bay, Yangtze River Estuary, and Subei Shoal in the continental shelf (Milliman and Meade 1983;Milliman et al. 1985). In fact, TSM concentrations can reach over 1000 g m −3 in some regions (Shen et al. 2010b;Zhang et al. 2010a). Farther offshore, the ocean is featured with clear open ocean waters with Kuroshio Current flowing from the southwest direction to the northeast direction. The CDOM absorption coefficient at the wavelength of 400 nm a CDOM (400) ranges approximately from greater than 1 m −1 in the coastal region to near 0 m −1 in the open ocean (Gong 2004).
For all the products of b bp (λ), a t (λ), a dg (λ), and a ph (λ), the noise level is low for both coastal and offshore waters. b bp (λ) reaches over~2 m −1 in highly turbid waters. In the Hangzhou Bay, Subei Shoal, and Yangtze River Estuary enhanced a dg (443) (Fig. 8c) and a t (443) (Fig. 8d) can be over~4 m −1 , while a ph (443) (Fig. 8b) only shows moderate increase in comparison to a ph (443) in the open ocean. This further suggests that high a t (λ) in coastal waters is attributed to the enhanced absorption of the dissolved and detrital matter, which are driven by the river runoff and sediment resuspensions in the Hangzhou Bay, Yangtze River Estuary, and Subei Shoal. Phytoplankton bloom is not a major cause for the features of enhanced a t (λ) along China's east coastal region. It is also noted that a t (671) in the open ocean is significantly higher than the a t (443) and a t (551) because a w (671) is much larger than a w (443) and a w (551). In the turbid China's east coastal region, b bp (λ) values generally become flat spectrally or slightly increase from b bp (443) (Fig. 8a) to b bp (671) (Fig. 8i). This implies that the power law exponential slope η of b bp (λ) is slightly negative or close to 0 (Shi and Wang 2014, 2017, 2019, and further suggests that the particle size in the water column of the highly turbid waters is generally larger than that of the other coastal regions. Similarly, Fig. 9 shows the climatology b bp (λ), a ph (λ), a dg (λ), and a t (λ) in the Mississippi River Estuary and tributaries from VIIRS-SNPP observations using the NIR-QAA combined IOP algorithm. Even though the IOP features in this region are not as pronounced as in China's east coastal region, enhanced a dg (443) (Fig. 9c) and a t (443) (Fig. 9d) can still be observed in the Mississippi River Estuary, Lake Pontchartrain, and Atchafalaya River Estuary. Similar to China's east coastal region, a dg (λ) (Fig. 9c,g,k) is the dominant component of a t (λ) (Fig. 9d,h,l). However, a ph (443) in this region (Fig. 9b) is generally higher than a ph (443) in China's east coastal region (Fig. 8b). Thus, its contribution to a t (443) in this region is not negligible. Higher a ph (λ) in this region also suggests that coastal and inland waters around the Mississippi River Estuary and its tributaries are more productive than China's east coastal region and La Plata River Estuary, even though the waters are less turbid than those two regions in terms of b bp (λ) and nL w (λ) magnitudes Wang 2010, 2017). Different from the b bp (λ) spectra in turbid waters along China's east coastal region, b bp (λ) values generally decrease with the increase of wavelength, e.g., from b bp (443) (Fig. 9a) to b bp (671) (Fig. 9i). This shows that the power law exponential slope η for b bp (λ) is positive in the Mississippi River Estuary, Lake Pontchartrain, and Atchafalaya River Estuary, while it is negative in China's east coastal region (Shi and Wang 2019). Since power law exponential slope η of b bp (λ) is related to the particle size distribution slope ξ in the water column (Kostadinov et al. 2009), the polynomial relationship between η and ξ implies that the particle size in the Mississippi River Estuary and tributaries is smaller than that in China's east coastal region (Shi and Wang 2019). Finally, climatology b bp (λ), a ph (λ), a dg (λ), and a t (λ) derived from the NIR-QAA combined IOP algorithm in the US east coastal region are shown in Fig. 10. Similar to the other regions discussed previously, these IOP products are in good quality with little noise for both the open ocean and turbid coastal waters. In comparison to the La Plata River Estuary region (Figs. 4-7), China's east coastal region (Fig. 8), and Mississippi River Estuary and tributaries (Fig. 9), the US east coastal region is less turbid with moderately high b bp (λ) in the Chesapeake Bay, Delaware Bay, and Pamlico Sound. Highly turbid waters at the tips of the Chesapeake Bay, Delaware Bay, and Albemarle Sound with b bp (443) over~2 m −1 can also be observed.
In the Chesapeake Bay and Delaware Bay, a dg (443) and a ph (443) are in the same order. This is different from the other highly turbid regions such as the La Plata River Estuary, China's , a ph (λ), a dg (λ), a t (λ) (m -1 ) Fig. 8. VIIRS-measured climatology b bp (λ), a ph (λ), a dg (λ), and a t (λ) in China's east coastal region derived using the NIR-QAA combined IOP algorithm for VIIRS bands at (a-d) 443 nm, (e-h) 551 nm, and (i-l) 671 nm. east coast region, and Mississippi River Estuary where a dg (443) is the dominant component in a t (443). In the US east coastal region, a dg (443) (Fig. 10c) is also significantly higher than a ph (443) (Fig. 10b) in the Pamlico Sound. This indicates that river runoff and sediment resuspension is the major process that drives a dg (443) in the Pamlico Sound, while high productivity and phytoplankton bloom in the Chesapeake Bay and Delaware Bay can lead to elevated a ph (443). Climatology a t (λ) images (Fig. 10d, h,l) show that a t (443) reaches over~3 m −1 in these turbid areas. Broadly enhanced a t (λ) can also be found in most parts of the Chesapeake Bay, Delaware Bay, and Pamlico Sound regions.
The IOP products in China's east coastal region, the Mississippi River Estuary, and the Chesapeake Bay as shown in Figs. 8-10 are consistent with the observations from various studies in these three regions (D'Sa et al. 2007;Werdell et al. 2009;Zhang et al. 2010bZhang et al. , 2012. In China's east coastal region, the IOP values in Fig. 8 qualitatively agree with the in situ measurements in the region (Zhang et al. 2010b) and also specifically in Lake Taihu (Zhang et al. 2012). In the Mississippi River Estuary, b bp (λ) slope from the blue to red bands as indicated in Fig. 9 is similar to the observations (D'Sa et al. 2007), and high a dg (λ) and b bp (λ) values are found in the Atchafalaya River Estuary since the sediment plume is trapped within the coastal current (Falcini et al. 2012). In the Chesapeake Bay, a ph (λ) trends lower from the upper Chesapeake Bay to the lower Chesapeake Bay as shown in Fig. 10. This agrees with the observations of the spatial variation of Chl a in the region (Werdell et al. 2009). On the other hand, Tzortziou et al. (2006) show that measured and model-estimated water-leaving radiances agree with each other very well in the Chesapeake Bay. This further suggests that the IOP retrievals from the bio-optical models as shown in this study are expected to match well with the sea truth.
Discussion nL w (λ) spectra are determined by the constituents of seawaters such as CDOM, phytoplankton particles, inorganic mineral sediments, and so on. The specific spectral features of each constituent such as absorption coefficients, e.g., a ph (λ), a dg (λ), and a t (λ), and backscattering coefficients, e.g., b bp (λ), determine the spectral variability of nL w (λ) or remote-sensing reflectance just beneath the surface r rs (λ).
Coastal waters with a diffuse attenuation coefficient K d (490) > 1 m −1 in the global ocean can account for about 1/7 of the global continental shelf (Shi and Wang 2010). Considering that coastal waters are often featured with complexity of the seawater constituents as well as the dynamic variability, it is critical to derive accurate IOP parameters in coastal waters in order to study climatic, hydrological, physical, biological, and biogeochemical processes and their interactions with each other. From this point of view, the new NIR-QAA IOP algorithm can broadly improve the global satellite IOP products in the global continental shelf region instead of just some highly turbid estuarine regions.
To derive all IOP parameters in the entire spectrum from the blue to the NIR bands in coastal and inland waters, it is a a t (λ) a dg (λ) a ph (λ) b bp (λ) b bp (λ), a ph (λ), a dg (λ), a t (λ) (m -1 ) Fig. 9. VIIRS-measured climatology b bp (λ), a ph (λ), a dg (λ), and a t (λ) in the Mississippi River Estuary and tributaries derived using the NIR-QAA combined IOP algorithm for VIIRS bands at (a-d) 443 nm, (e-h) 551 nm, and (i-l) 671 nm.
prerequisite to compute b bp (λ) first. Even though the power law of b bp (λ) has been widely used in the IOP modeling effort, it has been typically a challenge to accurately estimate reference b bp (λ 0 ) and b bp (λ) power law exponential slope η directly with the ocean reflectance in the visible bands. As shown previously, b bp (λ 0 ) and η in the QAA algorithm are computed using ocean reflectance in the visible bands empirically (Lee et al. 2002;IOCCG 2006), while η is set to a constant value in the GSM IOP algorithm (Garver and Siegel 1997;Maritorena et al. 2002). Sensitivity studies have shown that η can significantly affect the accuracy of IOP retrievals (Hoge and Lyon 1996). In addition, satellite-derived nL w (λ) spectra in the visible bands become less-related or even totally unrelated to the seawater constituent changes when nL w (λ) at the NIR wavelength is over a certain value Wang 2009a, 2014;Shen et al. 2010a). This suggests that IOP algorithms such as QAA and GSM using nL w (λ) at the visible bands have some significant limitations in deriving accurate water IOPs in turbid coastal and inland waters.
In comparison, the NIR-based IOP algorithm for turbid coastal and inland waters in this study is essentially analytical since the spectral features of the seawater constituents in the NIR wavelengths are used in this algorithm. Since it can compute the η values directly, this algorithm can effectively address the large variation of η values (Reynolds et al. 2001;Wozniak and Stramski 2004;Shi and Wang 2017). Otherwise, the η variation can lead to notable errors with empirical or semi-analytical IOP approaches.
A recent study (Shi and Wang 2017) has shown that the NIR-based b bp (λ) algorithm can be applied to highly turbid waters to derive reasonably accurate b bp (λ) in the visible and NIR bands. It can be safely used for highly turbid waters with nL w (745) and nL w (862) less than~6 mW cm −2 μm −1 sr −1 and 4 mW cm −2 μm −1 sr −1 , respectively. The Hydrolight simulations show that the NIR-based algorithm not only can be used to derive b bp (λ) spectra, but also other IOP parameters such as a ph (λ) and a dg (λ) with reasonably good accuracy over highly turbid waters. On the other hand, atmospheric correction uncertainty in the VIIRS NIR bands can be significant for nL w (745) and nL w (862)  a t (λ) a dg (λ) a ph (λ) b bp (λ) b bp (λ), a ph (λ), a dg (λ), a t (λ) (m -1 ) Fig. 10. VIIRS-measured climatology b bp (λ), a ph (λ), a dg (λ), and a t (λ) in US east coastal region derived using the NIR-QAA combined IOP algorithm for VIIRS bands at (a-d) 443 nm, (e-h) 551 nm, and (i-l) 671 nm.