Seasonal to interannual variability of Chlorophyll-a and sea surface 1 temperature in the Yellow Sea using MODIS satellite datasets

Abstract. The spatial and temporal variability of Chlorophyll-a concentration (CHL) and sea surface temperature (SST) in the Yellow Sea (YS) were examined using Empirical Orthogonal Function (EOF) analysis, which was based on the monthly, cloud-free Data INterpolating Empirical Orthogonal Function (DINEOF) reconstruction datasets for 2003–2015. The variability and oscillation periods on an inter-annual timescale were also confirmed using the Morlet wavelet transform and wavelet coherence analyses. At a seasonal time scale, the CHL EOF1 mode was dominated by a seasonal cycle of a spring and a fall bloom, with a spatial distribution that was modified by the strong mixing of the water column of the Yellow Sea Cold Warm Mass (YSCWM) that facilitated nutrient delivery from the ocean bottom. The EOF2 mode was likely associated with a winter bloom in the southern region, where it was affected by the Yellow Sea Warm Current (YSWC) that moved from southeast to north in winter. The SST EOF1 explained 99 % of the variance in total variabilities, which was dominated by an obvious seasonal cycle (in response to net surface heat flux) that was inversely proportional to the water depth. At the inter-annual scale, the wavelet power spectrum and global power spectrum of CHL and SST showed significant similar periods of variations. The dominant periods for both spectra were 2–4 years during 2003–2015. A significant negative cross-correlation existed between CHL and SST, with the largest correlation coefficient at time lags of 4 months. The wavelet coherence further identified a negative relationship that was significant statistically between CHL and SST during 2008–2015, with periods of 1.5–3 years. These results provided insight into how CHL might vary with SST in the future.


Introduction
Chlorophyll-a concentrations (CHL), as an index of phytoplankton pigment, are considered an important indicator of eutrophication in marine ecosystems, which is a process that may affect human life (Smith, 2006;Werdell et al., 2009).Additionally, it can be used to analyze the comprehensive dynamics of phytoplankton biomass (Muller-Karger et al., 2005).On the other hand, sea surface temperature (SST) anomalies indicate stratification of the water column, which is related closely to light and to nutrient loads of CHL (He et al., 2010).Certain studies have reported the spatio-temporal variability and relationship between CHL and SST (Gregg et al., 2005;Behrenfeld et al., 2006;Boyce et al., 2010).In the open ocean, Wilson and Coles (2005) analyzed global scale relationships between CHL and the monthly SST.
Similarly, the spatio-temporal variability of regional CHL and SST in the South Atlantic Bight and the Mediterranean Sea have been investigated using long-term satellite datasets (Miles and He, 2010;Volpe et al., 2012).Gao et al. (2013) examined the spatio-temporal distribution of CHL that was associated with SST in the western South China Sea using the Sea-viewing Wide Field-of-View Sensor (SeaWiFS) and National Oceanic and Atmospheric Administration Advanced Very High Resolution Radiometer (AVHRR) data.For coastal waters, Li and He (2014) examined spatio-temporal distribution of CHL that was associated with SST in the Gulf of Ocean Sci.Discuss., doi:10.5194/os-2017Discuss., doi:10.5194/os- -11, 2017 Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.Maine (GOM) using daily MODIS data.Moradi and Kabiri (2015) examined the spatio-temporal variability of CHL and SST in the Persian Gulf using MODIS Level-2 products.These studies found region-specific relationships between climate-driven SST and CHL.These findings also indicated that knowledge of the spatio-temporal variability in CHL and SST can assist scientists in developing a more comprehensive perspective of biological and physical oceanography of marine ecosystems in the global scale.
The Yellow Sea (YS) has an average water depth of only 44 m, and it is marginal seas surrounded by China, and Korea (Fig. 1a).It responds quickly to atmospheric climate change and, in turn, the YS influences local climate variability as a result of the air-sea feedback process.The Yellow Sea Warm Current (YSWC) moves from southeast to north in winter (Fig. 1b) (Teague and Jacobs, 2000;Lie et al., 2009;Yu et al., 2010), and the Yellow Sea Cold Water Mass (YSCWM; 122-125° E, 33-37° N) is entrenched at the bottom in summer (Fig. 1c) (Zhang et al., 2008).These water masses represent the two most important physical oceanographic features in the YS.
In addition, a southward coastal flow is present in winter along the eastern and western sides of the YS, which corresponds to the northward YSWC in the central sea area (Wei et al., 2016;Xu et al., 2016).These features affect the physical properties, water mass and circulation in the YS, and they are complicated both spatially and temporally (Chu et al., 2005).
To date, the importance of SST variability and the associated features such as thermal or tidal fronts, coastal waters, and currents in the YS have been addressed by Ocean Sci.Discuss., doi:10.5194/os-2017Discuss., doi:10.5194/os- -11, 2017 Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.numerous satellite-based studies (Tseng et al., 2000;Lin et al., 2005;Wei et al., 2010;Yeh and Kim, 2010;Shi and Wang, 2012), and the long term CHL trends and seasonal variations have been studied as well (Shi and Wang, 2012;Yamaguchi, et al., 2012;Liu and Wang, 2013).In recent years, warming signals of SST in the YS were reported by Yeh and Kim (2010) and Park et al. (2015), but few researchers have paid attention to how has the increasing SST affected the spatio-temporal pattern of CHL in the YS?What is the region-specific relationship between climate-driven SST and CHL?
To answer these questions, we combined remote sensing datasets and statistical analysis to investigate the patterns of variability of CHL and SST over seasonal and inter-annual periods at temporal scales during 2003-2015 in the YS.The present work provides a comprehensive description of the phytoplankton biomass and the physical conditions using 13 years of satellite-derived datasets.The objectives of the study were (i) to identify the seasonal spatial and temporal patterns of CHL and SST with the empirical orthogonal function (EOF) statistical model in the YS, (ii) to investigate the inter-annual trends of CHL and SST in a long-term time series with the continuous wavelet transform (CWT) analysis, and (iii) to explore the temporal correlations between CHL and SST using wavelet coherency analysis at a regional scale.

DINEOF
Due to the cloud coverage of the MODIS images over the YS, MODIS pixel values were missing for some months.The EOF and wavelet analyses generally require a complete time series of input maps without data voids.Therefore, a method to reconstruct missing data based on the Data Interpolating Empirical Orthogonal Functions (DINEOF) decomposition was applied to obtain complete CHL and SST data (Beckers and Rixen, 2003;Beckers et al. 2006).It is a self-consistent, parameter-free technique for gappy data reconstruction.Recently, DINEOF has been used widely to reconstruct SST (Miles and He, 2010;Huynh et al., 2016), CHL and winds (Miles and He, 2010;Volpe et al., 2012;Liu and Wang, 2013;Liu et al., 2014), total suspended matter (Sirjacobs et al., 2011;Alvera-Azcarate et al., 2015), and sea surface salinity (Alvera-Azcarate et al., 2016).This technique presents some advantages over more classical approaches (such as optimal interpolation), especially when working on CHL and SST datasets (Miles and He, 2010;Volpe et al., 2012).(Campbell, 1995), raw data were log-transformed prior to reconstruction 135 to homogenize the variance and to yield a near-normal distribution (Fig. 2).These

Empirical Orthogonal Function (EOF) analysis
After DINEOF reconstruction, cloud free CHL values were log-transformed before we included them in figures and before statistical analysis.In Section 4, to better discern the spatial heterogeneity and the degree of coherence and temporal evolution of the CHL and SST fields, a traditional EOF analysis was applied further to the monthly, cloud-free DINEOF CHL and SST datasets, which is an approach that is also used widely in other disciplines (Hu and Si, 2016a).Each data set was organized in an M×N matrix, where M and N represented the spatial and temporal elements, respectively.Taking CHL for instance, the matrix (, ) can be represented by (, ) = ∑ ()()  =1 , where () are the temporal evolution functions and (, ) are the spatial eigen-functions for each EOF mode.Prior to EOF analysis, the temporal means of each pixel were removed from the original data using: , where ′(, ) are the resulting residuals (anomalies).The first two modes were decomposed to analyze the major variability in CHL and SST.
To assess the significance of the EOF modes, we followed the methods described by North et al. (1982).The error produced in a given EOF (e j ) was calculated as:   =   ( 2/) 0.5 , where λ is the eigenvalue of that EOF, and n is the degrees of freedom.When the difference between neighboring eigenvalues satisfied   −  +1 ≥   , then the EOF modes represented by these two eigenvalues were significant statistically.
The continuous wavelet transform (CWT) was used to determine the inter-annual scales of variability and the oscillation periods of DINEOF CHL and SST.
Prior to the CWT analysis, the seasonal variation of each pixel was removed from the original data.The CWT is a tool for decomposing the non-stationary time series at different spatial or time scales into the time-frequency space by translation of the mother wavelet and by analyzing localized variations of power (Messié and Chavez, 2011).The mother wavelets used in this study were the "Morlet" wavelets, which is used commonly in geophysics, because it provides a good balance between time and frequency localization (Grinsted et al., 2004;Hu and Si, 2016b;She et al., 2016).The CWT can localize the signal in both the time and frequency domains, but the classical Fourier transform was able to localize the signal only in the frequency domain with no localization in time (Olita et al., 2011).In addition, cross-correlation functions (Venables and Ripley, 2002) were used to determine the degree of temporal correspondence between the CHL and SST time series datasets, after we removed the seasonal variations.Then, the wavelet coherence was used to show the local correlation between CHL and SST in time-frequency space (Ng and Chan, 2012) that was based on the cross-correlation result.We used the wavelet software provided by

Monthly Climatology of CHL and SST
The CHL monthly means during 2003-2015 followed a similar pattern from month to month, with more CHL in the shallow coastal waters and a decreased in the seaward direction (Fig. 3).Although the maximum CHL appeared to be fairly consistent seasonally, the spatial extent of blooms had significant seasonal fluctuations.Monthly mean imagery showed the largest spatial coverage of CHL in YS was in spring and the smallest coverage was in summer.The CHL in coastal waters was relatively high in spring during every year.Some portions of phytoplankton blooms occurred only in subsurface waters, which made it impossible to see using satellite imagery.Overall, the CHL was the greatest in coastal waters in spring or in regions with greater diluted water, such as near the Yangtze River, where the CHL was characterized by a long-lasting summer CHL maximum that started in April and ended in September.
The seasonal cycle was more evident in the SST field (Fig. 4).SST showed a sinusoidal seasonal cycle, with a persistent, seasonal, warming trend from winter (February) to summer (August).The SST in the YS during December-April were below 12 °C , increased to 15 °C in May, reached a maximum above 25 °C in August, and then decreased again in September-October.From December to May, there was a drastic temperature difference between northern waters and southern waters, but the SST during summer months was nearly uniform over the entire YS.Spatially, isotherms were generally parallel to the isobaths.There was a clear temperature contrast between coastal waters and offshore waters in winter and spring.Similarly, the thermal front in southeast waters was more visible during winter and spring.The  temporal means of monthly data during 2003-2015 (Fig. 7a and c) and variability associated with the standard deviations (STD) of monthly mean temporal values (Fig. 7b and d).One general evident spatial pattern was that mean CHL showed a sharp decrease from coastal waters to offshore regions (Fig. 7a).In our study, the highest CHL values (~ 6 mg m -3 ) and the lowest STD (~ 0.02 mg m -3 ) were observed in coastal waters (Fig. 7b) that were adjacent to the mouth of the Yangtze River where the water depth was less than 20 m.Compared with the coastal waters and the sea adjacent to large river mouths, central YS waters had lower CHL, but they displayed greater variability.In these regions, strong water mixing could make more deep-ocean nutrients available for utilization by phytoplankton in some months.The spatially-averaged time series showed clear inter-annual variability that was superimposed on the seasonal spring (April) and fall (August) blooms (Fig. 7e).A noticeable scenario was the increasing trend in CHL of ~ 0.03 mg m -3 year -1 throughout the YS during 2003-2015; this phenomenon would require more observations of subsurface nutrients to understand the underlying mechanisms.
The spatial distribution of the 13-year averaged SST in the YS showed the mean SST with a smooth transition from colder water along the coast to warmer water in offshore areas (Fig. 7c).Similarly, SST in the northern YS was colder compared to that in the southern YS.In contrast, the STD calculated from the 13-year time series of SST also revealed a high spatial distinction between the northern and southern regions (Fig. 7d).In the northern region, the STD reached its highest values in excess of 8 °C , in contrast to the STD of the southeastern YSWC region, which were

The dominant CHL EOF mode
The first EOF (EOF1) and second EOF (EOF2) modes accounted for 40 % of the total CHL variability in this study (Fig. 8), which are similar to those observed in previous studies on the regional and global CHL (Messié and Radenac, 2006;Thomas et al., 2012).The CHL EOF1 mode explained 34 % of the total variance.The anomalies were not distributed uniformly throughout the entire study area (Fig. 8a).
One of the centers of CHL anomalies was located in an area with the geographical coordinates of about 35-37° N and 122-126° E, which was affected mainly by the YSCWM in the central waters of the YS (Teague and Jacobs, 2000;Lie et al., 2009;Yu et al., 2010).In that location, the stronger mixing of the water column brought the deeper nutrients upward, in turn favoring a phytoplankton bloom (Liu and Wang, 2013).Another CHL positive center was in the southeast waters close to the YSWC, which indicated that the EOF1 mode could be explained by the influences of the currents in the YS.As such, EOF1 mode is a good representation of differences in the timing of the blooms.The temporal amplitude showed positive values from winter to spring (November to April) but negative values from summer and autumn (June to October) (Fig. 8b).The result was related to the seasonal cycles, with a high CHL

The dominant SST EOF mode
The SST EOF1 mode accounted for 99 % of the total variance (Fig. 9a).SST anomalies in the EOF1 mode for the entire study area were all positive, but they were not distributed uniformly throughout the entire study area.This indicated that SST exhibited a positive trend, which was consistent with the pattern that was depicted in As a result of water depth, the SST spatial EOF1 was highly correlated with the distribution of water depth in the YS.The magnitude of variability of EOF1 in shallow water was larger compared to that in the slope region, which suggested an inverse relationship between the general pattern of CHL and bathymetry (O'Reilly et al., 1987).The strong match between the mean CHL and SST patterns throughout the entire YS region can be explained in terms of lower primary production levels that corresponded to stronger stratification of the water column (Behrenfeld et al., 2006;Doney, 2006) and, thus, to warmer surface waters (Wilson and Coles, 2005).This is because the variation in SST, in the first order, one-dimensional sense, is inversely proportional to water depth (He and Weisberg, 2003).The shallow ocean waters overall have a larger seasonal cycle.In contrast, SST in the slope region remained  is linearly proportional to the water column in a shallow ocean (Yan et al., 1990; 328 Chen et al., 1994;Xie et al., 2002;Ichikawa and Beardsley, 2002;Xie et al., 2002; 329 Park et al., 2005).The temporal amplitude in SST was dominated strongly by a 330 seasonal periodicity that peaked in summer and winter (Fig. 9b).The influence of 331 EOF2 mode on SST variability could be omitted, since it only accounted for 1 % of 332 the total variance (Fig. 9c).EOF2 mode.

Scales of variability and oscillation periods on an inter-annual timescale
The CWT were applied to the long-term monthly CHL and SST datasets after removing seasonal variations.The wavelet power spectrum and the global power spectrum obtained through the Morlet wavelet transform highlighted the dominant scales of variability and oscillation periods of CHL and SST (Fig. 10).The global power spectra showed the multi-period for CHL and SST (right panels of Fig. 10a ).In our study, the statistically significant cross-correlation between the monthly 362 CHL and SST datasets after removing seasonal variations (Fig. 11a) also suggested

Fig. 1 .
Fig. 1.(a) Bathymetric and geographic map of the study area in the Yellow Sea,101 -Aqua data for CHL and SST during January

107 2003 -
December 2015 were used in this study.The CHL and SST data were level 3 108 fields provided by the NASA ocean color web page (http://oceancolor.gsfc.nasa.gov).109 The standard CHL product that was derived from the OC3Mv5 algorithm (OC3M 110 Ocean Sci.Discuss., doi:10.5194/os-2017-11,2017 Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.updated version after the 2009 reprocessing) and the daytime SST 11 μm product (which uses the 11 and 12 μm bands) were obtained.The level 3 product was collected in a 4 km spatial resolution from 30-40° N in latitude and 118-126° E in longitude for the YS region.
CHL and SST are characterized by different scales of variability in coastal or open ocean areas.This method identifies dominant spatial and temporal patterns in CHL and SST datasets, and it fills in missing data.Thus, DINEOF was applied to Ocean Sci.Discuss., doi:10.5194/os-2017-11,2017 Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.reconstruct the missing CHL and SST data in this study.Because the satellite CHL 133 values spanned three orders of magnitude and CHL retrievals are often distributed 134 log-normally 136 images show clearly the utility of the DINEOF method in reconstructing monthly, 137 high-resolution imagery from datasets with large amounts of cloud cover.For 138 example, in the CHL, DINEOF gives a low concentration of CHL in the southeast 139 regions of the YS (Fig. 2b).

2373. 3
Monthly mean and temporal variability of CHL and SST 238 The spatial patterns of CHL and SST concentration were produced by the 239 Ocean Sci.Discuss., doi:10.5194/os-2017-11,2017 Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.
Fig. 7. Long-term temporal mean and standard deviation maps of DINEOF CHL267 Fig. 8.The first two prominent EOF modes for CHL variability using DINEOF CHL

Fig. 7f .
Fig. 7f.Similar trends were observed by Liu and Wang (2013) and Park (2015) in OceanSci.Discuss., doi:10.5194/os-2017-11,2017   Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.fairly consistent especially in winter due to an increase in water depth (heat content) 326 and the persistent warm water supply from the southeast water because thermal inertia 327

363
that variability in CHL was slightly negatively correlated with variability in SST in 364 the YS during the study period.The negative correlation was confirmed also by the 365 scatter plot of CHL and SST (Fig. 11b).The cross-correlation between monthly CHL 366 Ocean Sci.Discuss., doi:10.5194/os-2017-11,2017 Manuscript under review for journal Ocean Sci. Discussion started: 28 April 2017 c Author(s) 2017.CC-BY 3.0 License.and SST showed a significant and negative cross-correlation (R = -0.21,p < 0.01) 367 with time lags of 4 months, which suggested that the CHL reached the maximum 368 value 4 months after the SST got the minimum value in the YS.Therefore, to further 369 examine the synchrony between CHL and SST, wavelet coherence analysis was used 370 to reveal the coherency between CHL and SST (Fig. 12).The wavelet squared 371 coherence values below the confidence level indicated that there were some randomly 372 distributed sections.The vectors indicated the phase difference between CHL and SST 373 at each time and period.The portion of figure 12 with significant correlation showed 374 the anti-phase relationship between CHL and SST, with arrows pointing left, which 375 suggested that the two series apparently have negative coherency in the 1.5-3 year 376 band during 2008-2015. 377