The mixed layer heat budget in the tropical Pacific is diagnosed using pentad (5 day) averaged outputs from the Global Ocean Data Assimilation System (GODAS), which is operational at the National Centers for Environmental Prediction (NCEP). The GODAS is currently used by the NCEP Climate Prediction Center (CPC) to monitor and to understand El Niño and La Niña in near real time. The purpose of this study is to assess the feasibility of using an operational ocean data assimilation system to understand SST variability.
The climatological mean and seasonal cycle of mixed layer heat budgets derived from GODAS agree reasonably well with previous observational and model-based estimates. However, significant differences and biases were noticed. Large biases were found in GODAS zonal and meridional currents, which contributed to biases in the annual cycle of zonal and meridional advective heat fluxes. The warming due to tropical instability waves in boreal fall is severely underestimated owing to use of a 4-week data assimilation window. On interannual time scales, the GODAS heat budget closure is good for weak-to-moderate El Niños. A composite for weak-to-moderate El Niños suggests that zonal and meridional temperature advection and vertical entrainment/diffusion all contributed to the onset of the event and that zonal advection played the dominant role during decay of the event and the transition to La Niña. The net surface heat flux acts as a damping during the development stage, but plays a critical role in the decay of El Niño and the transition to the following La Niña.
The GODAS heat budget closure is generally poor for strong La Niñas. Despite the biases, the GODAS heat budget analysis tool is useful in monitoring and understanding the physical processes controlling SST variability associated with ENSO. Therefore, it has been implemented operationally at CPC in support of NOAA’s ENSO forecasting.
Understanding changes in sea surface temperature is key to understanding the coupled atmosphere–ocean system. For example, for better understanding and ability to forecast the El Niño–Southern Oscillation (ENSO), which is the dominant mode of coupled ocean–atmosphere variability in the tropical Pacific, many studies have analyzed the physical mechanisms that govern the seasonal cycle and interannual variability of SST (Stevenson and Niiler 1983; Hayes et al. 1991; Chen et al. 1994; Kessler et al. 1998, hereafter KRC98; Wang and McPhaden 1999 (hereafter WM99), 2000, 2001a; Swenson and Hansen 1999; Vialard et al. 2001; Kim et al. 2007).
The near-surface ocean is forced by winds, downward shortwave and longwave radiation fluxes, and freshwater fluxes. The ocean then impacts the atmosphere via latent, sensible, and longwave radiative heat losses that are dependent on SST and near-surface atmospheric variables. Since SST is closely related to mixed layer temperature variability, SST variations are intimately connected with the heat budget of the mixed layer. Various approaches, differing in their use of input data, have been taken to analyze the heat budget of the mixed layer. One approach is the use of observational data. Because of the scarcity of the observational data, however, such analyses have difficulty in accurately calculating the necessary horizontal and vertical gradient terms in the heat budget equations (Hayes et al. 1991; WM99).
An alternate approach is the use of output from model simulations (Chen et al. 1994; KRC98; Vialard et al. 2001; Zhang et al. 2007; Zhang 2008, hereafter ZH08). Although the analysis based on model simulations can precisely calculate various terms in the budget equations, such analyses can deviate substantially from observed reality because of the uncertainty in atmospheric forcing and other model biases. Further, because of nonlinearities, heat budgets may close when data from daily model outputs are analyzed, but may not when only monthly outputs from the model simulations are available (Zhang et al. 2007). A third approach is the use of output from an ocean data assimilation system (Kim et al. 2007). A particular advantage of using ocean assimilation products is that the model solutions are partially constrained by observations so that departures from the observations, unlike for the model simulations, may not be as large.
Kim et al. (2007) used the data assimilation product called Estimating the Circulation and Climate of the Ocean (ECCO; available online at http://www.ecco-group.org) to analyze the mixed layer temperature variability in the Niño-3 region. ECCO is an adjoint-based estimation system that demands the estimated state satisfy the model equations exactly over a certain time interval while adjusting control variables, which are typically the initial state, surface forcing, and model parameters, so that the estimated states are as close to observations as possible. Kim et al. suggested that such systems ensure consistency of the estimated surface forcing with the estimated ocean state, thus guaranteeing the closure of heat budgets.
In this study, we use the pentad (5 day) averaged outputs from the Global Ocean Data Assimilation System (GODAS) (Behringer et al. 1998; Behringer and Xue 2004) produced at the National Centers for Environmental Prediction (NCEP). GODAS is a sequential estimation system that allows the estimated state to deviate from an exact solution of the underlying physical model by applying statistical corrections to the state. These corrections often make estimated states close to observations, but they imply internal sources and sinks of heat, salt, and momentum, et cetera. Therefore, the heat budgets derived from GODAS will not have a perfect closure as that in Kim et al. (2007). However, we will show that the heat budget derived from GODAS is approximately closed on seasonal to interannual time scales. In particular, this budget is useful in understanding and monitoring the physical processes controlling the SST variability associated with ENSO.
Previous model and observational studies have suggested that the mechanisms for mixed layer temperature variability are very complicated. As an example, for the seasonal cycle the net surface heat flux, subsurface entrainment/diffusion cooling, and tropical instability waves (TIWs) all play an important role (KRC98; WM99; Philander et al. 1986; Contreras 2002; Jochum and Murtugudde 2006). For the eastern Pacific on interannual time scales, vertical entrainment/diffusion is the most critical process controlling interannual SST variability (Harrison et al. 1990; Frankignoul et al. 1996; Wang and McPhaden 2000, 2001a; Zhang et al. 2007; Kim et al. 2007), while the surface heat fluxes act to damp interannual SST variations. For the central and western equatorial Pacific, studies have suggested that zonal advection by anomalous currents is the dominant mechanism for SST variation on interannual time scales (Kessler and McPhaden 1995).
For a heat budget analysis based on the output of an ocean data assimilation system, questions remain about how well earlier conclusions can be replicated and what new ones can be learned. In this study, we use the pentad (5 day averaged) outputs from GODAS to diagnose heat budgets of the mixed layer in the tropical Pacific. GODAS outputs have been extensively used at the Climate Prediction Center (CPC) of NCEP to monitor global ocean variability and its interaction with the atmosphere (see CPC’s monthly ocean briefing archive online at http://www.cpc.ncep.noaa.gov/products/GODAS). An advantage of using GODAS outputs for the mixed layer heat budget analysis is that, if realistic, it can be routinely updated in real time to monitor the mixed layer heat budget and to understand the sources of SST variability (particularly on ENSO time scales) in the tropical Pacific.
The purpose—and a unique aspect—of the paper is to demonstrate the feasibility of an ocean data assimilation product, that is, GODAS, for the analysis of the evolution of the mixed layer in the tropical Pacific. We will discuss the realistic and potentially problematic features of the analysis for the annual mean and seasonal cycle, as well as for interannual variability of the mixed layer temperature in the tropical Pacific. Based on the results and comparison with earlier studies, we demonstrate that the analysis of the mixed layer heat budget from an operational ocean assimilation system is an effective tool to monitor and understand SST variability on ENSO time scales. Special attention will be given to the issue of heat budget closure when the dynamical consistency of model solutions is not maintained due to the ingestion of data in the assimilation cycle.
We briefly describe the NCEP operational GODAS in section 2 and data and validation procedures in section 3. The methodology for the mixed layer heat budget calculations is discussed in section 4. Mixed layer heat budget governing the mean, seasonal cycle, and composite El Niño is presented in section 5.
2. A description of the NCEP GODAS
GODAS was implemented at NECP in 2004 (Behringer and Xue 2004) and is currently used to initialize the oceanic component of the NCEP Climate Forecast System (Saha et al. 2006). It replaced the Pacific Ocean Data Assimilation System (ODAS) version RA6 (Ji et al. 1995; Behringer et al. 1998). The major changes from the RA6 included 1) an extension to a quasi-global domain (75°S–65°N); 2) a replacement of the Geophysical Fluid Dynamics Laboratory’s Modular Ocean Model version 1 with version 3 (MOM3) (Pacanowski and Griffies 1999); 3) a change from momentum flux forcing only to momentum, heat, and freshwater flux forcings from the NCEP/Department of Energy Global Reanalysis 2 (R2 hereafter) (Kanamitsu et al. 2002); and 4) a change in the assimilation from temperature only to temperature and synthetic salinity that is constructed from temperature and a local temperature/salinity climatology.
The ocean model has a resolution of 1° × 1° that increases to ⅓° in the north–south direction within 10° of the equator and has 40 levels with a 10-m resolution in the upper 200 m. Other features of MOM3 include an explicit free surface, the Gent–McWilliams isoneutral mixing scheme (Gent and McWilliams 1990), and the K-profile parameterization (KPP) vertical mixing scheme (Large et al. 1994).
Temperature observations assimilated into GODAS include data from expendable bathythermographs (XBTs), Tropical Atmosphere Ocean (TAO) array in the tropical Pacific, Triangle Trans-Ocean Buoy Network (TRITON) in the tropical Indian Ocean, Prediction and Research Moored Array in the Tropical Atlantic (PIRATA), and Argo profiling floats (see references cited in Huang et al. 2008). In the assimilation cycle, the model state is corrected by observations within a 4-week window centered on the model time using a three-dimensional variational data assimilation (3DVAR) scheme (Behringer et al. 1998). The 4-week assimilation window is effective in eliminating unrealistic small-scale variations and improving large-scale structures, but it severely smoothes out variations associated with tropical instability waves (TIWs). Owing to the lack of direct salinity observations, synthetic salinity profiles constructed from temperature and a local T–S climatology are also assimilated into GODAS. During the assimilation cycle the surface fluxes from R2 are further corrected by restoring the model temperature of the first layer (5 m) to the optimal interpolation (OI) SST analysis version 2 (Reynolds et al. 2002) and restoring the model surface salinity to the annual sea surface salinity (SSS) climatology (Conkright et al. 1999). The restoring time scale is 5 days for temperature and 10 days for salinity. The strong restoration to observed SST is necessary so that the model SST is close to observations. The heat flux correction due to the SST relaxation is significant and has been included in our heat budget analysis.
GODAS has only pentad and monthly outputs. This study uses the pentad outputs of temperature, salinity, and three-dimensional ocean currents on a common 1° × 1° grid in the 1979–2008 period. The choice of pentad fields and 1° × 1° grid has little negative impact on the GODAS heat budget analysis since TIWs are severely underestimated in GODAS due to its use of a 4-week data assimilation window.
3. Data and GODAS validation
Observed data and analyses are used to validate GODAS, which include monthly and weekly OI SST, climatological salinity and temperature from the World Ocean Database 2001 (WOD01) (Conkright et al. 2002); pentad currents from Ocean Surface Current Analyses–Real Time (OSCAR) (Bonjean and Lagerloef 2002); daily temperature, salinity, and currents from TAO moorings (available online at http://www.pmel.noaa.gov/tao); surface heat fluxes from objectively analyzed air–sea fluxes (OAFlux) (Yu et al. 2008); and solar and longwave radiation heat fluxes from the International Satellite Cloud Climatology Project (ISCCP) (available online at http://isccp.giss.nasa.gov). Earlier validation of GODAS suggested that the temperature field is closer to observations than the Pacific ODAS and that the poor salinity field in ODAS is dramatically improved (Behringer and Xue 2004). Although this version of GODAS does not assimilate satellite altimetry, the sea surface height in GODAS is also reasonably consistent with altimetry and tide gauge records (Behringer and Xue 2004; Behringer 2007).
Before analyzing the mixed layer heat budget, we first quantify the accuracy of GODAS in representing the mixed layer temperature and ocean currents. The seasonal cycle of temperature along the equator is well simulated by GODAS, and differences from the observed seasonal cycle are generally less than 0.5°C (not shown). Zonal current along the equator (1°S–1°N) is compared with OSCAR (Fig. 1). The OSCAR currents (Fig. 1a) are based on an analysis of satellite altimeter and scatterometer measurements, and the seasonal cycle is based on the 1993–2007 analysis period (available online at http://www.oscar.noaa.gov/index.html). Compared to OSCAR currents, GODAS has a westward bias in the far western and eastern equatorial Pacific and an eastward bias in the central Pacific between 180° and 120°W (Fig. 1c). Biases in the western and central Pacific are likely associated with the assimilation of synthetic salinity, as these biases are dramatically reduced in an experimental GODAS assimilation run in which observed salinity from Argo floats is also assimilated (Behringer 2007).
Next, GODAS and OSCAR currents are compared with the measurements at four TAO mooring locations along the equator in the western (165°E), central (170°W), and eastern (140° and 110°W) Pacific. Figures 2a–h show the annual cycles of zonal and meridional currents from GODAS, OSCAR, and TAO observations, and Tables 1, 2 show the comparison statistics between OSCAR and TAO and between GODAS and TAO. In Tables 1, 2, the mean bias was calculated as the mean difference between model and TAO data for the common period of two datasets; rms errors (RMSEs) were calculated with the total currents; anomaly correlation coefficients (ACCs) and anomaly RMSEs (ARMSEs) were calculated from currents for which the means have been removed.
For zonal currents, OSCAR generally agrees with TAO better than GODAS does. The mean biases are respectively −18, 13, −4, −2 cm s−1 at 165°E, 170°W, 140°W, and 110°W in OSCAR; and −24, 23, 13, and −18 cm s−1 in GODAS (Table 1). The RMSEs are 19, 14, 8, and 5 cm s−1 at 165°E, 170°W, 140°W, and 110°W in OSCAR; and 26, 26, 15, and 18 cm s−1 in GODAS (Table 1). Interestingly, both OSCAR and GODAS have reasonably high ACCs (0.93, 0.94, 0.95, and 0.98 in OSCAR; 0.76, 0.80, 0.96, and 0.98 in GODAS) with TAO observations. The anomalous RMSEs (5, 6, 7, and 5 cm s−1 in OSCAR; 10, 11, 7, and 5 cm s−1 in GODAS) are much smaller than RMSEs that include the mean biases. In summary, both GODAS and OSCAR have large mean biases in zonal currents in the western (165°E) and central (170°W) Pacific, and GODAS has much larger mean biases than OSCAR in the eastern (140° and 110°W) Pacific. Once the mean biases are removed, both OSCAR and GODAS simulate TAO observations reasonably well. It will be shown in section 5 that GODAS is quite adequate in simulating anomalous zonal advective heat flux in the central–eastern tropical Pacific associated with ENSO.
For meridional currents, the OSCAR estimates are generally too weak and bear little resemblance to TAO observations (Figs. 2e–h). In contrast, GODAS currents have amplitudes comparable to those of observations. GODAS meridional currents are superior to OSCAR meridional currents in the western (165°E) and the central–eastern (170° and 140°W) Pacific. The RMSEs are respectively 4, 8, 4, and 2 cm s−1 at 165°E, 170°W, 140°W, and 110°W in OSCAR,; and 2, 4, 4, 3 cm s−1 in GODAS (Table 2). The ACCs are 0.59, 0.90, 0.56, and 0.44 in GODAS but near zero in OSCAR, except in the far eastern (110°W) Pacific.
It is very interesting that, for zonal currents, OSCAR is generally superior to GODAS in the central and eastern Pacific (170°, 140°, and 110°W) but both are poor in the western Pacific (165°E). An experimental GODAS run suggested that the large biases in GODAS zonal currents can be significantly reduced when the Argo salinity is assimilated (Behringer 2007). For the meridional currents, GODAS is generally superior to OSCAR in the western (165°E) and central–eastern (170° and 140°W) Pacific, but OSCAR is superior to GODAS in the eastern Pacific (110°W). However, the amplitude of GODAS meridional currents is more realistic than for OSCAR, which is too weak. The smaller RMSE in OSCAR than in GODAS in the eastern Pacific is probably due to the smaller amplitude in OSCAR.
4. Methodology for analyzing the mixed layer heat budget
a. Mixed layer depth
The criterion to calculate mixed layer depth (MLD) is often defined differently based on requirements of the analysis (You 1995; Sprintall and Tomczak 1992). We select the criterion to be a density difference of 0.125 kg m−3 between the surface and the bottom of the mixed layer. The results of the heat budget analysis, however, are not sensitive to the choice of the criterion. In fact, similar results were obtained when the criterion was chosen to be a temperature difference of 0.5 K.
The MLD of GODAS is calculated using the pentad fields of temperature and salinity. The seasonal cycle of the GODAS MLD is calculated based on the 1982–2004 period and is compared with the MLD of WOD01, calculated using monthly climatological fields of temperature and salinity in WOD01.
The seasonal cycle of MLD along the equator (1°S–1°N) from WOD01 and GODAS, and their differences, is shown in Fig. 3. The WOD01 MLD is relatively shallow (deep) in the western and eastern (central) tropical Pacific (Fig. 3a). The shallow MLD in the eastern tropical Pacific is associated with the shallow thermocline maintained by easterly trade winds in the central tropical Pacific. The shallow MLD in the western tropical Pacific is associated with excess precipitation over evaporation that forms a barrier layer (Sprintall and Tomczak 1992; Ando and McPhaden 1997). Compared to the WOD01, the MLD based on GODAS is about 20–30 m too deep in the west-central Pacific through the calendar year and ∼10–20 m too deep in the eastern Pacific during boreal fall.
b. Mixed layer temperature equation
where Tt = ∂Ta/∂t is the mixed layer temperature tendency and F is the combined forcing of zonal advection (Qu), meridional advection (Qυ), vertical entrainment (Qw), adjusted surface heat flux (Qq = Qadj/ρcph), and vertical diffusion (Qzz = Qdiff/ρcph); Qadj is the net surface heat flux plus heat flux correction minus the penetrative shortwave radiation [see Eq. (A4)]. Weak horizontal diffusion was ignored in our analysis.
To understand the physical processes that control the temperature variations in the mixed layer on different time scales, each variable associated with forcing F in Eq. (2) is decomposed into low frequency variation (≥75 day) and high frequency transients (hereafter referred to as eddy). Therefore, Eq. (1) becomes
where superscript L indicates the term calculated using low-pass filtered variables and E represents the combined terms from high frequency eddies (see details in appendix B). Equation (3) is further decomposed into climatology (bar) and its anomaly (prime). The equation for anomalous temperature is, omitting superscript L,
Details about each term in Eq. (4) are described in appendix B. The climatological mean and annual cycle of each term in Eq. (3) will be discussed in sections 5b and 5c. The anomalous heat budgets described by Eq. (4) will be used to construct a composite El Niño heat budget, and the characteristics of the anomalous heat budgets for a typical El Niño will be discussed in section 5d.
A cutoff period of 75 days is chosen to separate seasonal and longer time scale variability from that associated with TIWs, which exhibit a typical period of 20–30 days (Jochum and Murtugudde 2006). The annual climatology of the heat budgets did not change much when different cutoff periods between 30 and 90 days were selected (as also indicated by KRC1998). This suggests that 60–90 day period oceanic Kelvin waves, forced by the atmospheric Madden–Julian oscillation (MJO) and westerly wind bursts, do not make significant contributions to the climatological heat budgets. However, the cumulative effects of a sequence of oceanic Kelvin waves make a significant contribution to the anomalous heat budget on seasonal time scales and are believed to influence the onset and determination of El Niño (Seo and Xue 2005). In this paper, we will describe the heat budget using an El Niño composite, which tends to smear out the contributions of oceanic Kelvin waves that may be of importance in the analysis of specific El Niño events. The topic of how a sequence of oceanic Kelvin waves contributes to the anomalous heat budget on seasonal time scales during a specific El Niño event will be explored in a separate paper.
c. Closure of the temperature equation
To test the procedures used for computing mixed layer budgets, we first apply the proposed methodology to a control simulation (hereafter referred to as CNTRL) that is identical to GODAS except no observations are assimilated. The pentad fields from CNTRL are used in the calculation of all heat budget terms in Eq. (1). The pentad climatology is calculated for 1982–2004, and pentad anomalies are obtained by removing the pentad climatology from each heat budget term. The closure of the heat budgets is measured by the consistency between Tt and F of Eq. (1). Figures 4a,b show the time evolution of the Tt and F for the annual cycle and interannual variability in the Niño-3.4 region (5°S–5°N, 120°–170°W) during 1979–2007 in the CNTRL run. For both the annual cycle and the interannual variability, there is a close resemblance in the tendency and forcing term. Temporal correlations between Tt and F are above 0.95, and the RMSEs are less than 0.09°C month−1. However, in the annual cycle F is about 0.1°C month−1 cooler than Tt from June to December. The cold bias may be related to the underestimation of the eddy warming in CNTRL during summer/fall. In fact, the annual mean eddy warming averaged in the region 0°–4°N, 90°–140°W in CNTRL is 0.5°C month−1, significantly weaker than 0.8°C month−1 derived in the model study by Richards et al. (2009). The underestimation of the eddy warming might be due to the use of 5-day averaged fields and the 1° × 1° grid. The results nonetheless suggest that the temperature equation and its closure are approximately satisfied.
Figures 4c,d show the seasonal cycle and interannual variability of Tt and F for GODAS. Data assimilation is expected to introduce sources and sinks of heat and generate inconsistency between the forcing fields and the analyzed ocean states that will negatively impact the closure of the mixed layer heat budget. However, the mean seasonal cycle of Tt and F follow each other closely, suggesting that the heat sources and sinks due to data assimilation have only minor impact on the climatological heat budget. For the seasonal cycle, the ACC and RMSE between Tt and F are 0.97 and 0.06, very similar to those for CNTRL. The influence of data assimilation on the heat budget is more evident in the evolution of Tt and F anomalies, with an ACC (RMSE) of 0.70 (0.23), which is smaller (larger) than those for CNTRL. Note that a few factors contribute to the imbalance between Tt and F, including sources and sinks of heat due to data assimilation, and uncertainties in the parameterization of vertical entrainment and vertical diffusion and the use of pentad fields and the 1° × 1° grid. Therefore, we should not be surprised when Tt, and F are out of balance, as evident during the strong El Niño events (1982/83 and 1997/98; Fig. 4d). We will make a composite heat budget for those weak-to-moderate El Niño events where the closure is reasonably good and then analyze the heat budget during the strong El Niño events separately where the closure is poor. We will also diagnose errors in the heat budget through comparison with other observational and model heat budget analyses.
5. Analysis of the mixed layer heat budget
a. Surface heat fluxes
The net surface heat flux plays a critical role in forcing temperature changes in the mixed layer in the equatorial Pacific. The net surface heat flux forcings used in GODAS is specified from the NCEP R2 reanalysis, and the climatology of various components averaged in 1982–2004 is shown in Figs. 5a–d. Also shown are the climatologies of penetrative shortwave radiation (Fig. 5e), the surface heat flux correction (Fig. 5f) due to the SST relaxation in GODAS, and the adjusted net surface heat flux [Qadj, Fig. 5g; see Eq. (A4)], which is the net surface heat flux plus flux correction minus the penetrative radiation.
Shortwave radiation exhibits a clear semiannual cycle along the equatorial Pacific as the sun crosses the equator twice a year (Fig. 5a). It has two maxima near the date line with amplitude 240 W m−2 in October and 200 W m−2 in April. The latent heat flux also has a semiannual cycle with two maxima of 140 W m−2 in boreal winter and fall. A minimum of ∼100 W m−2 occurs in spring in the west-central tropical Pacific (Fig. 5d). The longwave radiation represents a net heat loss from the ocean to the atmosphere with an average amplitude of 50 W m−2 (Fig. 5b). The sensible heat flux is generally less than 5 W m−2 (Fig. 5c). Seasonality in longwave radiation and sensible heat flux is also small.
Penetrative shortwave radiation has a semiannual cycle with two maxima of about 40 W m−2 in boreal spring and fall in the far eastern tropical Pacific where the MLD is shallow (Fig. 5e). The heat flux correction due to SST relaxation is positive with a maximum of about 30 W m−2 in the eastern tropical Pacific in early spring (Fig. 5f) when the model cold bias is the largest (not shown). Because of the semiannual cycle in the shortwave radiation and latent heat fluxes, the adjusted net heat flux also has a clear semiannual cycle with two maxima in boreal spring and late fall in the central and eastern tropical Pacific (Fig. 5g). The analysis indicates that the longwave radiation, penetrative shortwave radiation, and heat flux corrections all make significant contributions to the closure of the heat budget, although their magnitudes are relatively small compared with shortwave radiation and latent heat fluxes.
The latent and sensible heat fluxes shown in Fig. 5 are similar to the OAFlux, and differences are generally less than 10 W m−2 (not shown). Compared with the net surface heat flux derived from the combination of ISCCP (shortwave and longwave) and OAFlux (latent and sensible) products, the net surface heat flux from R2 is 40–60 W m−2 too low in the equatorial Pacific (Fig. 5h), largely due to deficiencies in shortwave radiation. Shortwave radiation is 40–60 W m−2 too low in boreal spring compared to ISCCP and WM99. Satellite data from ISCCP suggests that mean shortwave radiation is as large as 280 W m−2 in boreal fall and 260 W m−2 in early boreal spring in the central equatorial Pacific (information available online at http://oaflux.whoi.edu).
To constrain the drift in surface temperature, GODAS includes a surface heat flux correction by relaxing model SST to observed SST. The mean flux correction is about 10–30 W m−2 (Fig. 5f). This correction partially compensates for biases in the net surface heat flux that, if not corrected, would lead to a cooling of upper ocean temperature during the assimilation cycle. However, the correction is not enough to compensate for all the deficiencies in the R2 net surface heat flux, which explains why GODAS surface temperature is still about 0.2°–0.4°C cooler than the observed SST (not shown).
b. Annual-mean mixed layer heat budget
Shown in Fig. 6 is the annual mean of the mixed layer heat budget calculated with the low-pass filtered GODAS data [Eq. (3)]. The mixed layer is heated on average by the adjusted surface heat flux (Qq) at a rate of 0.2–0.5°C month−1 in the central and western Pacific and 1–3°C month−1 in the eastern tropical Pacific (Fig. 6d). The heating is largely balanced by the cooling from meridional advection (Qυ, Fig. 6b), vertical entrainment (Qw, Fig. 6c), and vertical diffusion(Qzz, Fig. 6e). The maximum cooling by meridional advection is centered off the equator with a magnitude of 1°C month−1 near 2°N and 0.5°C month−1 near 3°S in the eastern tropical Pacific. The cooling from Qw is 0.2°C month−1 in the central tropical Pacific and 1°C month−1 in the far eastern tropical Pacific. The cooling is larger from vertical diffusion than from vertical entrainment, and the meridional extension is also broader because upwelling is mainly constrained within a narrow equatorial band. The annual mean zonal advection (Qu, Fig. 6a) contributes to a weak cooling (0.2°C month−1) across much of the tropical Pacific. TIW heating (Fig. 6f) is approximately 0.2°C month−1 in the eastern tropical Pacific east of 150°W, and is weaker than in CNTRL (0.5°C month−1), ZH08 (0.5°C month−1), Richards et al. (2009; 0.8°C month−1), and Jochum and Murtugudde (2006; 2°C month−1).
c. Seasonal cycle of the mixed layer heat budget
1) Seasonal cycle in GODAS
The seasonal cycle of the mixed layer heat budget in the equatorial Pacific (0.5°N, which is selected for the purpose of comparison in the following subsection) is discussed next. The mixed layer temperature tendency has a strong seasonal cycle in the equatorial eastern Pacific (Fig. 7e). The positive tendency in early boreal spring is largely due to excess heating by the adjusted net surface heat flux (Qq) (Fig. 7d) over cooling by vertical entrainment and diffusion (Qw + Qzz) (Fig. 7c). In contrast, the negative tendency during late boreal spring to summer is due to cooling by Qw + Qzz and Qυ dominating over heating by Qq.
The heating by Qq is dominated by a semiannual cycle (Fig. 7d). It is the largest in the equatorial eastern Pacific, with a primary maximum in boreal spring (4°C month−1) and a secondary maximum in boreal fall (2°C month−1). The seminannual variation of the heating by Qq is critically controlled by cloudiness and mixed layer depth, as was pointed out by KRC98. The magnitude of Qq is much weaker (0.5°C month−1) in the central and western tropical Pacific than in the east, because the MLD is relatively deep in the west. Heating in the west–central tropical Pacific is dominated by the semiannual signal in shortwave radiation (Fig. 5a), which largely governs the temperature tendency (Fig. 7e).
Cooling by Qw + Qzz remains confined to the central and eastern tropical Pacific, where the MLD is shallow and the vertical temperature gradient is the largest. The cooling has a primary maximum in boreal spring and a secondary maximum in boreal summer.
The contribution of Qu is the largest east of 130°W (Fig. 7a) and is dominated by the seasonal cycle: cooling from February to June and a heating from July to January. This cooling exhibits westward propagation, with cooling in the central and eastern Pacific in the fall (Fig. 7a). The cooling from Qυ is the largest from May to December when northerly currents are the strongest (Fig. 2h). The seasonal cycle of TIW heating (Fig. 7f) is about 0.5°C month−1 east of 130°W from June to December.
2) Comparison with other model simulations
The seasonal variation of Qq in the equatorial Pacific is similar to that in WM99. The seasonal variation of Qw + Qzz is close to that of ZH08. The pattern of TIW heating is almost the same as in ZH08.
A potential bias in GODAS is the zonal advective cooling in the eastern tropical Pacific in boreal spring (Fig. 7a), which is inconsistent with results reported in earlier studies. KRC98 and ZH08 suggested that Qu is weakly positive east of 120°W during boreal spring largely owing to the spring reversal of the South Equatorial Current (SEC). The erroneous cooling by Qu in GODAS is associated with errors in the surface zonal currents in GODAS for which the spring reversal of SEC does not extend as far eastward as in OSCAR (Fig. 1). The negative Qu in GODAS is generated by the westward surface zonal currents east of 105°W (Fig. 1b). Compared with the analysis of ZH08, the cooling in boreal fall and winter is confined too narrowly in the eastern Pacific and is related to the eastward biases in GODAS surface zonal currents in the region (Figs. 1c, 2a–d).
Cooling from Qυ in the eastern Pacific appears too strong in GODAS (0.5°–3°C, Fig. 7b) when compared with that in ZH08 (0.5°C). Considering that GODAS mean meridional currents do not agree well with observations in the eastern Pacific at 110°W (Table 2) and GODAS mean temperature has cold biases in the region (not shown), the Qυ climatology is likely problematic in the eastern Pacific. In addition, the TIW heating in GODAS (Fig. 7f) is only half of ZH08 (1°C month−1) and CNTRL (1°C month−1) and extends westward to only 125°W compared to 150°W in ZH08 and CNTRL.
3) Comparison with observational analyses
We use the observational analysis of WM99 to further validate the heat budgets of the mixed layer in GODAS. WM99 used observed winds, temperature, and ocean currents from the TAO moorings at four locations along the equator in the western (165°E), central (170°W), and eastern (140° and 110°W) Pacific. Changes in heat storage, horizontal heat advection, and heat fluxes at the surface in WM99 were estimated directly from data. In their estimates, vertical heat flux out of the base of the mixed layer was calculated as residual and surface heat fluxes were from COADS. The heat budgets in GODAS are calculated using unfiltered (total) data to facilitate comparison with WM99.
The temperature tendencies in the mixed layer of GODAS at the four TAO sites are shown in Fig. 8a. According to WM99, the temperature at 110°W warms from September to March and then cools from April to August. This seasonal variation is replicated by GODAS except the strong cooling tendency in June is underestimated and the weak warming tendency in November is missing. Because of the westward propagation of positive climatological SSTs, the peak warming tendency at 110°W, 140°W, 170°W, and 165°W subsequently progresses westward from February to April. This westward propagation of climatological SST is well simulated by GODAS. The secondary maximum warming at 165°E in boreal fall indicates a semiannual cycle in the western tropical Pacific (Yuan 2005).
Zonal advection, Qu, generally cools the eastern (110°W) tropical Pacific during August to February when the SEC is westward and warms it during the spring reversal of the SEC (Fig. 8b). The warming at 110°W in March–June (Fig. 6b in WM99) is simulated as cooling in GODAS (Fig. 8b), mainly because GODAS underestimates the spring reversal of the SEC in the far eastern Pacific. The cooling at 110°W in boreal fall is seriously overestimated by GODAS owing to westward biases in the surface zonal currents (Figs. 1c and 2d). Relative to 110°W, Qu at 140° and 170°W is much better simulated by GODAS, in good agreement with the observational analysis of WM99. Zonal advection at 165°E cools throughout the year and is largely consistent with WM99.
Meridional advection Qυ, at 110° and 140°W leads to warming along the equator with maximum amplitude in boreal summer and fall due to active TIWs (WM99; Chen et al. 1994; KRC98). The Qυ reaches a minimum during early boreal spring when TIWs are inactive (WM99). In GODAS, warming by Qυ at 110° and 140°W (Fig. 8c) is much weaker (0.2°–0.5°C month−1) than in WM99 (0.5°–2°C month−1) and has a secondary maximum warming during boreal spring (Fig. 8c). The weak warming in Qυ at 110° and 140°W is partly due to biases in the mean meridional currents (Figs. 2g,h) and partly to weak TIWs. Biases in the low-frequency meridional currents may be associated with biases in wind forcing and biases that result from assimilating synthetic salinity rather than the observed salinity (Huang et al. 2008).
Heating by Qq (Fig. 8d) is the largest in the eastern Pacific owing to a shallow MLD (Fig. 1b) and is dominated by a semiannual cycle at all four sites. The Qq in GODAS agrees well with that in WM99. One noticeable disagreement is that the second maximum at 110°W during boreal fall is about twice as large as in WM99. The differences may result from different methods in calculating Qq, although the adjusted net heat flux (Fig. 5g) and MLD (Fig. 1b) in GODAS agrees very well with those of WM99. The climatology of Qq in GODAS is calculated with the pentad Qq in the period 1982–2004 and therefore retains the nonlinear relationship between the adjusted net heat flux and MLD [Eq. (A9)]. In contrast, the climatology of Qq in WM99 is based on the climatology of the adjusted net heat flux divided by the climatology of MLD. It should also be noted that there is considerable uncertainty in the mean seasonal cycle of net heat flux itself (Wang and McPhaden 2001b).
The combined cooling by Qw + Qzz (Fig. 8e) is the largest in the eastern Pacific due to the shallow MLD and strong upwelling in the region. The Qw + Qzz generally agrees well with those in WM99. However, cooling at 110°W in April–August is significantly underestimated in GODAS. The cooling at 140°W is also weaker in GODAS than in WM99. The weaker cooling by entrainment and vertical diffusion in GODAS may be associated with biases in zonal wind stress in R2, which is too weak when compared with Quick Scatterometer (QuikSCAT) observations and other products (Josey et al. 2002). The treatment of the vertical entrainment and vertical diffusion as a residual in WM99 may also contribute to the difference between GODAS and WM99.
Mitchell and Wallace (1992) suggested that the seasonal cycle of SST and its westward propagation are driven by a weakened vertical mixing from December to March due to a weakened meridional wind. Xie (1994) further proposed that the weakened vertical mixing largely result from the coupling of SST, meridional wind, and evaporation. KRC98 suggested that both net heat fluxes and vertical mixing contributed to the warming tendency from late winter to early spring. They also pointed out that the solar radiation in spring is larger than that in fall owing to a minimum in cloud cover in spring. WM99 suggested that net surface heat fluxes and residual subsurface fluxes, equivalent to the combination of vertical entrainment and vertical diffusion, are the two dominant terms and tend to cancel out each other during spring. They emphasized a large warming due to TIWs during fall/winter and a correlation between the mean annual cycle of residual subsurface fluxes and zonal winds, which implied that the subsurface fluxes are the weakest during spring when zonal winds are the weakest. Vialard et al. (2001) suggested the annual cycle of SST in the Niño-3 region is largely controlled by net surface heat fluxes.
Compared to aforementioned results, the annual cycle heat budget in GODAS has several shortcomings. First, the TIW warming in boreal fall is severely underestimated in GODAS (Fig. 7f). Second, the low-frequency meridional advection, Qυ, is too strong in boreal summer/fall (Fig. 7b). The biases in Qυ are probably due to too strong cross-equatorial winds in the NCEP reanalysis during summer/fall (not shown) that are used to force the MOM3 model. Third, the Qw + Qzz in boreal summer/fall appears too weak compared to that in WM99, KRC98, and Vialard et al. (2001) (not shown).
Despite these shortcomings, GODAS analysis shows that a positive temperature tendency in the tropical eastern Pacific from December to March is associated with the strengthening of heating by the adjusted heat flux, consistent with WM99, KRC98, and Vialard et al. (2001). This term dominates over the intensification in cooling by vertical mixing, both of which are associated with a shoaling mixed layer depth from December to March (Fig. 3b). The shoaling mixed layer depth is also closely related to surface heat fluxes and surface winds. Therefore, a weakening surface winds, strengthening net surface heat fluxes, and shoaling mixed layer depth all contribute to the warming tendency of SST from December to March.
d. El Niño composite of mixed layer heat budget
In the previous section, the annual mean and the seasonal cycle of the mixed layer heat budget based on GODAS pentad output were discussed. Another important feature of coupled variability in the equatorial tropical Pacific is ENSO. The ENSO phenomenon (Philander 1990) is the most prominent interannual climate signal, so the ability of GODAS to capture the variability of the mixed layer heat budget during the ENSO cycle is assessed next.
1) El Niño composite
As defined by the NOAA official ENSO index, the Oceanic Niño Index (ONI), there were eight warm episodes from 1979 to the present that occurred in 1982–83, 1986–88, 1991–92, 1994–95, 1997–98, 2002–03, 2004–05, and 2006–07. Many studies have shown that ENSO episodes are phase locked to the seasonal cycle (Chang et al. 1995; Tziperman et al. 1998; Galanti and Tziperman 2000; McPhaden and Zhang 2009). The onset of ENSO events generally occurs early in the year and becomes mature late in the year or early in the next year (Rasmusson and Carpenter 1982) with the exception of the 1986–88 event, whose onset and mature phases were different from others.
An El Niño composite is constructed using five weak-to-moderate events (1991–/92, 1994–95, 2002–03, 2004–05, and 2006–07). The 1986–87 event was excluded because of its unique onset and decay phases. The 1982–83 and 1997–98 events are strong events in which nonlinear terms are likely very different from those during weak-to-moderate events (Jin et al. 2003). In addition, the analysis suggests that the closure of the heat budget during the two events is not as good as during other El Niño events. We therefore compare the heat budget of the two events with that discussed in literature separately and assess the feasibility of GODAS in simulating the two strongest events in the past 30 years.
2) Heat budget for the composite El Niño
Shown in Fig. 9 are the mixed layer heat budgets for the El Niño composite averaged between 1°S and 1°N in GODAS. The composite El Niño starts in January of an El Niño (year 0) and ends in a following La Niña (year +1, Fig. 9h). The analysis indicates that heating by entrainment and vertical diffusion (Qw + Qzz, Fig. 9c), meridional advection (Qυ, Fig. 9b), and zonal advevction (Qu, Fig. 9a) contribute to El Niño development during year 0: The heating by Qw + Qzz and Qυ is the strongest east of 140°W; and the heating by Qu is distributed over most of the tropical Pacific between 160°E and 110°W. The adjusted surface heat flkux (Qq) is mostly opposite to the composite SST anomalies (Fig. 9h) and strongly acts to damp El Niño in the eastern tropical Pacific in year 0 and early year +1. TIW heating (Fig. 9e) also acts weakly to damp El Niño development.
El Niño transits to La Niña around April of year +1. A moderate cooling (−0.5° to approximately −1°C month−1) by Qu and Qq plays a critical role during the period prior to the transition, while both Qw + Qzz and Qυ still contribute to a heating in the central-eastern Pacific. Therefore, zonal advection leads the transition from El Niño to La Niña.
After the decay of El Niño, cooling by Qw + Qzz, Qυ, and Qu becomes strong, leading to La Niña development during the summer/fall of year +1. Strong heating by Qq acts again to damp La Niña development later in year +1.
The balance between temperature tendency (Fig. 9f) and forcing (Fig. 9g) appears to be good in the central-eastern tropical Pacific, but it is worse in the far eastern tropical Pacific largely due to biases in zonal currents (Fig. 1c) and the underestimation of TIW heating (Fig. 6f), which is anomalously weak (strong) during El Niño (La Niña) events (Yu and Liu 2003). The underestimation of TIW heating in GODAS largely resulted from applying a 4-week data assimilation window that smoothes out variations shorter than 30 days. The TIW signal is further weakened by using pentad outputs on the 1° × 1° grid, which is available for the analysis.
The roles of different oceanic processes and surface heat fluxes in the El Niño composite in GODAS is further discussed for the Niño-3.4 SST index region (5°S–5°N, 120°–170°W). Figure 10a shows the spatial average of each flux term in the Niño-3.4 region. Lee et al. (2004) suggested that spatially integrated local temperature advection can be dominated by internal processes such as TIW that redistribute heat within the domain and they recommended a new boundary flux approach to analyze the mixed layer heat budget. Kim et al. (2007) further showed that the spatially integrated zonal advection tendency in the Niño-3 region derived from ECCO analysis is anticorrelated with the change of SST during the 1997–99 El Niño–La Niña cycle, which disagrees with the common understanding about the role of large-scale advection of warm pool water during El Niño. Ultimately, the two approaches must be consistent with one another when all terms in the heat balance are accounted for and properly interpreted, so our preference is to perform the analysis based on the more conventional local balance approach.
In this analysis the advective temperature terms are further decomposed into various components described in Eqs. (B5)–(B8). It is seen in Fig. 10a that the positive temperature tendency during the spring of year 0 is due to heating from Qw + Qzz, Qυ, and Qu. The temperature tendency changes sign from positive to negative at the end of year 0, which coincides with a rapid transition of Qu from heating to cooling. It is interesting that Qv does not have a transition until 4–5 months later. This suggests that cooling from Qu and Qq plays a dominant role during the decay phase of El Niño and transition to the following La Niña, consistent with ENSO recharge–discharge theory (Jin 1997). During the summer and fall of year +1, cooling from Qυ and Qu all contribute to the transition from El Niño to La Niña as well as the development of La Niña. However, Qw + Qzz does not appear to contribute to the development of La Niña in weak-to-moderate events, which is not the case for strong events such as in 1997–98 and 1982–83, which will be discussed later.
The decomposition of Qu in Fig. 10b shows that the heating by Qu during El Niño development in year 0 is dominated by −u′Tx (subscripts x, y, and z represent a partial derivative, hereafter), and the cooling by Qu during La Niña development in year +1 is also associated with −u′Tx. This suggests an important role of the zonal advection of climatological temperature by anomalous zonal currents in El Niño and La Niña development. Note that −uT ′x acts to weakly damp El Niño and TIW acts to weakly damp both El Niño and La Niña. For meridional advection, the heating (cooling) by Qυ (Fig. 10c) during the development of El Niño in year 0 (La Niña in year +1) in GODAS results from both −υT ′y and −υ′Ty that have comparable amplitude.
The decomposition of Qw + Qzz in Fig. 10d shows that both −wT ′z and −w′Tz contribute to the development of El Niño and their amplitudes are comparable. During the onset phase of El Niño, both −w′Tz and −wT ′z are important, suggesting that the role of upwelling anomalies (−w′Tz) is as important as subsurface temperature anomalies (−wT ′z) when El Niño is beginning to develop. During the decay phase of El Niño and development phase of La Niña, all components in Qw + Qzz act as a weak warming (Fig. 10d), which is quite different from those during the strong events of 1997–98 and 1982–83.
Base on GODAS, the onset phase of El Niño can be attributed to the heating by Qu, Qυ, and Qw + Qzz, which is largely consistent with observational analysis in Wang and McPhaden (2000). Since Qu is dominated by −u′Tx in GODAS (Fig. 10b), anomalous zonal currents acting on the large climatological SST gradient near the eastern edge of the warm pool plays a critical role during the onset phase of El Niño, a feature also noted in previous observational studies (Frankignoul et al. 1996; Picaut et al. 1996; Wang and McPhaden 2000). The results based on a coupled model simulation also suggested the importance of Qu in El Niño development (Zhang et al. 2007).
Meridional advection (Qυ) contributes to the development of El Niño east of 140°W in GODAS (Fig. 9b). The decomposition of Qυ indicates that both −υT ′y and −υ′Ty make important contributions to El Niño development in GODAS, while coupled model simulation (Zhang et al. 2007) suggested that −υ′Ty is not as important as −υT ′y. Wang and McPhaden argued that Qυ is a damping since TIW dominates this term at the equator and generally acts as a damping to SST, as indicated in model simulations (Yu and Liu 2003; ZH08; An 2008). This is also true in our analysis of GODAS in the equatorial Pacific (Figs. 10b,c), but the TIW term is underestimated by GODAS (Fig. 8c).
The decomposition of Qw + Qzz (Fig. 10d) shows that both wT ′z and w′Tz contribute to the onset phase of El Niño in GODAS. This suggests that both w′ and subsurface temperature anomalies (T ′) are key variables to determine El Niño onset, consistent with other observational studies (Wang and McPhaden 2000, 2001a; Zhang and McPhaden 2006, 2008). The separation of Qw + Qzz (not shown) indicates that Qzz plays an important role in the heat budget of the El Niño composite in GODAS. The observational analyses of 1987–88 (Hayes et al. 1991) and 1991–93 (Kessler and McPhaden 1995) El Niño events, model analysis of the 1997–98 event (Kim et al. 2007), and a coupled model simulation (Zhang et al. 2007) all suggested the important role of Qzz.
During the decay phase of El Niño, both Qu and Qυ play a leading role in GODAS. In fact, Qu leads Qυ by 4–5 months (Fig. 10a). Zonal advection (Qu) is dominated by u′Tx, which suggests that anomalous zonal currents are a precursor for the transition from El Niño to La Niña. These results are very consistent with previous studies based on observations (Wang and McPhaden 2000, 2001a; Vialard et al. 2001) and with coupled model simulations (An and Jin 2001; Zhang et al. 2007). Observational studies of Wang and McPhaden (2000, 2001a) suggested that Qq acts as a damping to El Niño because of the out-of-phase relationship between SST and Qq. This damping is particularly important during the decay phase of El Niño (Figs. 10a and 9d).
The Qq anomaly is largely due to latent heat flux (Qlh, Fig. 11a) and its dependence on anomalies of SST and other near-surface atmospheric variables. Shortwave radiation (Qshort) contributes secondarily to cooling during the peak phase of El Niño. Penetrative solar radiation (Qpen) contributes to heating (cooling) when mixed layer depth is anomalously deep (shallow). Longwave radiation and sensible heat flux (Qlong + Qsh) are generally weak. The flux correction (Qcorr) resulting from relaxation to observed SST is relatively large during the development phase of El Niño and La Niña, indicating biases in surface heat fluxes and errors in the model physics.
For the two strong events of 1997–98 and 1982–83, the mixed layer heat budget closure is generally poor (Fig. 2d). We pointed out earlier that the imbalance between F and Tt may result from multiple factors, including sources and sinks of heat due to data assimilation, and estimation of vertical diffusion. It should also be noted that the evolution of mixed layer heat budget during strong El Niño events may differ significantly from that for weak-to-moderate El Niño events.
Next, we assess if it is feasible for GODAS to simulate the mixed layer heat budget of the 1997–98 and 1982–83 events, which are the biggest of the past 30 years. Figure 12a shows the GODAS budget at 0°, 110°W from 1997 to 1998. It is clear that the development of the El Niño in 1997 is largely associated with the anomalous heating by subsurface Qw + Qzz and zonal advection Qu. The transition from El Niño to La Niña results from the anomalous cooling by Qu and net Qq. The development of La Niña in 1998 is largely associated with the anomalous cooling by Qw + Qzz and Qu. These features are largely consistent with the heat budget composite for weak-to-moderate events (Fig. 10a). Differences are that anomalous cooling by Qw + Qzz played the dominant role in the development of La Niña in 1998, while Qw + Qzz acts as a weak damping in the composite. Meridional advection Qυ did not contribute much to the development of El Niño and acted as a damping to the development of the La Niña in 1998. However, Qυ contributed significantly to the development of the composite El Niño and La Niña (Fig. 10a). The roles of Qw + Qzz and Qq in our analysis are consistent with the observational analysis of Wang and McPhaden (2001a, their Fig. 6). However, the magnitude of Qw + Qzz and Qq appears to be underestimated in our analysis. Some differences were also found in Qu and Qυ between GODAS and observations. Zonal advection (Qu) contributed to the development of El Niño in GODAS in 1997, which is consistent with observations. However, it contributed significantly to the development of La Niña in GODAS in 1998 instead, being a weak damping as in the observations. Further analysis showed that the anomalous cooling by Qu in GODAS in 1998 is largely associated with −u′Tx (Fig. 12b) because of westward biases in zonal current anomalies (Fig. 12d). Meridional advection (Qυ) contributed moderately to the development of El Niño in 1997 in the observation, but it was negligible in GODAS (Figs. 12a,c). TIW acts as a weak damping to the El Niño in 1997 and La Niña in 1998 (Figs. 12b,c). Differences at the other three TAO mooring sites (140°, 170°W; and 165°E) are also large (not shown).
For the 1982–83 event (Fig. 13a), anomalous heating and cooling by Qw + Qzz played a dominant role in the developments of El Niño and La Niña in GODAS, which is very consistent with the observational analysis of Wang and McPhaden (2000, their Fig. 5), although the magnitude of these terms was weaker in GODAS than in observations. Anomalous cooling by Qq damped the El Niño in 1982, consistent with observations. Anomalous heating by Qu contributed to the development of El Niño in GODAS in 1982, which is also consistent with observations; however, it became too strong in 1983 compared with the observations. This is largely due to the nonlinear heating (mainly −u′T ′x, Fig. 13b) that appears to be associated with strong westward biases of zonal current in GODAS (Fig. 13d). In observations, anomalous cooling by Qυ damped the development of the El Niño in 1982 but contributed to the development of the La Niña in 1983. However, in GODAS Qυ is very weak in 1982 and is a moderate damping to the development of the La Niña in 1983 because of heating from −υ′Ty and the nonlinear term (Fig. 13c). The TIW acts as a weak damping to the El Niño in 1982 and La Niña in 1983 (Figs. 13b,c).
The above results suggest that GODAS simulates Qq and Qw + Qzz reasonably well, but contain large biases in Qu and Qυ, most of which are due to biases in surface currents. It is found that surface current analysis near the equator can be significantly improved when observed salinity is assimilated into GODAS and the balance between density and current is properly maintained in the analysis cycle (Behringer 2007). Therefore, we hope that biases in Qu and Qυ will be significantly reduced in the next version of GODAS in which the Argo salinity will be assimilated.
6. Summary and discussion
The focus of this paper is to assess the adequacy of GODAS pentad data to provide a consistent picture of the mixed layer heat budget in the tropical Pacific from 1979 to present. We have demonstrated that the mean mixed layer depth was reasonably well simulated by GODAS east of the date line but was about 10–20 m too deep in the western Pacific, probably because the water in the western Pacific warm pool was too salty (Behringer 2007).
The climatological mean and seasonal cycle of mixed layer heat budgets derived from GODAS agree reasonably well with previous observational (WM99) and model-based estimates (KRC98; Vialard et al. 2001). However, significant differences and biases were noticed. The net surface heat fluxes into the ocean derived from the NCEP R2 and used to force GODAS are about 40–60 W m−2 too low in the eastern equatorial Pacific during boreal spring (Fig. 5h). Large biases were also found in GODAS zonal and meridional currents (Tables 1, 2), which contributed to biases in the annual cycle of zonal and meridional advective heat fluxes, Qu and Qυ (Figs. 7a,b). However, the interannual variability of the mixed layer heat budget has been simulated reasonably well by GODAS. This is demonstrated using the composite of anomalous heat budget for five weak-to-moderate El Niño events (Figs. 9, 10). The results, therefore, suggest that GODAS provides a reasonable estimate of both the seasonal cycle and interannual variability in the mixed layer heat budget. We conclude, therefore, that it is a useful tool for real-time monitoring of mixed layer variability and for further understanding coupled ocean–atmospheric interactions.
The mixed layer heat budget analysis tool derived from GODAS can also be derived from other operational ocean reanalyses (ORAs) routinely produced by operational centers (Balmaseda et al. 2008; Oke et al. 2005; Alves and Robert 2005). The effort of deriving multimodel heat budget analyses would likely enhance our confidence in the feasibility of ORAs for monitoring and understanding tropical Pacific SST variability. Monitoring the differences among multimodel analyses will also help us identify common biases and deficiencies in ORAs and the requirements those ORAs impose on the ocean observing system for operational data streams (Xue et al. 2010). As the ORAs also provide initial conditions for seasonal climate predictions, the mixed layer heat budget analysis presented here can also be applied to understand the evolution of the SST anomalies in the forecast, as well as to understand the initial adjustments (sometimes referred to as initial shock) and development of coupled model biases.
We thank Dr. Kelvin Richards, Dr. Shangping Xie, and an anonymous reviewer for their critical and constructive comments, which helped us improve the manuscript significantly. We also want to thank Dr. David Behringer, who provided us detailed information about the GODAS ocean analysis and for helpful discussions on our results. This work was partially supported by the Climate Test Bed of NOAA and also by NOAA’s Climate Program Office.
The temperature equation for the mixed layer, described by Stevenson and Niiler (1983), is adopted here:
where h is MLD; T and v are temperature and horizontal currents; subscript a represents quantities that are vertically averaged between the surface and h; v̂ and T̂ in the third term on the right-hand side are deviations from their vertically averaged quantities, and their contribution to the heat budgets is generally very small in the mixed layer (ZH08) and neglected in our discussions; we is the entrainment velocity across the bottom of the mixed layer, calculated by
where w−h is the vertical velocity at the bottom of the mixed layer. The last four terms in Eq. (A1) are various heat flux terms: QO is the net downward surface heat fluxes, composed of shortwave (Qshort) and longwave (Qlong) radiation, latent (Qlh) and sensible (Qsh) heat fluxes; Qpen is the shortwave radiation transmitted through the bottom of MLD and is parameterized (Pacanowski and Griffies 1999) as
Here Qcorr represents a heat flux correction introduced by a relaxation of model SST to observed SST in GODAS. Following WM99, we define the adjusted net surface heat flux as
in which ρ and cp are density and heat capacity of seawater;
in which Kc = 1 × 10−5 m2 s−1, β = 5, and Ri is the Richardson number defined as
where α = 8.75 × 10−6(T + 9)°C−1 and g = 9.8 m s−2; a minimum Ri is set to −0.1 to ensure valid diagnosis; and ν in Eq. (A6) is parameterized as
where νc = 1 × 10−4 and ν0 = 3.5 × 10−3 m2 s−1. We should point out that the diffusion parameterization used in our analysis in (A5)–(A8) is much simpler than originally implemented in GODAS [the K-profile parameterization (KPP) vertical mixing scheme, Large et al. (1994)], which may result in the imbalance of heat budgets. Horizontal diffusion has been neglected owing to its small magnitude.
We also note that vertical entrainment and diffusion are sometimes treated as a residual to keep the closure of mixed layer heat budgets (WM99). They were often combined and parameterized to be proportional to vertical differences between the temperatures in and below the mixed layer. One common parameterization is to assume a constant temperature difference of 1°–4°C (see reference cited in WM99). An alternative parameterization is to assume a depth from which colder water is entrained into the mixed layer. The entrainment depth of 0–20 m below the mixed layer is used in many studies (Stevenson and Niiler 1983; WM99; McPhaden et al. 2008).
Following Zhang et al. (2007), the vertical entrainment is rewritten as
The Eq. (A1) then becomes
Equation (A9) can be rewritten as
where Tt = ∂Ta/∂t is referred to as the temperature tendency and F as the forcing. The forcing is the sum of zonal advection: Qu = −u∂Ta/∂x, meridional advection: Qυ = −υ∂Ta/∂y, vertical entrainment: (Qw = −we∂T/∂z), adjusted surface heat flux: Qq = Qadj/(ρcph), and vertical diffusion: Qzz = Qdiff/(ρcph). The consistency between Tt and F can be used to check the closure of the temperature equation.
Decomposition of the Temperature Equation
Following KRC98, each variable associated with forcing F is decomposed into its low frequency (≥75 day) and high frequency components; for example, u = uL + uH. With this decomposition by omitting the subscript a, Eq. (A9) is written as
where superscript L indicates the term calculated using low-pass filtered variables:
, and ; R in Eq. (B1) indicates the residue terms that involve high-frequency variations:
The residual terms are combined as an “eddy” term:
therefore, Eq. (B1) becomes
Equation (B3) is further decomposed into climatology (bar) and its anomaly (prime). The equation for anomalous temperature is, omitting superscript L,
We combine vertical diffusion into the entrainment in Eq. (B4) because both Qw and Qzz are parameterized to be proportional to the vertical gradient of temperature. The vertical diffusion
is rewritten as Qzz = −ω∂T/∂z. Here ω = Kz/h represents an equivalent entrainment velocity and can be decomposed into its climatology and anomaly: ω = ω + ω′. Therefore, Q′zz = −ω∂T ′/∂z − ω′∂T/∂z − ω′∂T ′/∂z + ω′∂T ′/∂z. The combination of anomalous vertical entrainment and vertical diffusive heat fluxes becomes
Lastly, the anomalous heat by eddy is calculated as
Corresponding author address: Dr. Boyin Huang, Wyle Information Systems and NOAA/Climate Prediction Center, 5200 Auth Road, Room 605-A WWB, Camp Springs, MD 20746. Email: email@example.com