Global upper ocean heat content and climate variability

Abstract

Observational data from the Global Temperature and Salinity Profile Program were used to calculate the upper ocean heat content (OHC) anomaly. The thickness of the upper layer is taken as 300 m for the Pacific/Atlantic Ocean and 150 m for the Indian Ocean since the Indian Ocean has shallower thermoclines. First, the optimal spectral decomposition scheme was used to build up monthly synoptic temperature and salinity dataset for January 1990 to December 2009 on 1° × 1° grids and the same 33 vertical levels as the World Ocean Atlas. Then, the monthly varying upper layer OHC field (H) was obtained. Second, a composite analysis was conducted to obtain the total-time mean OHC field ($\bar{\bar{H}}$) and the monthly mean OHC variability ($\widetilde{\text{H}}$), which is found an order of magnitude smaller than $\widehat{\text{H}}$. Third, an empirical orthogonal function (EOF) method is conducted on the residue data ($\widehat{\text{H}}$), deviating from $\widetilde{\text{H}}$ + $\widetilde{\text{H}}$, in order to obtain interannual variations of the OHC fields for the three oceans. In the Pacific Ocean, the first two EOF modes account for 51.46% and 13.71% of the variance, representing canonical El Nino/La Nina (EOF-1) and pseudo-El Nino/La Nina (i.e., El Nino Modoki; EOF-2) events. In the Indian Ocean, the first two EOF modes account for 24.27% and 20.94% of the variance, representing basin-scale cooling/warming (EOF-1) and Indian Ocean Dipole (EOF-2) events. In the Atlantic Ocean, the first EOF mode accounts for 49.26% of the variance, representing a basin-scale cooling/warming (EOF-1) event. The second EOF mode accounts for 8.83% of the variance. Different from the Pacific and Indian Oceans, there is no zonal dipole mode in the tropical Atlantic Ocean. Fourth, evident lag correlation coefficients are found between the first principal component of the Pacific Ocean and the Southern Oscillation Index with a maximum correlation coefficient (0.68) at 1-month lead of the EOF-1 and between the second principal component of the Indian Ocean and the Dipole Mode Index with maximum values (around 0.53) at 1–2-month advance of the EOF-2. It implies that OHC anomaly contains climate variability signals.

Appendix A: Optimal spectral decomposition

A new data analysis/assimilation scheme, the optimal spectral decomposition (OSD), has been recently developed to analyze fields from noisy and sparse oceanographic data. Let (x, z) be horizontal and vertical coordinates and t be time. A physical variable $c\left( {{\mathbf{x}},z,t} \right)$ at depth ${z_k}$ is decomposed using the generalized Fourier series (Chu 1999; Chu et al. 2003a, b; Chu et al. 2005a, b)

$$ c\left( {{\mathbf{x}},{z_k},t} \right) = {c_0}\left( {{z_k},t} \right) + \sum\limits_{{m = 1}}^M {{a_m}\left( {{z_k},t} \right){ }} {\Phi_m}\left( {{\mathbf{x}},{z_k}} \right),{ }{\mathbf{x}} \in R\left( {{z_k}} \right) $$

(A1)

where c ₀ is the horizontal mean of the variable c. M is the truncated mode number, and ${\Phi_m}\left( {{\mathbf{x}},{z_k}} \right)$ and ${a_m}\left( {{z_k},t} \right)$ are the orthogonal basis functions (or called modes) and the spectral coefficients, respectively; $R\left( {{z_k}} \right)$ is the area bounded by the lateral boundary $\Gamma \left( {{z_k}} \right)$ at depth ${z_k}$. The basis functions $\left\{ {{\Phi_m}\left( {{\mathbf{x}},{z_k}} \right)} \right\}$ are eigen-functions of the horizontal Laplace operator $\left( {\nabla_h^2 \equiv {\partial^2}/\partial {x^2} + {\partial^2}/\partial {y^2}} \right)$ with the basin geometry and certain physical boundary conditions. For temperature and salinity, the homogeneous Neumann boundary condition is taken at the solid boundary $\Gamma (z)$ (i.e., no heat and salt fluxes),

$$ \nabla_h^2{\Phi_m} = - {\lambda_m}{\Phi_m},{ }{\mathbf{n}} \bullet {\nabla_h}{\Phi_m}{|_{\Gamma }} = 0,{ }m = 1,2,...,M, $$

(A2)

where n is the unit vector normal to $\Gamma (z)$. The basis functions $\left\{ {{\Phi_m}} \right\}$ are independent of the data and therefore available prior to the data analysis. The OSD method has two important procedures: optimal mode truncation and determination of spectral coefficients {a _m }. The optimal mode truncation (M) is determined using Vapnik (1983) variational method. Application of the generalized Fourier series expansion (1) to the observational points with P as the total number of observations leads to an algebraic equation

$$ {\mathbf{Aa}} = {\mathbf{QY}} $$

(A3)

where a = (a ₁, a ₂, …, a _M) is the state vector (M-dimensional); A is a P × M matrix; Q is a P × P square matrix (P > M); Y is a P-dimensional observation vector, consisting of a signal ${\mathbf{\bar{Y}}}$ and a noise Y’. Due to high level of noise contained in the observations, the algebraic equation (A3) is ill-posed and needs to be solved by a rotation matrix regularization method (Chu et al. 2004) that provides (a) stability (robustness) even for data with high noise and (b) the ability to filter out errors with a priori unknown statistics. Interested readers are referred to Chu et al. (2003a, b, 2004).