Analysis of the Monthly and Spring-Neap Tidal Variability of Satellite Chlorophyll-a and Total Suspended Matter in a Turbid Coastal Ocean Using the DINEOF Method

: Missing spatial data is one of the major concerns associated with the application of satellite data. The Data INterpolating Empirical Orthogonal Functions (DINEOF) method has been proven to be an effective tool for ﬁlling spatial gaps in various satellite data products. The Ariake Sea, which is a turbid coastal sea, shows the large spatial and temporal variability of chlorophyll-a (Chl-a) and total suspended matter (TSM). However, ocean color satellite data for this region usually have large gaps, which affects the accurate analysis of Chl-a and TSM variability. In this study, we applied the DINEOF method to ﬁll the missing pixels from the regionally tuned Moderate Resolution Imaging Spectroradiometer (MODIS)-Aqua (hereafter, MODIS) Chl-a and MODIS-derived TSM datasets for the period 2002–2017. The validation results showed that the DINEOF reconstructed data were accurate and reliable. Furthermore, the Empirical Orthogonal Functions (EOF) analysis based on the reconstructed data was used to quantitatively analyze the spatial and temporal variability of Chl-a and TSM at both monthly and individual events of spring-neap tidal scales. The ﬁrst three EOF modes of Chl-a showed seasonal variability mainly caused by precipitation, the sea surface temperature (SST), and river discharge for the ﬁrst EOF mode and the sea level amplitude for the second. The ﬁrst three EOF modes of TSM exhibited both seasonal and spring-neap tidal variability. The ﬁrst and second EOF modes of TSM displayed spring-neap tidal variability caused by the sea level amplitude. The second EOF mode of TSM also showed seasonal variability caused by the sea level amplitude. In this study, we ﬁrst applied the DINEOF method to reconstruct the satellite data and to capture the major spatial and temporal variability of Chl-a and TSM for the Ariake Sea. Our results demonstrate that the DINEOF method can reconstruct patchy oceanic color datasets and improve spatio-temporal variability analysis. Additionally, the r values were calculated for correlations between the sea level amplitude, wind speed, precipitation, sea surface temperature, and river discharge. Significant correlations were found between precipitation and river discharge (r = 0.93), precipitation and SST (r = 0.88), and river discharge and SST (r = 0.83).


Introduction
Chlorophyll-a (Chl-a) from ocean color remote sensing has been effectively used to investigate phytoplankton dynamics at seasonal and inter-annual scales, in order to obtain a better understanding of the role of phytoplankton in marine biogeochemistry; the global carbon cycle; and the response of marine ecosystems to climate variability, change, and feedback processes. Similarly, the total suspended sediment (TSM), which is a major constituent in coastal waters, is involved in processes such as light attenuation, pollutant propagation, and waterway blockages. The spatial distribution of TSM is an indicator of deposition and erosion patterns in estuaries and coastal zones and a necessary input deposition and erosion patterns in estuaries and coastal zones and a necessary input for estimating the material fluxes from the land through rivers to the sea [1][2][3][4][5]. However, satellite ocean color datasets are fraught with missing data caused by improper conditions, such as clouds and rain [3,6]. These missing data may cause bias in the analysis of satellite-derived time series [6], especially for short-term scales. Therefore, it is essential that missing data are filled to better utilize the datasets.
The Ariake Sea ( Figure 1) is a semi-enclosed embayment (~20 km wide and 100 km long) located in Kyushu Island of Japan. The waters are shallow, with an average depth of ~15 m, and ~5 m near the coasts. There are seven main river systems discharging into the Ariake Sea, namely the Rokkaku, Kase, Chikugo, Yabe, Kikuchi, Shira, Midori, and Kuma rivers, and they supply a large amount of nutrients and suspended sediments to the Ariake Sea [16]. Additionally, the Ariake Sea is a tidally-driven embayment and the spring-neap tidal range is the largest among Japanese coastal waters [17]. This environment makes the Ariake Sea a productive, turbid, and dynamic coastal water body. Figure 1. Location of the Ariake Sea, Japan (a). The water depth of the bay is shown in light to dark blue (b). The seven main rivers of Rokkaku, Kase, Chikugo, Yabe, Kikuchi, Shira, Midori, and Kuma are indicated by arrows. The three regional areas, i.e., Off Saga, Isahaya Bay, and off Kumamoto, are highlighted by red boxes. The observation station for tidal level data of Ariake Sea, named Oura, is represented by a magenta filled circle. The four observation stations for sea level amplitude, wind speed, precipitation, and sea surface temperature, are represented by blue filled circles.
Chl-a and TSM are two common indicators of water quality. There have been several efforts to estimate the spatial and temporal variability of Chl-a and TSM in the Ariake Sea using satellite ocean color data [1,18]. Using significantly improved, regionally-tuned Moderate Resolution Imaging Spectroradiometer (MODIS) datasets from 2002 to 2017, Yang et al. [1] reported that the Chl-a variability during the spring-neap tidal cycle in the region off Saga (Figure 1) was generally controlled by TSM in seasons other than summer, Figure 1. Location of the Ariake Sea, Japan (a). The water depth of the bay is shown in light to dark blue (b). The seven main rivers of Rokkaku, Kase, Chikugo, Yabe, Kikuchi, Shira, Midori, and Kuma are indicated by arrows. The three regional areas, i.e., Off Saga, Isahaya Bay, and off Kumamoto, are highlighted by red boxes. The observation station for tidal level data of Ariake Sea, named Oura, is represented by a magenta filled circle. The four observation stations for sea level amplitude, wind speed, precipitation, and sea surface temperature, are represented by blue filled circles.
Chl-a and TSM are two common indicators of water quality. There have been several efforts to estimate the spatial and temporal variability of Chl-a and TSM in the Ariake Sea using satellite ocean color data [1,18]. Using significantly improved, regionally-tuned Moderate Resolution Imaging Spectroradiometer (MODIS) datasets from 2002 to 2017, Yang et al. [1] reported that the Chl-a variability during the spring-neap tidal cycle in the region off Saga (Figure 1) was generally controlled by TSM in seasons other than summer, when the Chl-a variability was mostly controlled by river discharge. They also found that the influence of TSM on the Chl-a variability was small in Isahaya Bay and off Kumamoto (Figure 1) due to the low TSM concentrations, based on significantly improved, regionally-tuned MODIS datasets from 2002 to 2017. Based on the monthly average of Sea-viewing Wide Field-of-view Sensor (SeaWiFS) data from May 1998 to December 2001, Ishizaka et al. [18] reported that Chl-a peaks occurred in summer and fall following peaks in rain and river discharge based on the monthly average of Sea-viewing Wide Field-ofview Sensor (SeaWiFS) data from May 1998 to December 2001. However, the satellite ocean color data used in previous studies had many spatial gaps, which made it difficult to capture the short-term variability of Chl-a and TSM [1]. Gaps in data, both spatial and temporal, can also affect the analysis of climatology data.
Additionally, other environmental variables, including the sea level amplitude, sea surface temperature (SST), and wind, are also key factors responsible for Chl-a and/or TSM dynamics in the Ariake Sea. For the sea level amplitude, the individual events of spring-neap tidal cycles have often been used [1,19]. For example, the sea level amplitude is the highest and lowest during spring and neap tide, respectively. However, the seasonal pattern of the sea level amplitude and its effect on the Chl-a and TSM variability in the Ariake Sea has not yet been investigated. It was proposed that for estuaries with sufficient nutrients, light availability and SST are key factors for phytoplankton biomass indicated by the Chl-a concentration [20][21][22]. Moreover, Yamaguchi et al. [23] proposed that the SST in winter possibly has a strong influence on phytoplankton biomass dynamics in the waters. As far as wind is concerned, Arai et al. [24] found that the MODIS-derived Chl-a in the winters of 2010 and 2015 for the Ariake Sea often increased during calm ocean wind conditions, but not when the seawater was turbulent, due to the convection caused by strong northern winds.
In this study, the DINEOF method was applied to reconstruct the regionally-tuned MODIS Chl-a and TSM datasets for the period from 2002 to 2017. Then, the major spatial and temporal variability of the reconstructed Chl-a and TSM was investigated based on the EOF analysis. Furthermore, the impact of the TSM, sea level amplitude, wind speed, precipitation, river discharge, and SST on the variability of Chl-a was quantitatively estimated using correlation analysis. The impact of the environmental variables on the variability of TSM was also quantitatively analyzed. We placed emphasis on the seasonal and spring-neap tidal variability of Chl-a and TSM using the monthly and individual events of spring-neap tidal cycle data. The main objective was to improve our understanding of the spatio-temporal variability of Chl-a and TSM and the mechanisms behind it with the DINEOF reconstructed data for the Ariake Sea.

Satellite Data and Pre-Processing
MODIS Aqua level 2 remote sensing reflectance (Rrs) and Chl-a data products of the reprocessed (2018.1) version, provided by the NASA Ocean Biology Processing Group (https://oceancolor.gsfc.nasa.gov/ (accessed on 21 December 2021)), were used in this study. The daily data from July 2002 to December 2017 with a 1 km and daily resolution and in a 169 (latitude) × 120 (longitude) pixel grid bounded by 32 • N to 33.33 • N latitude and 130 • E to 130.89 • E longitude were extracted. The domain covered included the Ariake Sea, which was the target area, as well as the Yatsushiro and Tachibana Bays, and a part of the East China Sea (Figure 1). There were only 1582 available daily data because of the problematic weather conditions, such as cloudy and rainy days. Flags of LAND, HIGLINT, HILT, HISATZEN, CLDICE, HISOLZEN, LOWLW, MAXAERITER, and NAVFAIL (Table 1) were used to mask the data in the pixel grid to exclude questionable data [25]. Finally, daily data with a spatial coverage of less than 20% were abandoned as those data were likely contaminated by the cloud edge. This reduced the dataset to only 899 daily data for further data processing [1]. The regionally-tuned Chl-a and TSM algorithms of Yang et al. [25] and Yang et al.
[1], respectively, were applied to the recalculated MODIS Rrs using the technique in Hayashi et al. [26], in order to obtain improved satellite retrievals of Chl-a and TSM. Both Chl-a and TSM data were comprised of composites of all the individual events of the four tidal stages, namely, spring to neap (SN), neap (N), neap to spring (NS), and spring (S) tides. The coastal area with a total number of daily data less than 450 (50.61%) was masked to avoid the bias of the low data in the analysis of the spatial and temporal variability of Chl-a and TSM [1]. This also ensures that the quality of data used for reconstruction meets the requirement of DINEOF, which requires the missing pixel time frequency to be less than 95% [7,8].

In Situ Data and Pre-Processing
Chl-a and other environmental data, including the sea level amplitude, wind speed, precipitation, SST, and river discharge, were used for validation of the DINEOF reconstruction and for the forcing factors of Chl-a and TSM, respectively. The in situ Chl-a data were collected from Nagoya University (2002-2017), Nagasaki University (2002-2011), and the three Fisheries Research institutes located around Ariake Sea (2011)(2012)(2013)(2014)(2015). In situ Chl-a data from the Nagoya and Nagasaki universities were measured using a pre-calibrated Turner Designs Fluorometer 10-AU following the method of Welschmeyer [27], while the Chl-a data from the Fisheries Research institutes were measured using the fluorometric method with 90% acetone extraction [28]. The spatial distribution of the data is shown in Figure 1 [25].
The tide data from the Oura observation site ( Figure 1) were downloaded from the Japan Oceanographic Data Center (https://www.jodc.go.jp/jodcweb/JDOSS/index_j.html (accessed on 21 December 2021)) (July 2002 to December 2017) and the Japan Meteorological Agency (https://www.data.jma.go.jp/gmd/kaiyou/db/tide/suisan/index.php (accessed on 21 December 2021)) (January 2011 to December 2017). In Yang et al. (2020), the time period of each tidal stage from July 2002 to December 2017 was separated by the tide data from the Japan Oceanographic Data Center (July 2002 to December 2017) and from the Japan Meteorological Agency (January 2011 to December 2017) as the data sources did not cover the entire time period. However, the Japan Oceanographic Data Center has extended its database and now provides data from January 2011 to December 2017. Therefore, the tide data during this time period were also downloaded to validate the consistency of the tide data from the two data sources. The wind speed, precipitation, and SST from January 2002 to December 2017 were also downloaded from the Japan Meteorological Agency. The observation sites were Shiraishi (No. 1), Shimabara (No. All of the in situ Chl-a data from the same location during the same tidal stage were averaged to compare them with the MODIS Chl-a from the corresponding tidal stage. The in situ Chl-a data selected were from the time period of 10:00 to 17:00 because MODIS usually passes over the Ariake Sea at local times from 13:00 to 14:00. Therefore, we chose in situ data collected within 3 h of the MODIS passing to better control the quality of the in situ MODIS matchups.
Daily tide data from the Japan Oceanographic Data Center were used to calculate the daily sea level amplitude, which is the difference between the maximum and minimum tidal levels. The daily sea level amplitude was then averaged for each month of each year. Furthermore, average data from all years were again averaged to obtain the monthly sea level amplitude data.
For the environmental variables, the average was calculated for all of the observation sites to obtain the daily wind speed and SST data, while the sum was calculated to obtain the daily precipitation and river discharge data. Then, monthly data were obtained in the same way as that for the monthly sea level amplitude data.

DINEOF
DINEOF was used to reconstruct missing data in the regionally refined MODISderived and improved MODIS Chl-a and TSM datasets (July 2002-December 2017), in order to extract the main patterns of the spatial and temporal variability of Chl-a and TSM. We used the DINEOF 3.0 package [7,8], which can be downloaded from the website of GeoHydrodynamics and Environment Research (GHER). The main procedure includes the following steps: (1) Each data set is stored in a 3-D matrix (y × x × t) with the condition y × x > t. The y and x represent the dimensions of the latitude and longitude of each image, respectively, and the t represents the total number of images. Note that the natural logarithm of Chl-a and TSM data were used as DINEOF inputs to avoid negative values in the reconstructed data; (2) the mean value of all the data over space and time is subtracted from the matrix, and the missing data are set to zero to avoid a biased initial guess; (3) iterative singular value decomposition [29] and cross-validation of 3% of the valid data (randomly selected) are performed to compute the optimal EOF modes; (4) the optimal number of EOF modes is used to reconstruct the whole dataset. A detailed description of the DINEOF method can be referred to in [7,8].
To validate the accuracy of the DINEOF reconstruction method, we randomly selected 1% of the valid pixels from the original datasets of Chl-a and TSM and intentionally regarded them as 'missing values'. Additionally, we chose to retain the values in the valid pixels of the original datasets, so the reconstructed data only involved the invalid pixels. After the DINEOF reconstruction, the reconstructed data were compared with the original data for 1% of the randomly selected data set.
The composite temporal functions, i.e., principal components (PCs), were constructed to investigate the monthly variability of Chl-a and TSM. The procedure is as follows: (1) Average the temporal functions of the same tidal stages to obtain the monthly springneap tidal cycle of temporal functions for each year; (2) average the temporal functions in (1) of the same year to obtain the monthly spring-neap tidal cycle of temporal functions; and (3) average the temporal functions in (2) of the same month to obtain the monthly temporal functions. As a result, 147 and 146 monthly temporal functions were obtained for Chl-a and TSM for 2002-2017, respectively.

Statistical Analysis
The correlations among the Chl-a, TSM, sea level amplitude, wind speed, precipitation, river discharge, and SST were quantified by calculating the Pearson's correlation coefficients (i.e., r) [30]. Then, the significance of the correlations was estimated by the student's t test to evaluate the hypothesis of no correlation. For monthly and individual events of spring-neap tidal cycle data, the degrees of freedom were 10 and ≤54 (14 tidal cycles × 4 tidal stages), respectively. Under these numbers, the absolute r values larger than 0.57 and 0.27, respectively, were significantly different from zero (5% significance level or p < 0.05) [30]. In addition, regressions between two variables and the statistical errors of the regressions in terms of R 2 , the root mean square error (RMSE), and bias, were also evaluated.
The procedure of the methods is described in the flow chart ( Figure 2) below.
Remote Sens. 2021, 13, x FOR PEER REVIEW 7 of 21 cycles × 4 tidal stages), respectively. Under these numbers, the absolute r values larger than 0.57 and 0.27, respectively, were significantly different from zero (5% significance level or p < 0.05) [30]. In addition, regressions between two variables and the statistical errors of the regressions in terms of R 2 , the root mean square error (RMSE), and bias, were also evaluated. The procedure of the methods is described in the flow chart ( Figure 2) below.

Figure 2.
Flow chart of the procedure of the methods.

DINEOF reconstruction
Although the spatial coverage of the original data in each image was more than 20%, the spatial coverage in the Ariake Sea area may be less than 5%, which is the minimum value required for DINEOF [7,8]. In such cases, we discarded the original Chl-a and TSM images. Eventually, there were 547 and 564 images for Chl-a and TSM, respectively. The dimensions, percent of missing data, number of retained EOFs, and total variance of the retained EOF modes are shown in Table 2. The percent of missing data and total variance of the retained EOF modes were similar for both Chl-a and TSM, but the number of retained EOFs for Chl-a was much smaller than for TSM. The number (547) of Chl-a images was also much smaller, which may have caused the difference in the number of retained EOFs between Chl-a and TSM. An example of the comparison between reconstructed and original Chl-a and TSM data during the spring-neap tidal cycle is shown in Figure 3. For the original Chl-a (b) and TSM (d), there were many spatial gaps during the spring-neap tidal cycle in the Ariake Sea, especially for TSM. In contrast, the reconstructed Chl-a (a) and TSM (c) data were

DINEOF reconstruction
Although the spatial coverage of the original data in each image was more than 20%, the spatial coverage in the Ariake Sea area may be less than 5%, which is the minimum value required for DINEOF [7,8]. In such cases, we discarded the original Chl-a and TSM images. Eventually, there were 547 and 564 images for Chl-a and TSM, respectively. The dimensions, percent of missing data, number of retained EOFs, and total variance of the retained EOF modes are shown in Table 2. The percent of missing data and total variance of the retained EOF modes were similar for both Chl-a and TSM, but the number of retained EOFs for Chl-a was much smaller than for TSM. The number (547) of Chl-a images was also much smaller, which may have caused the difference in the number of retained EOFs between Chl-a and TSM. An example of the comparison between reconstructed and original Chl-a and TSM data during the spring-neap tidal cycle is shown in Figure 3. For the original Chl-a (b) and TSM (d), there were many spatial gaps during the spring-neap tidal cycle in the Ariake Sea, especially for TSM. In contrast, the reconstructed Chl-a (a) and TSM (c) data were smooth and showed a reasonable spatial distribution, i.e., the reconstructed Chl-a and TSM concentrations were generally higher near the coastal areas than in the offshore areas. smooth and showed a reasonable spatial distribution, i.e., the reconstructed Chl-a and TSM concentrations were generally higher near the coastal areas than in the offshore areas.

Validation of Reconstructed Data
In addition to the DINEOF internal cross-validation, which was designed to obtain the optimal EOF modes for final data reconstruction, we also randomly selected 1% (depends on the size of the original datasets) of the valid pixel values from the original data set as 'missing values' to compare the reconstructed and original data in these pixels. The comparisons for Chl-a and TSM are shown in the density plots in Figure 4. For both Chla and TSM, the reconstructed and original data showed good correlations in terms of the slope, coefficient of determination (R 2 ), bias, and root mean square error (RMSE). Moreover, the density, i.e., number of data in each grid bin (1 km × 1 km), increased towards the 1:1 line. Therefore, the data reconstructed by the DINEOF method were accurate and reliable. Validation of the corresponding original and reconstructed Chl-a data with in situ data ( Figure 5) yielded statistically significant (p < 0.05) correlations. Moreover, the valid and invalid data points almost overlapped, which further confirmed that the reconstructed data were close to the original data. It should be noted that the original data were validated with in situ data in Yang et al. [25].

Validation of Reconstructed Data
In addition to the DINEOF internal cross-validation, which was designed to obtain the optimal EOF modes for final data reconstruction, we also randomly selected 1% (depends on the size of the original datasets) of the valid pixel values from the original data set as 'missing values' to compare the reconstructed and original data in these pixels. The comparisons for Chl-a and TSM are shown in the density plots in Figure 4. For both Chl-a and TSM, the reconstructed and original data showed good correlations in terms of the slope, coefficient of determination (R 2 ), bias, and root mean square error (RMSE). Moreover, the density, i.e., number of data in each grid bin (1 km × 1 km), increased towards the 1:1 line. Therefore, the data reconstructed by the DINEOF method were accurate and reliable. smooth and showed a reasonable spatial distribution, i.e., the reconstructed Chl-a and TSM concentrations were generally higher near the coastal areas than in the offshore areas.

Validation of Reconstructed Data
In addition to the DINEOF internal cross-validation, which was designed to obtain the optimal EOF modes for final data reconstruction, we also randomly selected 1% (depends on the size of the original datasets) of the valid pixel values from the original data set as 'missing values' to compare the reconstructed and original data in these pixels. The comparisons for Chl-a and TSM are shown in the density plots in Figure 4. For both Chla and TSM, the reconstructed and original data showed good correlations in terms of the slope, coefficient of determination (R 2 ), bias, and root mean square error (RMSE). Moreover, the density, i.e., number of data in each grid bin (1 km × 1 km), increased towards the 1:1 line. Therefore, the data reconstructed by the DINEOF method were accurate and reliable. Validation of the corresponding original and reconstructed Chl-a data with in situ data ( Figure 5) yielded statistically significant (p < 0.05) correlations. Moreover, the valid and invalid data points almost overlapped, which further confirmed that the reconstructed data were close to the original data. It should be noted that the original data were validated with in situ data in Yang et al. [25]. Validation of the corresponding original and reconstructed Chl-a data with in situ data ( Figure 5) yielded statistically significant (p < 0.05) correlations. Moreover, the valid and invalid data points almost overlapped, which further confirmed that the reconstructed data were close to the original data. It should be noted that the original data were validated with in situ data in Yang et al. [25].

EOF Analysis of Chl-a and TSM
The optimal EOF modes from the DINEOF reconstruction process were 36 and 49 for Chl-a and TSM, respectively (Table 2). Generally, the major spatial and temporal variability of geophysical variables is revealed in the first three EOF modes. Therefore, our EOF analysis of Chl-a and TSM is based on the first three EOF modes, which explained 75.9% and 75.5% of the total variance in the reconstructed Chl-a and TSM datasets, respectively. The three spatial EOF modes of Chl-a and TSM are shown in Figure 6a,b.

EOF Analysis of Chl-a and TSM
The optimal EOF modes from the DINEOF reconstruction process were 36 and 49 for Chl-a and TSM, respectively (Table 2). Generally, the major spatial and temporal variability of geophysical variables is revealed in the first three EOF modes. Therefore, our EOF analysis of Chl-a and TSM is based on the first three EOF modes, which explained 75.9% and 75.5% of the total variance in the reconstructed Chl-a and TSM datasets, respectively. The three spatial EOF modes of Chl-a and TSM are shown in Figure 6a

EOF Analysis of Chl-a and TSM
The optimal EOF modes from the DINEOF reconstruction process were 36 and 49 for Chl-a and TSM, respectively (Table 2). Generally, the major spatial and temporal variability of geophysical variables is revealed in the first three EOF modes. Therefore, our EOF analysis of Chl-a and TSM is based on the first three EOF modes, which explained 75.9% and 75.5% of the total variance in the reconstructed Chl-a and TSM datasets, respectively. The three spatial EOF modes of Chl-a and TSM are shown in Figure 6a,b.  The time-series of the monthly temporal functions (i.e., PCs) of Chl-a and TSM and of the monthly environmental factors, including the sea level amplitude, wind speed, precipitation, river discharge, and sea surface temperature, are shown in Figures 7a-f and 8a- The time-series of the monthly temporal functions (i.e., PCs) of Chl-a and TSM and of the monthly environmental factors, including the sea level amplitude, wind speed, precipitation, river discharge, and sea surface temperature, are shown in Figure 7a    The time-series of the monthly temporal functions (i.e., PCs) of Chl-a and TSM and of the monthly environmental factors, including the sea level amplitude, wind speed, precipitation, river discharge, and sea surface temperature, are shown in Figure 7a

Seasonal Variability of Chl-a, TSM, and Factors
For the first EOF mode of Chl-a, the spatial EOF (Figure 6a), which explained 53.6% of the total variance, generally showed negative values for the whole Ariake Sea. The variance was slightly larger off Saga, in Isahaya Bay, and off Kumamoto than that in other areas; the corresponding temporal EOF (Figure 6a) generally decreased from winter to summer (July), and then increased in autumn. Therefore, the first EOF mode of Chl-a increased from winter to summer (peaked in July), and then decreased in autumn, for the whole Ariake Sea.
For the second EOF mode of Chl-a, the spatial EOF (17.8%) (Figure 6a) showed all positive values, and the variance was homogeneous; the corresponding temporal EOF (Figure 7b) generally decreased from winter to spring (December to April), and then increased in summer (July), before decreasing again in autumn (October). Therefore, the second EOF mode of Chl-a showed a similar temporal variability to that of the temporal EOF.
For the third EOF mode of Chl-a, the spatial EOF (4.5%) (Figure 6a) showed opposing patterns between the region off Saga and Kumamoto, Isahaya Bay, and the southern part of Ariake Sea; the corresponding temporal EOF generally increased slightly from winter to summer (June), and then slightly decreased in autumn. Therefore, the third EOF mode of Chl-a slightly decreased from winter to summer (June), and then slightly increased in autumn for the region off Saga and Kumamoto, and Isahaya Bay, whereas it exhibited an inverse temporal variability in the southern part of the Ariake Sea.
For the first EOF mode of TSM, the spatial EOF (63.1%) (Figure 6b) displayed opposing patterns between the region off Saga and Kumamoto, Isahaya Bay, and the southern part of the Ariake Sea. Moreover, the variance off Saga was much larger than that in the other regions of the Ariake Sea. The corresponding temporal EOF (Figure 7d) generally increased from winter to spring, and then decreased in autumn. Therefore, the first EOF mode of TSM generally increased from winter to spring, and then decreased to autumn for off Saga, Isahaya Bay, and off Kumamoto, whereas the temporal variability was the inverse of this for the southern part of the Ariake Sea.
For the second EOF mode of TSM, the spatial EOF (9.7%) (Figure 6b) showed all negative values, and the variance was also homogenous; the corresponding temporal EOF (Figure 7e) generally increased from winter to spring, and then decreased to autumn. Therefore, the second EOF mode of TSM showed an inverse temporal variability to that of the temporal EOF for the whole Ariake Sea.
For the third EOF mode of TSM, the spatial EOF (2.7%) showed opposing patterns between the region off Saga and the other areas of Ariake Sea; the corresponding temporal EOF generally increased slightly from winter to summer (June), and then decreased in autumn. Therefore, the monthly TSM of the third EOF mode exhibited a similar temporal variability to that of the temporal EOF for the whole Ariake Sea, except for the region off Saga.
Furthermore, the environmental factors also displayed seasonal variability (Figure 8). The monthly sea level amplitude (Figure 8a) generally decreased from winter to spring (March), increased in summer (August), and then decreased in autumn. The monthly wind speed (Figure 8b) increased in winter, decreased in spring, and showed more variability in summer and autumn. The monthly precipitation (Figure 8c), river discharge (Figure 8d), and sea surface temperature (Figure 8e) generally increased from winter to summer, and then decreased in autumn. Moreover, the precipitation peaked in June, whereas the river discharge and sea surface temperature peaked in July.

Correlations among Monthly Chl-a, TSM, and the Environmental Factors
Pearson's correlation coefficients (r) were calculated between the monthly PCs of Chl-a and the monthly PCs of TSM and the environmental factors (Table 3), and between the monthly PCs of TSM and the environmental factors (Table 4). Significant (p < 0.05) correlations are highlighted in the tables. Additionally, the r values were calculated for correlations between the sea level amplitude, wind speed, precipitation, sea surface temperature, and river discharge. Significant correlations were found between precipitation and river discharge (r = 0.93), precipitation and SST (r = 0.88), and river discharge and SST (r = 0.83).  Table 3 shows the correlations between the monthly PCs of Chl-a, TSM, and environmental factors. The monthly PC1 of Chl-a was significantly correlated with the monthly PC3 of TSM, precipitation, SST, and river discharge, and the r values were higher for the latter three factors. The monthly PC2 of Chl-a was significantly correlated with the monthly PC1 and PC2 of TSM, sea level amplitude, and wind speed, and the r value was higher for the monthly PC2 of TSM. The monthly PC3 of Chl-a was significantly correlated with the monthly PC3 of TSM, precipitation, SST, and river discharge, and the r value was higher for the monthly PC3 of TSM. Table 4 shows the correlations between the monthly PCs of TSM and environmental factors. The monthly PC1 of TSM was not correlated with any environmental factors. The monthly PC2 of TSM was significantly correlated with the monthly sea level amplitude and wind speed, with a higher r value for the wind speed. The monthly PC3 of TSM was significantly correlated with the monthly sea level amplitude, precipitation, SST, and river discharge.

Individual Events of the Spring-Neap Tidal Variability of Chl-a and TSM
The fourteen individual events of spring-neap tidal cycles (Table 5) with continuous tidal stage data and corresponding Chl-a and TSM data were selected from 2002 to 2017. There were four, two, two, and six individual events of spring-neap tidal cycles from winter, spring, summer, and autumn, respectively. Table 5. Statistics of the fourteen selected individual events of spring-neap tidal cycles. SN, N, NS, and S represent spring to neap, neap, neap to spring, and spring tides, respectively. The fourteen tidal cycles were separated by the four seasons, i.e., winter (December-February), spring (March-May), summer (June-August), and autumn (September-November) (Date: M/D/Y). The time series of the PCs of Chl-a and TSM and the environmental factors are shown in Figures 9 and 10, respectively. The PCs of Chl-a showed no clear patterns during the spring-neap tidal cycle, whereas the PC1 and PC2 showed a clear seasonal pattern, where Chl-a was generally higher in summer than in winter (Figure 9a-c). In contrast, the PC1 and PC2 of TSM exhibited a clear pattern indicating that TSM was generally higher during the NS and S tides than during the SN and N tides, except for the PC3 of TSM, whereas the PCs of TSM showed no clear seasonal patterns (Figure 9d-f).

Correlations among Individual Events of Spring-Neap Tidal Chl-a, TSM, and Environmental Factors
Pearson's correlation coefficients (r) were also calculated between the individual events of spring-neap tidal PCs of Chl-a and the environmental factors (Table 6), and between the individual events of spring-neap tidal PCs of TSM and the environmental factors ( Table 7). The significant r values are highlighted in the tables.
As shown in Table 6, the PC1 of Chl-a was significantly correlated with SST. The PC2 of Chl-a was significantly correlated with the PC1 and PC3 of TSM and SST, and a higher r value was obtained for the PC1 of TSM. The PC3 of Chl-a was significantly correlated with the PC3 of TSM and precipitation, and a higher r value was obtained for the precipitation.
As demonstrated in Table 7, both the PC1 and PC2 of TSM were significantly correlated with the sea level amplitude. The PC3 of TSM was significantly correlated with the precipitation and river discharge. The time series of the PCs of Chl-a and TSM and the environmental factors are shown in Figures 9 and 10, respectively. The PCs of Chl-a showed no clear patterns during the spring-neap tidal cycle, whereas the PC1 and PC2 showed a clear seasonal pattern, where Chl-a was generally higher in summer than in winter (Figure 9a-c). In contrast, the PC1 and PC2 of TSM exhibited a clear pattern indicating that TSM was generally higher during the NS and S tides than during the SN and N tides, except for the PC3 of TSM, whereas the PCs of TSM showed no clear seasonal patterns (Figure 9d-f).

Correlations among Individual Events of Spring-Neap Tidal Chl-a, TSM, and Environmental Factors
Pearson's correlation coefficients (r) were also calculated between the individual events of spring-neap tidal PCs of Chl-a and the environmental factors (Table 6), and between the individual events of spring-neap tidal PCs of TSM and the environmental factors ( Table 7). The significant r values are highlighted in the tables.
As shown in Table 6, the PC1 of Chl-a was significantly correlated with SST. The PC2 of Chl-a was significantly correlated with the PC1 and PC3 of TSM and SST, and a higher r value was obtained for the PC1 of TSM. The PC3 of Chl-a was significantly correlated with the PC3 of TSM and precipitation, and a higher r value was obtained for the precipitation. As demonstrated in Table 7, both the PC1 and PC2 of TSM were significantly correlated with the sea level amplitude. The PC3 of TSM was significantly correlated with the precipitation and river discharge.   Furthermore, the significant correlations in Tables 6 and 7 were analyzed with the scatter plots in Figures 11 and 12, respectively. In Figure 11, the PC1 of Chl-a was negatively correlated with the SST and the values were generally lower in summer. The PC2 of Chl-a was negatively correlated with the PC1 of TSM, whereas it was positively correlated with the PC3 of TSM and the SST. The correlation between the PC2 of Chl-a and the PC1 of TSM was stronger than the other two correlations in terms of r and R 2 . The PC2 of Chl-a was generally lower in winter than that in other seasons. The PC3 of Chl-a was positively correlated with the PC3 of TSM and the precipitation. The correlations were weak due to the small variability of the PC3 of Chl-a values.
tively correlated with the SST and the values were generally lower in summer. The PC2 of Chl-a was negatively correlated with the PC1 of TSM, whereas it was positively correlated with the PC3 of TSM and the SST. The correlation between the PC2 of Chl-a and the PC1 of TSM was stronger than the other two correlations in terms of r and R 2 . The PC2 of Chl-a was generally lower in winter than that in other seasons. The PC3 of Chl-a was positively correlated with the PC3 of TSM and the precipitation. The correlations were weak due to the small variability of the PC3 of Chl-a values. As shown in Figure 12, the PC1 of TSM was positively correlated with the sea level amplitude and the values were generally lower during SN and N tides than during NS and S tides. The PC2 of TSM was negatively correlated with the sea level amplitude and the values were also generally lower during SN and N tides than during NS and S tides. The PC3 of TSM was positively correlated with the precipitation and river discharge, whereas the data were scattered.  Table 6, i.e., PC1 of TSM vs. sea level amplitude (a), PC2 of TSM vs. sea level amplitude (b), PC3 of TSM vs. precipitation (c), and PC3 of TSM vs. river discharge (d). The data for winter, spring, summer, and autumn are represented by filled circle, diamond, triangle, and plus symbols, and the data for spring to neap (SN), neap (N), neap to spring (NS), and spring (S) tides are represented by red, green, blue, and magenta colors. The black lines in each plot are regression lines.

DINEOF Reconstruction of MODIS Chl-a and TSM for the Ariake Sea
There are often spatial gaps in satellite ocean color data due to various issues, such as high sensor-zenith angles, high sun glint contamination, and significant cloud cover [2,31]. Filling the missing pixels before using the data is highly advantageous, especially in small water bodies like the Ariake Sea, in order to better monitor coastal eutrophication, harmful algal blooms, and sediment plumes. The DINEOF method has been extensively used to fill gaps in ocean color data, including Chl-a and TSM, for various satellite sensors, such as MODIS, Geostationary Ocean Color Imager (GOCI), and merged Visible Infrared Imaging Radiometer Suite (VIIRS) Suomi National Polar-orbiting Partnership (SNPP)/  Table 6, i.e., PC1 of TSM vs. sea level amplitude (a), PC2 of TSM vs. sea level amplitude (b), PC3 of TSM vs. precipitation (c), and PC3 of TSM vs. river discharge (d). The data for winter, spring, summer, and autumn are represented by filled circle, diamond, triangle, and plus symbols, and the data for spring to neap (SN), neap (N), neap to spring (NS), and spring (S) tides are represented by red, green, blue, and magenta colors. The black lines in each plot are regression lines.
As shown in Figure 12, the PC1 of TSM was positively correlated with the sea level amplitude and the values were generally lower during SN and N tides than during NS and S tides. The PC2 of TSM was negatively correlated with the sea level amplitude and the values were also generally lower during SN and N tides than during NS and S tides. The PC3 of TSM was positively correlated with the precipitation and river discharge, whereas the data were scattered.

DINEOF Reconstruction of MODIS Chl-a and TSM for the Ariake Sea
There are often spatial gaps in satellite ocean color data due to various issues, such as high sensor-zenith angles, high sun glint contamination, and significant cloud cover [2,31]. Filling the missing pixels before using the data is highly advantageous, especially in small water bodies like the Ariake Sea, in order to better monitor coastal eutrophication, harmful algal blooms, and sediment plumes. The DINEOF method has been extensively used to fill gaps in ocean color data, including Chl-a and TSM, for various satellite sensors, such as MODIS, Geostationary Ocean Color Imager (GOCI), and merged Visible Infrared Imaging Radiometer Suite (VIIRS) Suomi National Polar-orbiting Partnership (SNPP)/ National Oceanic and Atmospheric Administration (NOAA)-20 [3,[11][12][13][14][15]31].
A previous study [1] on the spring-neap tidal variability of Chl-a and TSM in the Ariake Sea reported that a significant amount of spatial data were missing during the individual tidal events, which hampered a detailed analysis of the spatial and temporal variability of Chl-a and TSM. This was a big issue for the detailed analysis of spatial and temporal variability of Chl-a and TSM during the individual events of the spring-neap tidal cycles. Moreover, missing data also affect the validity of the climatology data, and thus the interpretation of the long-term Chl-a and TSM variability.
In spite of its importance in the Ariake Sea, the major spatial and temporal variability of Chl-a and TSM and the mechanisms behind it have not been comprehensively investigated. Using the EOF method, which has been proved to be an effective tool for capturing the major spatial and temporal variability in a large geophysical dataset [9,32,33], we could gather considerable information on the variability of Chl-a and TSM for the Ariake Sea. Applying the DINEOF method to reconstruct the spatial gaps in the individual events of spring-neap tidal stages over a sixteen-year period (2002-2017), MODIS Chl-a and TSM produced generally smooth and reasonable data compared to the original images.
Validation of the reconstructed data demonstrated that the reconstructed data were accurate and reliable and could be used for EOF analysis of the Chl-a and TSM variability in the Ariake Sea. This allowed us to examine the correlations between Chl-a and TSM, sea level amplitude, wind speed, precipitation, SST, and river discharge, and between TSM and the other factors on the temporal scales of monthly and spring-neap tidal cycles.

EOF Analysis: Monthly Variability of Chl-a
The first EOF mode of Chl-a showed that the variance of Chl-a was higher off Saga, in Isahaya Bay, and off Kumamoto than in the middle and southern regions of the Ariake Sea, and that Chl-a increased from winter to summer and then decreased to autumn. The correlation analysis demonstrated that the temporal variability of Chl-a was positively correlated with precipitation, SST, and river discharge. This positive correlation between Chl-a and precipitation and river discharge was consistent with the previous studies. For example, Yang et al. [1] reported that the monthly Chl-a based on the regionally-tuned MODIS/Aqua Chl-a dataset during 2002 to 2017 was strongly positively correlated with the monthly river discharge, and Ishizaka et al. [18] reported that the monthly SeaWiFS Chl-a from May 1998 to December 2001 showed two peaks after the peaks of rain and river discharge, for the Ariake Sea. Moreover, Le et al. [34] reported that the long-term Chl-a time series of SeaWiFS and MODIS observations from 1998 to 2011 showed that the Chl-a was higher and lower in the wet and dry season, respectively, and that river discharge could explain~60% of the seasonal changes. The main reason for this is that increased precipitation enhanced the river discharge, which supplies large amounts of nutrients to the seawater.
The positive correlation between SST and Chl-a may be explained by the fact that increased SST enhances the water stratification, which in turn increases the light availability in turbid coastal waters, and thus promotes the growth of phytoplankton [19]. Similar results were obtained by Ji et al. [35] for the coastal waters of East China Sea (north coast of Jiangsu Province, Yangtze River estuary region, and the Hangzhou Bay) using the DINEOF reconstructed monthly MODIS/Aqua Chl-a and SST datasets during 2003-2016. In the Bay of Bengal, where phytoplankton productivity is highly influenced by river discharge, Jutla et al. [36] reported positive correlations between monthly SeaWiFS Chl-a and SST from 1997 to 2009. In this study, the first EOF mode of Chl-a was positively correlated with the third EOF mode of TSM for all of the Ariake Sea except for the region off Saga, indicating that the spatial and temporal variability of Chl-a of the first EOF mode was controlled by precipitation, SST, and river discharge.
The spatial variance of the second EOF mode of Chl-a was small, and similar to that of the second EOF mode of TSM. Additionally, the second EOF mode of Chl-a generally decreased from winter to spring, increased in summer, and then decreased to autumn. This temporal variability was very similar to that of the second EOF mode of TSM and sea level amplitude. The second EOF mode of TSM was also strongly and positively correlated with the sea level amplitude. It is possible that both the second EOF mode of Chl-a and TSM was mainly controlled by the sea level amplitude. As far as we know, there have been few studies investigating the impact of sea level amplitude on the monthly Chl-a and/or TSM. Yang et al. [1] reported that monthly Chl-a generally increased during N and NS tides when the sea level amplitude was low and that monthly TSM generally increased during NS and S tides when the sea level amplitude was high based on regionally-tuned MODIS/Aqua Chl-a and TSM datasets from 2002 to 2017. However, the seasonal variability of the sea level amplitude and the response of the monthly Chl-a and TSM to it was not investigated.
For the third EOF mode of Chl-a, the temporal variability of Chl-a was the inverse between off Saga, Isahaya Bay, and off Kumamoto, and the southern part of the Ariake Sea. In other words, Chl-a decreased from winter to summer, and then increased to autumn, for off Saga, Isahaya Bay, and off Kumamoto, and the Chl-a increased from winter to summer, and then decreased to autumn, for the southern part of the Ariake Sea. Furthermore, the third EOF mode of Chl-a was positively correlated with the third EOF mode of TSM for off Saga and the southern part of the Ariake Sea, whereas it was negatively correlated for Isahaya Bay and off Kumamoto. The third EOF mode of Chl-a was also positively correlated with precipitation, SST, and river discharge, which increased from winter to summer and then decreased to autumn. It was noticed that the third EOF mode of TSM also decreased from winter to summer, and then increased to autumn, for off Saga, whereas it increased from winter to summer, and then decreased to autumn, for the other areas of the Ariake Sea. Therefore, the temporal variability of Chl-a and precipitation, SST, and river discharge was the same for the southern part of the Ariake Sea, whereas it was the inverse for off Saga, Isahaya Bay, and off Kumamoto. In addition, the temporal variability of TSM and precipitation, SST, and river discharge was the same for the whole Ariake Sea, except for off Saga.
Previously, Yang et al.
[1] reported that Chl-a was negatively correlated with TSM during the spring-neap tidal cycle for all seasons except summer off Saga, whereas TSM was possibly dominated by Chl-a during the tidal cycle for all seasons in Isahaya Bay and off Kumamoto. In this study, the correlations between the monthly PC3 of Chl-a and the monthly PC3 of TSM for off Saga, Isahaya Bay, and off Kumamoto were the inverse of that in Yang et al. [1]. The possible reason for this may be that the monthly PCs of Chl-a and TSM were the average values during the spring-neap tidal cycle or that the similar regional correlations of Chl-a and TSM were represented in the higher EOF modes.

EOF Analysis: Monthly Variability of TSM
The first EOF mode of TSM showed the general distribution of TSM in the Ariake Sea [1], with a higher variance of TSM in coastal areas, including off Saga, Isahaya Bay, and off Kumamoto, which may be due to the larger spring-neap tidal range off Saga and the higher river discharges in coastal areas [16,17]. The time series of TSM was only negatively correlated with that of the second EOF mode of Chl-a, suggesting that TSM and Chl-a are negatively correlated off Saga, in Isahaya Bay, and off Kumamoto. However, the second EOF mode of Chl-a cannot be the controller for the seasonal variability of the first EOF mode of TSM due to the negative correlation of Chl-a and TSM.
For the second EOF mode of TSM, the homogeneous spatial distribution of the variance was probably caused by wind. However, this cannot be proven because high resolution satellite or modeled wind data are not available for the Ariake Sea. The temporal variability of the second EOF mode of TSM resembled that of the sea level amplitude. It is possible that the second EOF mode of TSM was controlled by the combined force of wind and sea level amplitude. This is the first study that demonstrates a correlation between the monthly TSM and sea level amplitude for the Ariake Sea.
For the third EOF mode of TSM, the temporal variability was very similar to that of precipitation, SST, and river discharge, for the whole Ariake Sea, except for off Saga. It is possible that the third EOF mode of TSM was mainly controlled by precipitation, SST, and river discharge, for the whole Ariake Sea, except for off Saga, where the controller was not clear.

EOF Analysis: Spring-Neap Tidal Variability of Chl-a
The correlations between Chl-a and environmental factors, and between TSM and environmental factors, during the spring-neap tidal cycles were firstly quantitatively calculated. Previously, Yang et al.
[1] compared the spatially-averaged Chl-a and TSM during the spring-neap tidal cycle using the monthly and individual events of spring-neap tidal cycle data. However, the correlation between Chl-a and TSM was only qualitatively analyzed, and the impact factors only included TSM and river discharge. In this study, we also considered the sea level amplitude, wind speed, precipitation, and SST. Moreover, we calculated the correlation coefficients among the PCs of Chl-a and TSM, and the other impact factors for the fourteen selected individual events of spring-neap tidal cycle data.
The first EOF mode of Chl-a was positively correlated with SST for the whole Ariake Sea. The scatter plot of the PC1 of Chl-a versus SST indicated that Chl-a was generally higher in summer when SST was higher, and that the spring-neap tidal variability of the first EOF mode of Chl-a was not clear.
For the second EOF mode of Chl-a, it was negatively correlated with the first EOF mode of TSM for off Saga, Isahaya Bay, and off Kumamoto, whereas they were positively correlated for the southern part of the Ariake Sea. The first EOF mode of TSM was generally much higher during NS and S tides than that during SN and N tides for all the seasons except for summer, whereas the second EOF mode of Chl-a showed no clear patterns during the spring-neap tidal cycle, but it was generally lower in winter than in other seasons due to the lower temperature in winter.
For the third EOF mode of Chl-a, it was significantly correlated with the third EOF mode of TSM and precipitation. However, the correlations were weak because the range of the PC3 of Chl-a was narrow after excluding the two lowest and highest values. In addition, the seasonal and spring-neap tidal variability of the third EOF mode of Chl-a was not clear.

EOF Analysis: Spring-Neap Tidal Variability of TSM
For the first EOF mode of TSM, it was positively correlated with the sea level amplitude for off Saga, Isahaya Bay, and off Kumamoto, whereas it was negatively correlated for the southern part of the Ariake Sea. The scatter plot of the PC1 of TSM versus the sea level amplitude (Figure 11a) indicated that TSM was generally much higher during NS and S tides, when the sea level amplitude was higher than that during the SN and N tides. Yang et al.
[1] reported a similar variability of TSM during the spring-neap tidal cycle. This is due to the increased tidal resuspended sediment during the NS and S tides [37].
The second EOF mode of TSM was also positively correlated with the sea level amplitude (Table 6) for the whole Ariake Sea, indicating that the second EOF mode of TSM was generally higher during NS and S tides, when the sea level amplitude was much higher than that during SN and N tides.
The third EOF mode of TSM was positively correlated with precipitation and river discharge. The positive correlation of TSM and river discharge was supported by the fact that river discharge supplies a large amount of nutrients and suspended sediment to the Ariake Sea [16]. However, the positive correlation of TSM and precipitation may be due to the co-correlation of precipitation and river discharge.

Conclusions
The DINEOF method was successfully applied to fill the spatial gaps in all individual events of spring-neap tidal stages of the regionally-tuned MODIS Chl-a and MODISderived TSM datasets from 2002 to 2017. Then, we conducted EOF analysis of the reconstructed data to capture the major spatial and temporal variability of Chl-a and TSM associated with the sea level amplitude, wind speed, precipitation, SST, and river discharge, for the Ariake Sea.
The first three EOF modes of Chl-a mainly revealed the seasonal variability of Chl-a. For the first EOF mode of Chl-a, it showed the general spatial distribution of Chl-a, and the variance was larger in the coastal areas. The monthly PC1 of Chl-a was positively correlated with the monthly precipitation, SST, and river discharge, and the Chl-a concentrations increased from winter to summer, and then decreased to autumn. For the second EOF mode of Chl-a, it showed a homogeneous spatial distribution similar to that of the second EOF mode of TSM, which may have been caused by the wind speed. Both second EOF modes of Chl-a and TSM were positively correlated with the sea level amplitude, and their seasonal variability was similar to that of the sea level amplitude. For the third EOF mode of Chl-a, it showed regional differences between off Saga, Isahaya Bay, and off Kumamoto, but the spatial and temporal variability were complex.
The first three EOF modes of TSM exhibited both seasonal and spring-neap tidal variability. The first EOF mode of TSM also demonstrated the general spatial distribution, and the variance was larger off Saga, in Isahaya Bay, and off Kumamoto. However, the temporal variability was not clear. For the second EOF mode of TSM, it had a similar spatial and temporal variability to that of the second EOF mode of Chl-a. The controller may be a combination of the sea level amplitude and wind speed. The third EOF mode of TSM displayed regional variability between off Saga and the other areas. The PC3 of TSM was mainly positively correlated with precipitation, SST, and river discharge. The third EOF mode of TSM may be controlled by these factors in the other areas, but the controller was not clear for off Saga.
Our study showed that the DINEOF-reconstructed Chl-a and TSM data were reliable for further analysis. Moreover, we clearly captured the major spatial and temporal variability of Chl-a and TSM for the Ariake Sea, and we comprehensively and quantitatively analyzed the linkages among the Chl-a, TSM, sea level amplitude, wind speed, precipitation, SST, and river discharge at monthly and spring-neap tidal cycle scales.