Variations of Water Transparency and Impact Factors in the Bohai and Yellow Seas from Satellite Observations

Variations of


Introduction
In marine and lake ecosystems, the penetration and availability of light underwater is a controlling factor of the physical (sediment suspension and heat transfer), biological (phytoplankton photosynthesis), and chemical (cycling of nutrients) processes of the water body [1].Therefore, water transparency is a key parameter for describing the optical properties in the water body.Water transparency is closely related to many environmental parameters, such as chlorophyll-a (Chl-a) [2][3][4], phytoplankton biomass [5], total suspended matter (TSM) [6], nutritional status [7], sea surface salinity [8], total nitrogen and total phosphorus [9], and widely used as evaluation indicators of marine environmental quality [10].
Remote Sens. 2021, 13, 514 2 of 19 The Bohai and Yellow Seas (BYS) are located in the northwestern margin of the Pacific Ocean [11].The sea area has a complex hydrodynamic environment, and also affected by a large amount of freshwater and sediments from the surrounding two world-famous rivers (the Yellow River and the Yangtze River).The marine community is rich and diverse, with clear seasonal changes of water optical properties.Due to the complexity of the optical properties and the importance of the marine ecological environment, it is desired to obtain the long-term monitoring data of water transparency and establish its relationship with environmental factors in the BYS, to improve the understanding of the dynamic variations of marine ecosystems in the region.
For the in situ water transparency measurement, the Secchi disk, which is the oldest "optical instrument" and dates from 1865, is a simple, fast, and low-cost commonly used tool [12].The Secchi disk depth (SDD) is the depth when a black and white Secchi disk disappears from the observer's sight in the human eye when it is deployed into the water [13] and is considered to be a reasonable estimate of water transparency at visible light range [14].Due to the limitations of a ship survey, the measured data can be discrete and sparse in time and space, and poorly synchronized, which can affect the detailed and objective evaluation of the variation of water transparency.As proved in the past few decades, the satellite sensors with repeated measurements of the same sea area can offer a wealth of ocean information, and improve the spatial and temporal coverage of environmental monitoring.It has been an important source of information for assessing global and regional aquatic environment conditions and trends [15][16][17].On-board ocean color sensors have been successfully used to monitor SDD, such as the Moderate Resolution Image Spectroradiometer (MODIS) [18][19][20].
In previous studies, the water transparency estimation using remote sensing approaches were mainly applied to inland water in China [21][22][23] with few reports on the coastal and marine regions of China seas.For example, the Bohai Sea (BHS) was found to have the lowest water transparency among the East China Seas from in situ collected data over a 21-year period (1959-1979) [24].Considering the residual error in satellite remote sensing reflectance (R rs ) data, an SDD empirical retrieval model coupled with a band difference method was developed for China Eastern Coastal Zone and shows that increasingly more waters within 30-35 and 20-25 m isobaths were deteriorating from 2002 to 2014 [25].For the central BHS, accompanying rapid economic growth at the end of the 1980s, there was a significant decrease in the water quality, but no significant variations were detected between 2003 and 2014 [26].Using a regionally tuned model from the Geostationary Ocean Color Imager (GOCI) data over the BYS, the water transparency shows obvious temporal variations, the daily change is primarily driven by the changes in the solar zenith angle, and the monthly variation is mainly driven by the stability of the water body stratification [27].However, the main driving force of the long-term SDD changes in the BYS has not yet been evaluated and analyzed in detail.
Therefore, this study applied the new mechanistic model from Lee et al. (2015) [28] to MODIS R rs data from 2003 to 2019 in order to explore the SDD changes in the BYS.Marine environmental factors, including near-surface Chl-a concentration, TSM concentration, chromophoric dissolved organic matter (CDOM) concentration, sea surface temperature (SST), photosynthetically available radiation (PAR), sea surface salinity (SSS), and wind vectors were assessed for the study region using satellite observations.Therefore, the aims of this study are to (1) use long-term MODIS R rs monthly mean products to estimate the SDD values of the BYS from SDD algorithm of Lee et al. (2015), and analyze the spatial distribution and temporal variation trend of SDD from 2003 to 2019; (2) carry out statistical analysis between the main components from the EOF decomposition and the environmental factors that affect SDD changes, and clarify the main impact factors of SDD changes in the BYS; and (3) elucidate the water components and marine environmental factors for the driving mechanism of SDD changes.

Study Area
The Bohai and Yellow Seas (31 • -41 • N, 117 • -127 • E) are semiclosed continental shelf shallow seas, with depths generally less than 100 m [29], as shown in Figure 1.The Bohai Sea is a shallow shelf basin covering 7.70×10 4 km 2 and an average water depth of 18 m [30].The Yellow Sea (YS) is contiguous to the Bohai Sea through the Bohai Strait with an average water depth of about 44 m and coverage of 3.8×10 5 km 2 [31].Affected by the East Asian monsoon, a north wind prevails in winter, and the average wind speed is about 8-9 m s −1 .The high-salt tongue penetrates the northern part of the Yellow Sea and enters the Bohai Sea via the Bohai Strait.The circulation is largely affected by the sea surface wind and brings strong vertical mixing in the water body [32].The south wind prevails in summer with an average wind speed of 4-6 m s −1 , resulting in low surface salinity in the YS.Most areas of the Yellow Sea are regulated by the semidiurnal tide, which is the main driving force of the Bohai Sea environment [33].

Study Area
The Bohai and Yellow Seas (31°-41°N, 117°-127°E) are semiclosed continental shelf shallow seas, with depths generally less than 100 m [29], as shown in Figure 1.The Bohai Sea is a shallow shelf basin covering 7.70×10 4 km 2 and an average water depth of 18 m [30].The Yellow Sea (YS) is contiguous to the Bohai Sea through the Bohai Strait with an average water depth of about 44 m and coverage of 3.8×10 5 km 2 [31].Affected by the East Asian monsoon, a north wind prevails in winter, and the average wind speed is about 8-9 m s -1 .The high-salt tongue penetrates the northern part of the Yellow Sea and enters the Bohai Sea via the Bohai Strait.The circulation is largely affected by the sea surface wind and brings strong vertical mixing in the water body [32].The south wind prevails in summer with an average wind speed of 4-6 m s -1 , resulting in low surface salinity in the YS.Most areas of the Yellow Sea are regulated by the semidiurnal tide, which is the main driving force of the Bohai Sea environment [33].

NASA Ocean Color Data
Monthly Level 3 standard mapped image of Rrs at 412, 443, 469, 488, 531, 555, 660, 680, and 745 nm, taken at a spatial resolution of 4 km × 4 km, were acquired from MODIS-Aqua between January 2003 and December 2019.The data were downloaded from the National Aeronautics and Space Administration (NASA) Ocean Biology Distributed Active Archive Center (https://oceancolor.gsfc.nasa.gov/l3/).Meanwhile, sea surface Chl-a concentration, TSM concentration, and photosynthetically available radiation (PAR) calculated from monthly mean Rrs data were obtained for establishing the relationship with water transparency variations.Chl-a concentrations were calculated using the standard OC3/OC4 (OCx) blue-to-green band ratio algorithm [34] combined with the Hu color in-

NASA Ocean Color Data
Monthly Level 3 standard mapped image of R rs at 412, 443, 469, 488, 531, 555, 660, 680, and 745 nm, taken at a spatial resolution of 4 km × 4 km, were acquired from MODIS-Aqua between January 2003 and December 2019.The data were downloaded from the National Aeronautics and Space Administration (NASA) Ocean Biology Distributed Active Archive Center (https://oceancolor.gsfc.nasa.gov/l3/).Meanwhile, sea surface Chl-a concentration, TSM concentration, and photosynthetically available radiation (PAR) calculated from monthly mean R rs data were obtained for establishing the relationship with water transparency variations.Chl-a concentrations were calculated using the standard OC3/OC4 (OCx) blue-to-green band ratio algorithm [34] combined with the Hu color index algorithm [35].The TSM concentration was taken by applying the algorithm constructed in Zhang et al. (2010) [36].The PAR is derived from the solar energy at the top of the atmosphere that removed the ocean-atmosphere system reflection and the atmosphere absorption, and commonly related to marine primary productivity [37].CDOM concentration represented by CDOM absorption coefficients at 440 nm (a CDOM (440)) can be separated from absorption coefficients of CDOM and nonalgal particles based on the QAA_CDOM algorithm [38,39].

GHRSST Products
The daily gridded AVHRR SST data at 4 km × 4 km spatial resolution were provided by the Group for High-Resolution Sea Surface Temperature (GHRSST, https://data.nodc.noaa.gov/ghrsst/L3C/)[40,41].Monthly mean SST time series data were composited from daily averaged SST data from January 2003 to December 2019.

RMESS Products
The wind vector data near the ocean surface (including wind speed and wind direction) were derived from Quick Scatterometer (QuikSCAT) and the first Advanced Scatterometer (ASCAT) of Remote Sensing Systems (RSS, http://www.remss.com/).Compared with QuikSCAT, ASCAT data shows a high level of consistency at all wind speed regimes by using a similar methodology and wind algorithm [42,43].Monthly gridded data of wind vectors from January 2003 to November 2009 were derived from QuikSCAT, and data from March 2007 and December 2019 were obtained from ASCAT.For products with the same month in both datasets, the specified variables were averaged in the data grid to generate monthly average fusion data.

CMEMS Dataset
SSS data were taken from the European Copernicus Marine Environment Monitoring Service (CMEMS, https://marine.copernicus.eu/).Global ocean multiobservation reanalysis products in the framework of CMEMS operationally provide dynamics information of the ocean and marine ecosystems and rely on in situ measured networks to calibrate and validate the data from satellites [44].The SSS L4 dataset was obtained by a multidimensional optimal interpolation technique from interpolating in situ and satellite-derived SSS data [45].The grid data were resampled to 4 km × 4 km resolution from the raw spatial resolution of 0.25 • × 0.25 • (about 25 km × 25 km).

Algorithm to Retrieve SDD
Based on classical underwater visibility theory, an analytical model for SDD retrieval has been developed [46], but the modeled data are not in agreement with in situ measurements.To exactly interpret the physical processes of sighting a Secchi disk in water by the human eye, a new underwater visibility theory described by Lee et al. (2015) was then proposed and an innovative mechanistic model for SDD was established.Using this model to retrieve water transparency, it mainly includes three main steps.First, the quasianalytical algorithm (QAA) is used to estimate the total absorption coefficient a(λ) and the total scattering coefficient b b (λ) of the water from R rs (λ) data [47].If R rs (670) is less than 0.0015 sr −1 , the reference band is taken 550 nm, otherwise, 670 nm is selected.Second, according to the model developed by Lee et al. (2005Lee et al. ( , 2013) ) [48,49], the diffuse attenuation coefficient K d (λ) is evaluated as follows: Here, the four parameters m 0-3 that vary with solar altitudes and depth owing to the change of light distribution are queried through a lookup table.The variable θ s is the solar zenith angle and η w is the ratio of the pure seawater backscattering coefficient b bw (λ) to b b (λ).The pure seawater coefficients of absorption and backscattering at each wavelength were tabulated by Morel [50] and Smith and Baker [51].Finally, following this new underwater visibility theory [28], SDD is given by SDD where Min(K d (λ)) is the minimum K d value in the visible domain wavelength of MODIS, R rs PC is the R rs value matched to the wavelength with minimum K d , and C t r is the contrast threshold for human eyes to sight a Secchi disk and taken as an average of 0.013 sr −1 based on the measurements of Blackwell (1946) [52].For global combined in situ datasets of 338 R rs and SDD matchup with SDD range of 0.1 m to 30 m, which contains 197 pairs of matching points from field measurements off China Seas, the accuracy of the new model was statistically verified.The mean absolute relative difference between measured and model retrieved SDD is 18.2%, and this difference includes 10% or more uncertainties in the observation of Secchi disk depth.

Trend Analysis
To evaluate the temporal variation of water transparency, the Mann-Kendall (MK) test was carried out to evaluate the temporal changes and significance of trends.Meanwhile, Sen's slope (SS) was used to obtain the magnitude of variation [53].The combination of both methods has been extensively used for trend analysis in hydrometeorological time series.In the MK test, the Z score obeys the standard normal distribution.When its absolute value is greater than 1.96, the evaluated trends are significant at 5% significance level [54]; it can be used to estimate the significance of the variation trend.
Regression analysis was used to evaluate the quantitative relationship between dependent and independent variables.For the single variable linear regression of the observed factors for time, the linear trend can be obtained from the least square fitting [55].Before the regression was applied, all assumptions of the analysis (i.e., normality, homoscedasticity, and independence of predictor variables) were tested and met.Specifically, we used this test to verify that there was no autocorrelation between predictor variables.

Empirical Orthogonal Function Decomposition
The empirical orthogonal function (EOF) is a principal component analysis (PCA) method that was first introduced into meteorology and climate research by Lorenz in the 1950s [56], and now it has been widely used in the field of oceanography [57].In general, the eigenvector corresponds to the spatial sample, so called the spatial mode (hereafter referred to as EOFs); the principal components (PCs) correspond to the time change.It can identify key variability patterns from the time series data of the spatial map.When a certain pattern can be statistically related to different physical processes, this allows for the underlying driving mechanisms to be evaluated.
When EOF decomposition is applied to the special parameters observation, the long time series data over spatial gridded points should be taken for anomaly processing first.Then, the eigenvalues and eigenvectors can be calculated easily using singular value decomposition (SVD) techniques to substitute for the classical covariance matrix approach to save the entire computational time.Arranging the eigenvalues in descending order, a list of eigenvectors corresponding to each nonzero eigenvalue is taken as a spatial distribution mode.The temporal mode can be found by multiplying the eigenvectors with the transpose of the original matrix.The size of the eigenvalue determines the contribution of the corresponding eigenvector to the total variance.The first EOF mode represents the long-term trend of time series changes, while higher-order EOF harmonics are detrended.At last, the North test is chosen to test the significance of the error range of the eigenvalue and judge whether the result of decomposition is meaningless noise or a signal with physical significance [58].It is worth noting that since each independent variability pattern is spatially dependent, the EOF analysis results would rely on the spatial area of the data grid involved.

Interannual Variations in SDD Dynamics
For the MODIS L3 monthly mean R rs images of BYS, the percentage of valid pixels is about 80.4%.The missing data of satellite remote sensing appeared mostly in winter and was mainly caused by cloud cover, contamination of sun glint, atmosphere correction, or adjacency from land.With the confidence in the model description of Lee et al. (2015) in the China seas, monthly mean SDD images of BYS were estimated from MODIS monthly mean R rs .The spatial distribution of climatological annual mean SDD images from 2003 to 2019 was calculated (Figure 2), and lower/higher-value area is displayed.Overall, it shows a large SDD variation range in the study area (0.4-15.6 m).The lower-value area is found in the Bohai Sea (0.5-4.3 m), and the higher-value area appeared in the center South Yellow Sea, which is consistent with the variation trend of the water depth contour marked in Figure 1.With the influence of the shallow water depth, the terrestrial input by the Yellow River and the Yangtze River, and the coastal currents, the transparency inshore where the water depth is below 50 m is relatively low, and in contrast, the water transparency is high in the deep ocean far from the coastline, which is less affected by the continent.There are tongue-shaped low-transparency areas near Chengshanjiao and southeast of Subei Shoal.The standard deviation of the SDD of the entire sea area is less than 6 m.The lack of data in the northern Jiangsu Shoal and the coastal waters of Hangzhou Bay are largely due to the failure of the atmospheric correction algorithm in the coastal turbid waters.
spatially dependent, the EOF analysis results would rely on the spatial area of the data grid involved.

Interannual Variations in SDD Dynamics
For the MODIS L3 monthly mean Rrs images of BYS, the percentage of valid pixels is about 80.4%.The missing data of satellite remote sensing appeared mostly in winter and was mainly caused by cloud cover, contamination of sun glint, atmosphere correction, or adjacency from land.With the confidence in the model description of Lee et al. (2015) in the China seas, monthly mean SDD images of BYS were estimated from MODIS monthly mean Rrs.The spatial distribution of climatological annual mean SDD images from 2003 to 2019 was calculated (Figure 2), and lower/higher-value area is displayed.Overall, it shows a large SDD variation range in the study area (0.4-15.6 m).The lower-value area is found in the Bohai Sea (0.5-4.3 m), and the higher-value area appeared in the center South Yellow Sea, which is consistent with the variation trend of the water depth contour marked in Figure 1.With the influence of the shallow water depth, the terrestrial input by the Yellow River and the Yangtze River, and the coastal currents, the transparency inshore where the water depth is below 50 m is relatively low, and in contrast, the water transparency is high in the deep ocean far from the coastline, which is less affected by the continent.There are tongue-shaped low-transparency areas near Chengshanjiao and southeast of Subei Shoal.The standard deviation of the SDD of the entire sea area is less than 6 m.The lack of data in the northern Jiangsu Shoal and the coastal waters of Hangzhou Bay are largely due to the failure of the atmospheric correction algorithm in the coastal turbid waters.The MK and SS analysis results for each region in Figure 2 show that the water transparency of most areas of the BYS did not obtain statistical clustering characteristics (1.96 > Z-value > −1.96); insignificant increasing or decreasing trend was observed except for the coastal waters of Haizhou Bay, which showed a significant negative trend (Z-value < −1.96).Although the transparency of the coastal waters in the BHS was also significantly increased, the three bays (Laizhou Bay, Liaodong Bay, and Bohai Bay) have about 3 to 4 months of the ice in winter, which will cause most of the SDD data missing in winter.Therefore, trend analysis is not statistically significant.The SDD in the central BHS has increased, but the time trend is not significant.The NYS area showed an increase in SDD in the northern water depth of less than 50 m, and a negative trend covers the spatial extent of other areas.SDD in the SYS showed an increasing temporal trend in the central part of the YS.
Taking January and July as the representative months of winter and summer, respectively, based on the spatial average of the monthly data of each sea area, linear regression was performed on the data of SDD time anomalies, and the annual average time series were also statistically analyzed (Figure 3).Due to the existence of the winter ice period, the three bays of the BHS were excluded on average and only consider the central region.Regarding the trend of SDD anomalies, the temporal variations of the three sea areas all show a certain degree of fluctuation.In winter, both NYS and SYS areas decreased significantly, and BHS showed an insignificant decrease.In summer, NYS showed an increasing trend, SYS and BHS showed a decrease in variation, but none of the trends were significant.For the insignificant trend, on the one hand, when the new transparency algorithm was established and compared with the measured data, 10% of the uncertainties came from the visual deviation of the in situ measurements; on the other hand, it can be attributed to the use of a unified atmospheric correction algorithm when acquiring MODIS R rs product data, which leads to a large deviation in the SDD data obtained through inversion.The influence of these uncertain errors was ignored in this study.In terms of annual changes, the SDD of the entire study sea area shows a decreased trend, and the variation rate is lowest in the BHS at 0.003 m y −1 and highest in the SYS at 0.015 m y −1 .The trend analysis of the SDD anomaly also shows that before 2009, the water transparency was higher than the annual mean, and lower from 2009 to 2017, which also means that the water quality of the entire sea area may have deteriorated during this period.Since 2017, it shows a weak increasing trend, which implies that the water quality in the study sea area seems to have improved.The variations in water transparency are most likely related to changes in human activities and weather conditions, so more efforts are needed to maintain a healthy marine ecosystem.

Seasonal Patterns of SDD and Related Marine Environmental Factors
Climatological SDD time series data also indicated seasonal patterns of low in spring and winter and high in summer and autumn, as displayed in Figure 4. Related marine environmental factors, including Chl-a, TSM, CDOM, SST, PAR, SSS, and wind vectors (wind speed and direction) are all given in similar patterns.
In spring, the water transparency is low, less than 4 m in the Bohai Sea, the coastline of the Yellow Sea, Subei Shoal, and Yangtze River estuary.Values are higher in the central South Yellow Sea, reaching 14 m.The low-value area near Chengshanjiao is related to the strong mixing caused by the flow around the peninsula, and also to the sediment of the Yellow River carried by the North Shandong coastal currents.When the SST gradually rises in summer, the increasing thickness of the thermocline strengthens the vertical stability of coastal waters, and the upper and lower layers of seawater are not easily mixed [59], making this the most transparent time of the year.The water transparency is 1.6-8 m in the Bohai Sea and 2.4-19.4m in the Yellow Sea, which is significantly higher than that in spring.The low-value area in the east of the Subei Shoal decreases to the smallest area of the year.In autumn, as the northwest wind weakens and the northwest wind strengthens, the convective mixing of seawater gradually increases, and the concentration of suspended matter increases, which leads to a general decrease in water transparency.The high-value area in the center of the Bohai Sea no longer exists, and water transparency on the entire sea area is relatively uniform, with its values <5 m.The water transparency range of the Yellow Sea is 0.2-16 m.The water transparency in winter is the lowest among the four seasons.The water transparency in the Bohai Sea is 0.1-2 m, and the spatial difference is the smallest.The transparency of the coastal waters of the Yellow Sea is less than 3 m.Most of the central sea areas in the Yellow Sea are still in high transparency (6-8 m).From the climatological image of the optical components of the seawater, the regional distribution trends of Chl-a, TSM, and CDOM in the Bohai Sea are relatively consistent: from the high-value area of Bohai Bay, Laizhou Bay, and Liaodong Bay to the low-value area of Bohai Strait.All factors show regular gradient changes, and the isoline is basically consistent with the trend of the coastline.Near the Subei Shoal in the Yellow Sea, coastal currents and terrestrial input leads to the high content of suspended particles and chromophoric dissolved organic matter throughout the whole year.The center of the North Yellow Sea and the South Yellow Sea are mainly affected by the Yellow Sea warm current.The seawater of this area contains a certain portion of suspended sediment and nutrient matter.The Chl-a concentration is relatively high in spring and winter.The coastal waters of the North Yellow Sea and the Shandong Peninsula are affected by the combined action of density currents and coastal currents, and the contours of each factor are almost parallel to the coastline.
As seen from the series of climatological changes of marine environmental factors, the spatial patterns of PAR, SST, and wind vectors are nearly coherent in the seasonal cycles.The PAR is mainly related to solar radiation intensity.Accounting for the difference in solar altitude angles of different seasons, SST presents the characteristics of low in winter and high in summer, and its spatial distribution roughly corresponds to the latitude zone.On account of the East Asian monsoon [60], northerly winds prevail in winter, the average wind speed is about 7-8 m s −1 ; spring monsoons alternate and wind directions become unstable; in summer, southerly and southeast winds prevail and the average wind speed is about 4-6 m s −1 .In terms of the spatial distribution of SSS, the salinity in the near-shore is low, resulting from the input of freshwater from the river, and gradually increases from inshore to offshore.Winter and spring are mainly under the influence of the Yellow Sea warm water mass, and the high-salt water body carried is tongue-shaped, extending to the central Yellow Sea.In summer and autumn, two low-salt areas formed near the mouths of the Yellow River and the Yangtze River are more obvious, and the freshwater input by the rivers extends from the coastal areas to the deep ocean.

Variations of Marine Environmental Factors and Relevance for SDD
The principal modes obtained by the EOF decomposition can be used to achieve the dominant space and temporal patterns and main impact factors from the time series satellite images [61][62][63].A strong environmental factor can appear in several principal modes, and there many environmental factors can take effect in one mode.Figure 5 shows spatial patterns (a-c) and temporal amplitudes (d) for the first three SDD EOF modes.
The first SDD mode explains 58.75% of the total variance.The change of spatial patterns is not obviously apparent.The amplitude value is positive everywhere, but the time amplitude changes drastically, indicating that the SDD temporal trends are consistent throughout the study area.The amplitude of the spatial pattern increases from the coast to the open sea, with the largest changes found in the center of the Yellow Sea.It is similar to the annual average SDD image.The corresponding time mode reveals a clear seasonal cycle change signal, which is positive from May to September, indicating that the water transparency is higher than the annual mean, and the other months are negative, indicating that the transparency is lower.It also shows that seasonal ocean and atmospheric processes, such as PAR, SST, and wind vectors, are the important factors that cause periodic changes of the mode.Since the optical components of water also have seasonal cycle changes, they must be included when analyzing the physical processes related to leading patterns.The second and third EOF modes of SDD play a small role in SDD changes.Higher-order EOFs are not considered in this study, and only explain the minor component of the total variance.modes, and there many environmental factors can take effect in one mode.Figure 5 shows spatial patterns (a-c) and temporal amplitudes (d) for the first three SDD EOF modes.The first SDD mode explains 58.75% of the total variance.The change of spatial patterns is not obviously apparent.The amplitude value is positive everywhere, but the time amplitude changes drastically, indicating that the SDD temporal trends are consistent throughout the study area.The amplitude of the spatial pattern increases from the coast to the open sea, with the largest changes found in the center of the Yellow Sea.It is similar to the annual average SDD image.The corresponding time mode reveals a clear seasonal cycle change signal, which is positive from May to September, indicating that the water transparency is higher than the annual mean, and the other months are negative, indicating that the transparency is lower.It also shows that seasonal ocean and atmospheric processes, such as PAR, SST, and wind vectors, are the important factors that cause periodic changes of the mode.Since the optical components of water also have seasonal cycle changes, they must be included when analyzing the physical processes related to leading patterns.The second and third EOF modes of SDD play a small role in SDD changes.Higher-order EOFs are not considered in this study, and only explain the minor component of the total variance.
The second mode of EOF accounts for 4.14% of the total variance.Combined with the spatial patterns and temporal amplitudes, it decreases in summer and increases in autumn.The transparency changes between the center Yellow Sea and the region near the Yangtze River estuary are opposite.The presence of a thermocline generated by the increase in surface water temperature in the central Yellow Sea in summer prevents the suspended particles at the bottom from mixing with the upper water body.In offshore wa- The second mode of EOF accounts for 4.14% of the total variance.Combined with the spatial patterns and temporal amplitudes, it decreases in summer and increases in autumn.The transparency changes between the center Yellow Sea and the region near the Yangtze River estuary are opposite.The presence of a thermocline generated by the increase in surface water temperature in the central Yellow Sea in summer prevents the suspended particles at the bottom from mixing with the upper water body.In offshore waters, the shallower water depth and the shallower thermocline weaken the effect of hindering the resuspension of particles.The change of transparency in autumn may be affected by the decrease of freshwater drained down from the Yangtze River estuary and the high-salt water body brought by the Yellow Sea warm water mass.Therefore, this mode is closely related to SST and SSS.The third mode of EOF contains 2.97% of the total variance, and the explained SDD feature is weaker.It shows that the SDD in the central Yellow Sea increases in autumn and winter, whereas the SDD near the coastal area decreases.This may be owing to the increased sea surface wind intensity in autumn and the effect of suspended sediment on the bottom of the coastal waters, which leads to a greater degree of transparency than the central South Yellow Sea.The large spatial amplitude of the southern branch of the Shandong Peninsula and the Bohai Sea is perhaps caused by seasonal blooms from the changes of Chl-a.
Since the overall transparency of the Bohai Sea and the North Yellow Sea are relatively low, the water properties can be regarded to be similar.It should be noted that since the spatial area is relatively small, the characteristics of the signal in this area can be easily covered by the strong signal of the South Yellow Sea when EOF is decomposed.We also performed EOF analysis in the BNYS and SYS areas at the same time.The spatial and temporal modes from the EOF decomposition of long-term data are not shown in this study.
To clarify what environmental factors cause the time changes of SDD in the area of BNYS, SYS, and BYS respectively, Table 1 displays the correlation coefficient (r) between SDD EOF temporal modes and the marine environmental factors in each region.The first EOF mode in the BYS shows a strong negative correlation with Chl-a, TSM, a CDOM (440), SSS, wind speed, and a significant positive correlation with PAR and SST.This mode has the highest correlation with TSM (r = −0.910).Although the weakest correlation with Chl-a (r = −0.508),strong correlation with TSM and CDOM (r = −0.905)illustrates that the change of SDD is determined by the optical components of the water body.The correlation between the second mode and PAR, SST, and SSS is significant, but the correlation is weaker.The third mode has the highest correlation with wind speed, but not strong.It can be taken that PAR, SST, SSS, and wind speed are covariates of SDD changes.It is consistent with the analysis above for each SDD EOF mode on BYS.The contribution rate of PC1 in the BNYS area to the total variance is about 72.64%, which is higher than that of SYS and BYS.The first mode also has a higher correlation with TSM and CDOM, whereas the second and third modes are significantly correlated with SSS and Chl-a, respectively.The correlation between SDD changes and environmental factors in SYS is similar to that of BYS overall, but the first mode has the highest correlation with CDOM.As the SDD changes in the entire sea area of the BYS have the highest correlation with TSM, and the relationship between the water optical properties and the water components (such as Chl-a concentration and TSM concentration) strongly depends on the composition of suspended particles in water, so it is necessary to characterize the dominant suspended particulate matter types in the water.The ratio of the Chl-a concentration to TSM concentration (Chl-a/TSM) can be considered as a good indicator to quantify the composition of suspended particulates.When the ratio is greater than 0.5 × 10 −3 , it means that the optical properties of the water body are mainly dominated by phytoplankton [64][65][66].Analyzing the Chl-a/TSM ratio of the six stations (Figure 6), the values are concentrated in the range of 0.25 × 10 −3 -3.0 × 10 −3 , most of which are in the range of 0.5 × 10 −3 -2.5 × 10 −3 .In terms of seasonal changes, S1 and S6 indicate a trend of being high in summer and autumn, and show that the Chl-a/TSM ratio of the sea surface at the S6 station is lower than 0.5 × 10 −3 in winter and spring, which means the suspended particles are mainly inorganic and mainly resuspended sedimentary particles.Station S1 is strongly mixed in both the upper and lower water column layers, and only shows characteristics of inorganic particles in winter.Other stations show more or less phytoplankton-dominated features throughout the year, which is that the remote sensing reflectance spectrum of the water in these stations shows a higher absorption coefficient in the red band.Stations S2, S4, and S5 showed double peaks in spring and autumn.The water at all stations has a low Chl-a/SPM ratio in spring but a high TSM concentration, which may be carrying more nutrients and more likely to cause algal blooms [67].
particles in winter.Other stations show more or less phytoplankton-dominated features throughout the year, which is that the remote sensing reflectance spectrum of the water in these stations shows a higher absorption coefficient in the red band.Stations S2, S4, and S5 showed double peaks in spring and autumn.The water at all stations has a low Chla/SPM ratio in spring but a high TSM concentration, which may be carrying more nutrients and more likely to cause algal blooms [67].Further analysis of the role played by marine environmental factors, SDD changes, and the corresponding seasonal changes of marine environmental factors at the six stations is presented in Figure 7. Overall, the water transparency is the greatest at S4 and S6, and lowest at S1 and S5.The seasonal variation trends of SST, PAR, SSS, and wind speed at each station are roughly the same, but the changes in the optical components of the water body are quite different.The detailed contrastive analysis shows the following features: (1) The SDD changes of S1 and S5 are only different in July; S1 continues to increase whereas S5 decreases.Station S1 is located near the mouth of the Yellow River in the Bohai Sea, and station S5 is near the mouth of the Yangtze River.Comparing the changes in related marine environmental factors from June to August, this trend of change is caused by the optical components of the water at station S5 in July, and can be attributed to the huge amount of freshwater from the Yangtze River estuary.(2) The variation trends of S2 and S3 are similar in spring, but their difference increases afterward.Different from S2, S3 has a higher concentration of Chl-a and CDOM that leads to lower transparency.It may be caused by the large-scale algal blooms that occurred every summer in the western part of the Yellow Sea since 2008, which led to the increase of Chl-a concentration and the decrease of water transparency [68,69].(3) The water depth at S4 and S6 are greater than 50 m, and the water transparency varies from July to October.Although the TSM at each station is relatively low, S6 is adjacent to the diffusion area of the Yangtze River estuary freshwater.The increase in nutrients in seawater drives the photosynthesis of phytoplankton, and Chl-a increases correspondingly and reduces the water clarity.
of the Yellow Sea since 2008, which led to the increase of Chl-a concentration and the decrease of water transparency [68,69].(3) The water depth at S4 and S6 are greater than 50 meters, and the water transparency varies from July to October.Although the TSM at each station is relatively low, S6 is adjacent to the diffusion area of the Yangtze River estuary freshwater.The increase in nutrients in seawater drives the photosynthesis of phytoplankton, and Chl-a increases correspondingly and reduces the water clarity.

Discussion
In natural waters, the downward spectral radiation within water bodies decreases with increasing water depth; the attenuation process underwater is closely related to the absorption and scattering of the optical components (water molecules, phytoplankton, total suspended particulate matter, and CDOM) [70,71].In water dominated by inorganic particulate matter, scattering is relatively more significant to SDD than absorption [1].In the visible light range, water molecules strongly absorb at red wavelength, whereas scattering is very weak.Phytoplankton has two strong absorption peaks at the blue and red wavelength, and strong scattering at the same time [72].The total suspended particulate matter is composed of the inorganic part and organic part.Inorganic particulate matter has a high refractive index relative to water molecules, so the light underwater would be strongly scattered, but the optical properties of debris-like organic particulate matter are similar to CDOM [73].CDOM can absorb the light underwater at the blue band and scatter at the yellow band, thus making the sea surface water present as yellow [74].The Bohai and the Yellow Seas have typical water characterized as Case II.Apart from water molecules, the three main optical components jointly determine the optical characteristics of the water body.Their absorption and scattering characteristics overlap each other during the underwater light penetration process, and each contribution to SDD is nonlinear [75,76].The long-term SDD changes in the BYS are highly correlated with TSM and CDOM, indicating that the underwater light attenuation of the entire sea area is mainly due to the combined effect of these two optical components.The fact that chlorophyll-a is not highly correlated with SDD also shows that the contribution of phytoplankton in the study area is masked by other water components.Generally speaking, the influence between phytoplankton and water transparency is reciprocal.On the one hand, the absorption of underwater light by phytoplankton biomass will reduce the amount of light that can penetrate the water body.On the other hand, the reduction of light underwater results in lower carbohydrate content in phytoplankton cells and restricts the growth of phytoplankton [77].Although the SDD changes in the BYS are mainly caused by TSM and CDOM, there are significant differences in the weights of the three optical components for different types of water [78], and the driving factors for SDD changes are also different [79].
From the perspective of ocean optics and ecology, the effects of marine environmental factors on the changes in water transparency cannot be ignored.As analyzed in Table 1, PAR, SST, SSS, and wind speed are covariates that cause SDD changes in the BYS.These factors directly vary the distribution of underwater light and the process of underwater ecosystems, thereby affecting the primary production of seaweed, microalgae, and benthic organisms [80].The intensity of PAR, corresponding to the intensity of solar radiation, determines the amount of input that may enter the underwater light field.Due to its regular seasonal changes, it cannot be a decisive variable in the SDD variation.There is also a significant positive correlation between SST and SDD (Figure 7).In summer, the sea surface absorbs more solar radiation, and the sea surface temperature rises.The upper-ocean water body gradually forms a stable thermocline, which prevents upward movements of the bottom particles and nutrients.The material is mixed with the upper water body, so the water transparency is generally higher during this period.At the same time, the increase in water transparency promotes the penetration of solar irradiance, increases the light absorption in the deeper part of the water, and then varies the heat transmission of the upper and middle water bodies of the ocean [81].Since the global SST has increased year-by-year since 2000 [82,83], it may cause an increase in interannual SDD in the deep-sea area of the Yellow Sea (stations S4 and S6).On the spatial distribution of SSS, low-salinity areas are formed near the coast due to runoff input.A large amount of suspended sediment and CDOM accumulate in the water body, and the CDOM concentration is significantly higher in low-salinity water bodies, but the sediment input is usually trapped there with the high sedimentation rate caused by regional circulation [9,84].The high-salinity areas in the central Yellow Sea are the mixing of coastal low-salt water bodies with high-salt water bodies carried by the Yellow Sea warm water mass.The significant positive correlation between SDD and SSS (Table 1) further highlights the powerful role of CDOM in light attenuation in water bodies.The wind field is the most important factor that affects the TSM changes, especially in the coastal estuary area.The central Yellow Sea is deeper.The stable thermocline in summer prevents the suspended particles driven by the tide and circulation from reaching the upper water body, and the thermal layer gradually disappeared in autumn, the strong wind intensified the mixing of the upper and lower water bodies, resulting in more bottom sediments suspended.
Furthermore, in the analysis of long-term SDD changes, it is essential to take account of the influences of water depth changes caused by topography and tides [85,86].As shown in Figure 2, the depth of the coastal waters on the west side of the station S3 and station S5 increases with the distance from the coastline.The SDD changes also show a similar increasing trend, and the contour is approximately parallel to the coastline.For deep-sea areas, especially water areas with fast currents and large tidal height changes, for example, the transparency of seawater at station S6 is lower than that at station S5; it is also considered to be affected by the mixing of tidal currents.

Conclusions
In this study, climatological SDD images from 2003 to 2019 are estimated using the innovative mechanistic model of Lee et al. (2015) over MODIS L3 monthly mean R rs data of BYS.The spatial distribution of annual mean SDD shows lower SDD values (0.5-4.3 m) in the Bohai Sea, but higher values in the center of the South Yellow Sea with a large variation range (0.4-15.6 m).The terrestrial input from the Yellow River and the Yangtze River, and coastal currents, lead to water transparency that is relatively low in the inshore area where the water depth is below 50 m but high in the deep ocean.There are tongueshaped low-transparency areas near Chengshanjiao and southeast of Subei Shoal.The standard deviation of the SDD of the entire sea area is less than 6 m.The seasonal SDD variations also indicate seasonal patterns of low in spring and winter and high in summer and autumn.On the SDD annual changes analysis of the six stations, water transparency at S1 and S3 both decline and other stations increase.
To determine the main impact factors from the time series satellite images, EOF decomposition was carried out.The dominant space and temporal patterns obtained from EOF analysis show that the first SDD EOF mode in the BYS is the highest, and closely correlated with TSM.The high correlation coefficients of CDOM also illustrate that the SDD variation is mainly dominated by the optical components in the seawater, although correlation with Chl-a is the weakest.The contribution of each optical component to the changes of SDD is found to be different in the seawater with different content of optical components.Meanwhile, the dominant suspended particulate matter composition is characterized by Chl-a/TSM.From the second and third EOF modes, PAR, SST, SSS, and wind speed are found to be the covariates that cause SDD changes in the BYS.Furthermore, the impact of changes in water depth caused by topography and tides also need to be taken into consideration.Future research of SDD changes in the marine ecosystem should consider other influencing factors, such as phytoplankton biomass, algal bloom, and river discharge.

Figure 1 .
Figure 1.Study area over the Bohai and Yellow Seas.Contours show 20, 50, and 100 m isobaths from ETOPO1-1 Arc-Minute Global Relief model data.

Figure 1 .
Figure 1.Study area over the Bohai and Yellow Seas.Contours show 20, 50, and 100 m isobaths from ETOPO1-1 Arc-Minute Global Relief model data.

Figure 2 .
Figure 2. Spatial distribution of the climatological annual mean of Secchi disk depth (SDD) data derived from MODIS for 2003-2019 (left panel), and variation trend of each sea area using the Mann-Kendall (MK) test (middle panel) and Sen's

Figure 5 .
Figure 5. Spatial patterns (a-c) and temporal amplitudes (d) for the first three SDD empirical orthogonal function (EOF) modes.The percentage of total variance explained by the principal component 1 to 3 is 58.75%, 4.14%, and 2.97%, respectively.

Figure 5 .
Figure 5. Spatial patterns (a-c) and temporal amplitudes (d) for the first three SDD empirical orthogonal function (EOF) modes.The percentage of total variance explained by the principal component 1 to 3 is 58.75%, 4.14%, and 2.97%, respectively.

Figure 6 .
Figure 6.Monthly variations of Chl-a/TSM ratios from the six stations (a-f).Each boxplot chart includes the minimum, maximum, median, first quartile, third quartile, and outliers.

Figure 6 .
Figure 6.Monthly variations of Chl-a/TSM ratios from the six stations (a-f).Each boxplot chart includes the minimum, maximum, median, first quartile, third quartile, and outliers.

Table 1 .
Correlation coefficient between SDD EOF temporal modes and marine environmental factors in BNYS, SYS, and BYS.