@prefix vivo: . @prefix edm: . @prefix ns0: . @prefix dcterms: . @prefix skos: . vivo:departmentOrSchool "Science, Faculty of"@en, "Earth and Ocean Sciences, Department of"@en ; edm:dataProvider "DSpace"@en ; ns0:identifierCitation "Tang, Youmin, Hsieh, WIlliam W. 2003. Nonlinear modes of decadal and interannual variability of the subsurface thermal structure in the Pacific Ocean. Journal of Geophysical Research Oceans. 108 C3 3084"@en ; ns0:rightsCopyright "Hsieh, Wiliam W."@en ; dcterms:creator "Tang, Youmin"@en, "Hsieh, William W."@en ; dcterms:issued "2016-11-10T16:18:07"@en, "2003-03-18"@en ; dcterms:description "The nonlinear principal component analysis, a neural network technique, is applied to the observed upper ocean heat content anomalies (HCA) in the Pacific basin from 1961 to 2000. By applying the analysis to high-passed and low-passed data, nonlinear interannual and decadal modes are extracted separately. The first nonlinear interannual mode is mainly characterized by the El Niño-Southern Oscillation (ENSO) structure in the tropical Pacific, with considerable asymmetry between warm El Niño and cool La Niña episodes; for example, during strong El Niño, the negative HCA in the western tropical Pacific is much stronger than the corresponding positive HCA during strong La Niña. The first nonlinear decadal mode goes through several notable phases. Two of the phases are related to decadal changes in the La Niña and El Niño characteristics, revealing that the decadal changes for La Niña episodes are much weaker than the changes for El Niño episodes. Other phases of the decadal mode show a possible anomaly link from the middle latitudes to the western tropical Pacific via the subtropical gyre. The decadal changes in the HCA around 1980 and around 1990 were compared and contrasted. An edited version of this paper was published by AGU. Copyright 2003 American Geophysical Union."@en ; edm:aggregatedCHO "https://circle.library.ubc.ca/rest/handle/2429/32827?expand=metadata"@en ; skos:note "Nonlinear modes of decadal and interannual variability of the subsurface thermal structure in the Pacific Ocean Youmin Tang Courant Institute of Mathematical Sciences, New York University, New York, New York, USA William W. Hsieh Department of Earth and Ocean Sciences, University of British Columbia, Vancouver, British Columbia, Canada Received 27 November 2001; revised 19 November 2002; accepted 17 December 2002; published 18 March 2003. [1] The nonlinear principal component analysis, a neural network technique, is applied to the observed upper ocean heat content anomalies (HCA) in the Pacific basin from 1961 to 2000. By applying the analysis to high-passed and low-passed data, nonlinear interannual and decadal modes are extracted separately. The first nonlinear interannual mode is mainly characterized by the El Niño-Southern Oscillation (ENSO) structure in the tropical Pacific, with considerable asymmetry between warm El Niño and cool La Niña episodes; for example, during strong El Niño, the negative HCA in the western tropical Pacific is much stronger than the corresponding positive HCA during strong La Niña. The first nonlinear decadal mode goes through several notable phases. Two of the phases are related to decadal changes in the La Niña and El Niño characteristics, revealing that the decadal changes for La Niña episodes are much weaker than the changes for El Niño episodes. Other phases of the decadal mode show a possible anomaly link from the middle latitudes to the western tropical Pacific via the subtropical gyre. The decadal changes in the HCA around 1980 and around 1990 were compared and contrasted. INDEX TERMS: 3339 Meteorology and Atmospheric Dynamics: Ocean/atmosphere interactions (0312, 4504); 4215 Oceanography: General: Climate and interannual variability (3309); 4522 Oceanography: Physical: El Niño; 1620 Global Change: Climate dynamics (3309); KEYWORDS: nonlinear principal component analysis, heat content, El Niño Citation: Tang, Y., and W. W. Hsieh, Nonlinear modes of decadal and interannual variability of the subsurface thermal structure in the Pacific Ocean, J. Geophys. Res., 108(C3), 3084, doi:10.1029/2001JC001236, 2003. 1. Introduction [2] Among the low-frequency variability of the thermal fields in the Pacific Ocean, interannual variability and decadal variability are the two most interesting [e.g., Wal- lace et al., 1998; Trenberth and Hurrell, 1994]. While these two well-defined variabilities reside in the whole Pacific basin within at least the upper 400-m ocean, they also show strong regional features. The interannual variability, domi- nated by the El Niño-Southern Oscillation (ENSO) phe- nomenon, is centered in the equatorial Pacific, whereas the decadal variability is most strongly manifested in the mid- latitude North Pacific, as characterized by an elliptical anomaly located in the subtropic gyre [Zhang et al., 1999]. Understanding and interpreting the interannual and decadal variabilities have long been of interest [e.g., Klee- man et al., 1996, 1999], not only for their major impacts on the regional and global climates and ecologies, but also for assessing possibly forced climate variability, such as anthro- pogenic global warming [Latif et al., 1997]. [3] An important aspect of studying the low-frequency variability in the Pacific Ocean is to characterize the major spatial and temporal characteristics in a low-dimensional space. Until very recently, this has been implemented by principal component analysis (PCA, also called EOF anal- ysis), and by related techniques, for example, singular spectrum analysis (SSA, also called extended EOF analy- sis), and principal oscillation pattern (POP) analysis, with either observed data [Zhang et al., 1999] or modeled data [Miller et al., 1998]. The interannual and decadal modes are described by the first few leading eigenvectors, giving the spatial patterns, and by the corresponding time series. To focus on a specific timescale, the data are usually filtered prior to applying PCA. For instance, for detecting decadal variability, we used a filter which removes signals with periods under 5 years, while for studying interannual variability, we filtered out periods above 5 years. The leading interannual and decadal PCA modes (Figures 1a and 1b) characterize the spatial anomaly patterns at different frequency oscillations (Figure 2a and 2b). [4] In this paper, a nonlinear algorithm to extract low- dimensional structure from multivariate data sets, i.e., nonlinear principal component analysis (NLPCA), is applied to the oceanic heat content anomalies in the upper JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 108, NO. C3, 3084, doi:10.1029/2001JC001236, 2003 Copyright 2003 by the American Geophysical Union. 0148-0227/03/2001JC001236$09.00 29 - 1 400 m (HCA) over the Pacific basin to detect nonlinear modes of decadal-scale and interannual variability. There is no a priori reason to believe that the thermal structures in Pacific Ocean are linear. As the data contain nonlinear lower-dimensional structure, the PCA will miss the non- linearity. Compared with the sea surface temperature, the upper ocean heat content is better for describing and understanding interannual and decadal variability [Zhang et al., 1999], as it reflects the thermocline displacement and contains the ocean’s ‘‘memory.’’ NLPCA was developed originally by Kramer [1991] in the chemical engineering literature, was applied to the Lorenz three-component chaos system by Monahan [2000], and to several meteoro- logical and oceanographic data sets [Monahan, 2001; Monahan et al., 2001; Hsieh, 2001; Hamilton and Hsieh, 2002]. [5] This paper is structured as follows: Section 2 briefly describes the methodology and the data. Section 3 presents the nonlinear interannual mode, section 4 presents the nonlinear decadal mode, section 5 presents the decadal changes in the 1980s and the 1990s, and Section 6 is the summary and conclusion. 2. Method and Data 2.1. NLCPA [6] If the data are in the form x(t) = [x1,. . ., xl], where each variable xi, (i = 1,. . .,l ), is a time series containing n observations, the PCA method looks for u, a linear combi- nation of the xi, and an associated vector a, with u tð Þ ¼ a  x tð Þ; ð1Þ so that hkx tð Þ  au tð Þk2i is minimized; ð2Þ where h...i denotes a sample or time mean. Here u, called the first principal component (PC), is a time series, while a, the first eigenvector of the data covariance matrix (also called an empirical orthogonal function, EOF), often describes a spatial pattern. [7] The fundamental difference between NLPCA and PCA is that NLPCA allows a nonlinear mapping from x to u whereas PCA only allows a linear mapping. To perform NLPCA, a nonlinear mapping is made; that is, u tð Þ ¼ f x tð Þ;wð Þ; ð3Þ where f denotes the nonlinear mapping function from the data space to the u (the nonlinear PC) space, and w denotes the parameters determining the f structure inherent to the data set. Denoting g as the inverse mapping function from u to the data space, we have x0 tð Þ ¼ g u; ~wð Þ; ð4Þ where g is the f-adjoint operator. For linear PCA, g is simply the transpose of f. Here x0(t) is the approximation to data set x(t), when the 1-D PC space is used to describe the data set. Figure 1. EOF1 of the HCA data for (a) the high-passed data (i.e., with the 61-month running mean subtracted from the original data) and (b) the low-passed data (i.e., the 61-month running mean). The value in percentage is the explained variance by each mode. Contour interval is 0.2C, with dashed contours for negative anomalies. Figure 2. First mode PC associated with the EOF spatial patterns in Figure 1. For better legibility, the PCs for different data sets have been shifted vertically by 0.25. The tick marks along the abscissa indicate the start of the year. 29 - 2 TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN As in linear PCA, the cost function defined by the error between x(t) and x0(t) is used to determine the parameters w and ~w; that is, hkx tð Þ  x0 tð Þk2i is minimized: ð5Þ [8] An important issue in NLPCA is how to derive the nonlinear operators f and g from the inherent structure of the data set. This has been implemented by neural networks (NN) [Kramer, 1991], since NN can simulate any nonlinear continuous functions [Cybenko, 1989]. Figure 3a shows the architecture of the NLPCA, which is capable of extracting a 1-D open curve approximation to the data. However, this algorithm cannot be used to extract closed curve solutions, as the bottleneck neuron u is not an angular variable. Kirby and Miranda [1996] introduced a circular node or neuron, and showed that the NLPCA with a circular node (henceforth abbreviated as NLPCA.cir) at the bottleneck is capable of extracting closed curve solutions. The algorithm of the NLPCA.cir is identical to the architecture of the NLPCA of Kramer, except at the bottleneck layer, where instead of a single neuron u, there are now two neurons p and q, con- strained to lie on a unit circle in the p-q plane (Figure 3b), so there is only 1 angular degree of freedom (q) to present the nonlinear PC (NLPC). In this paper, both NLPCA and NLPCA.cir algorithms are used. When we discuss the decadal mode, we use NLPCA.cir, since the analyzed data, obtained by smoothing the original data set with a low-pass filter, are well characterized by closed curve solutions. [9] In contrast to PCA, as the mapping function g from the PC space to the data space is nonlinear, there is not a single spatial pattern associated with an NLPCAmode. The approx- imation x0(t), however, corresponds to a sequence of different patterns that can be visualized cinematographically. For linear PCA, the approximation au (equation (2)) produces a standing wave pattern as the PC varies, whereas with NLPCA the spatial pattern generally changes as the NLPC varies. We will use the x0(t) corresponding to a few u (q) values to explore the changing spatial structures of the NLPCAmodes. [10] An important aspect of the NLPCA is the size of the network, i.e., the number of hidden neurons m in the encoding (and also in the decoding layer) for representing the nonlinear functions f and g. A larger m increases the nonlinear modeling capability of the network, but could also lead to overfitted solutions (i.e., wiggly solutions which fit to the noise in the data). Based on a general principle of parsimony, the m values were varied from 2 to 4 and the weight penalty parameters [Hsieh, 2001] were varied from 0.01 to 0.05 for smoothing. For a given m, an ensemble of 30 NNs with random initial weights and bias parameters was run. Also, 20% of the data was randomly selected as test data and withheld from the training of the NNs. Runs where the mean square error (MSE) was larger for the test data set than for the training data set were rejected to avoid overfitted solutions. The NN with the smallest MSE was selected as the solution for the given m. The solutions from different m were further compared with respect to their MSE to get the optimal NN structure. 2.2. Data [11] The data used are the monthly 400-m depth-averaged heat content anomalies (HCA) during 1961–2000, from the data set of subsurface temperature and heat content pro- vided by the Joint Environmental Data Analysis Center at the Scripps Institution of Oceanography. This data set consists of all available XBT, CTD, MBT and hydrographic observations, optimally interpolated by White [1995] to a three-dimensional grid of 2 latitude by 5 longitude, and 11 standard depth levels between the surface and 400 m. This data set has recently been successfully assimilated into a Figure 3. (a) A schematic diagram of the NN model for calculating nonlinear PCA (NLPCA). There are three ‘‘hidden’’ layers of variables or ‘‘neurons’’ (denoted by circles) sandwiched between the input layer x on the left and the output layer x0 on the right. Next to the input layer is the encoding layer, followed by the ‘‘bottleneck’’ layer (with one neuron u), which is then followed by the decoding layer. A nonlinear function maps from the higher dimension input space to the lower dimension bottleneck space, followed by an inverse transform mapping from the bottleneck space back to the original space represented by the outputs, which are to be as close to the inputs as possible by minimizing the cost function J = hkxx0k2i. Data compression is achieved by the bottleneck, with the bottleneck neuron giving u, the nonlinear principal component (NLPC). (b) A schematic diagram of the NN model for calculating the NLPCAwith a circular node at the bottleneck (NLPCA.cir). Instead of having one bottleneck neuron u, there are now two neurons p and q constrained to lie on a unit circle in the p-q plane, so there is only one free angular variable q, the NLPC. This network is suited for extracting a closed curve solution. TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN 29 - 3 hybrid coupled model for ENSO prediction [Tang and Hsieh, 2003], and used in the study of decadal oscillations [e.g., Miller et al., 1997, 1998; Schneider et al., 1999]. [12] To study the decadal mode, the data were first smoothed by a 61-month running mean (referred to as the low-passed data hereinafter). The residual field between the original data and the low-passed data (referred to as the high- passed data) will be used to extract the interannual mode. To reduce the large number of spatial variables, the HCA data were preprocessed by retaining only the first six EOF modes, which account for 41% and 93% of the variance for the high- passed data and low-passed data, respectively. 3. Interannual Mode [13] The six leading PCs from the high-passed HCA are input to the NLPCA network to extract the NLPCA mode 1 (NLPCA1). Figure 4 shows the projection of the NLPCA1 solution in the PC1-PC2-PC3 space. The NLPCA1 accounts for 26% of the total variance versus 22% by the PCA mode 1. The trajectory of the NLPCA1 describes a curve in the PC space, indicating nonlinearity as compared to the PCA (straight line). The NLPC, u, time series is shown in Figure 5a, well characterized by irregular oscillations at 2- to 5-year timescale, while Figure 5b is the frequency distribu- tion curve (FDC) for u. We next examine the spatial anomaly patterns associated with some specific u values, namely those marked in Figure 5b. The neural network maps from u to the output PCs (x0), which when individu- ally multiplied to the associated EOF spatial pattern, and summed over the six modes, yield the spatial anomaly pattern of the NLPCA1 for the given u. As shown in Figure 6, the spatial structures of this nonlinear interannual mode are mainly characterized by ENSO features in the tropical Pacific, i.e., a seesaw oscillation along the equator. The most probable spatial pattern, corresponding to C in Figure 5b, describes a neutral state, i.e., negligible anomalies in the tropical Pacific (not shown). Patterns A and B depict extreme and typical La Niña episodes, respectively, while D and E represent typical and extreme El Niño, respectively (Figure 6). In the middle latitude, the interannual variability Figure 4. The first NLPCA mode for the high-pass filtered HCA plotted as (overlapping) squares in the PC1-PC2-PC3 3-D space. The linear (PCA) mode is shown as a dashed line. The NLPCA mode and the PCA mode are also projected onto the PC1-PC2 plane, the PC1-PC3 plane, and the PC2-PC3 plane, where the projected NLPCA is indicated by (overlapping) circles, the PCA is indicated by thin solid lines, and the projected data points are indicated by dots. One end of the NLPCA curve with maximum PC1 value is associated with the minimum value of the NLPC u and an extreme La Niña situation, while the opposite end of the curve corresponds to maximum u and extreme El Niño. The plotted PCs have been scaled up by a factor of 10. Figure 5. (a) NLPC1, u, and (b) the frequency distribution curve (FDC) for the NLPC1. The data have been high- passed prior to the NLPCA. 29 - 4 TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN is weak, particularly during the cool episodes of ENSO, in contrast to the interannual variability in the sea surface temperature (SST), where there are significant anomalies in the midlatitudes [Giese and Carton, 1999]. [14] Asymmetries between El Niño and La Niña spatial anomaly patterns, which are absent in the linear mode, are readily manifested in NLPCA1 (Figure 6). One notices much stronger anomalies occurring in the western tropical Pacific during extreme El Niño ( pattern E) than during extreme La Niña ( pattern A), even though in the eastern tropical Pacific, the anomalies are of similar magnitude. Furthermore, north of 30N, the anomalies are considerably stronger during El Niño than during La Niña (from comparing the amount of shaded area in pattern D with that in B, and between E and A). A useful way to characterize the asymmetry between El Niño and La Niña is by the spatial correlation coefficient. Between pattern A and E, the correlation is0.75, departing considerably from the correlation of 1 for the linear PCA mode. Another interesting nonlinear behavior is seen between typical El Niño ( pattern D) and extreme El Niño ( pattern E); as one proceeds from D to E, the cool anomalies in the western equatorial Pacific intensifies as expected, but the warm anomalies in the eastern equatorial Pacific weak- ens; that is, E is obtained from D by adding cool HCA in both the western and eastern equatorial Pacific. [15] We can compare our NLPCA results with the con- ventional composite method. Composites of HCA for 5 typical La Niña years (1971/1972, 1975/1976, 1984/1985, 1988/1989, 1995/1996) and 5 typical El Niño years (1972/ 1973, 1982/1983, 1986/1987, 1991/1992, 1997/1998) are shown in Figure 7, where the warm episodes have stronger heat content anomalies in the equatorial Pacific, especially in the western equatorial Pacific, than the cool episodes, in agreement with our NLPCA results. Of course, the averag- ing process in the composite method does not allow a distinction between typical and extreme El Niño conditions as in the NLPCA results. Also with the composite approach, one has to somewhat subjectively decide which ENSO episodes to include in the composite. [16] One reviewer cautioned that the data had unrealisti- cally small amplitudes in the southwestern tropical Pacific before the early 1980s [Lysne and Deser, 2002], compared to other data sources, and could affect our NLPCA calcu- lations. Fortunately, the extreme u values were attained after the earlier defective period, as seen in Figure 5a. We also recomputed the NLPCA excluding the earlier defective period, and the new extreme patterns A and E (not shown) are not very different from those in Figure 6. [17] Figure 8 is the Hovmöller diagrams showing the time evolution of the HCA along the equator from the NLPCA1, the linear PCA mode 1 and the leading six linear PCA modes. As in Figure 8a, the NLPCA1 rather well reflects observed features such as the eastward propagation of HCA, the oscillatory periods of 2–5 years, and the asymmetry of anomalies between El Niño and La Niña episodes. These features are absent or not obvious in the PCA mode 1 (Figure 8b), indicating that the NLPCA1 approximates the data set better than the PCA mode 1. 4. Decadal Mode [18] The first nonlinear decadal mode for the low-passed HCA data extracted from the NLPCA.cir network (Figure 3b) [Hsieh, 2001] is shown in the PC space (Figure 9). This Figure 6. Spatial anomaly patterns associated with the NLPC at A, B, D and E in Figure 5b. The contour interval is 0.4C, and areas with absolute values over 0.2C are shaded. TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN 29 - 5 mode explains 72% of the HCA variance, versus only 38% by the first PCA mode. The NLPC q in Figure 10a shows that the decadal variations are characterized by two jumps in q. The first jump, occurring in the early 1980s as detected also by linear PCA [Zhang et al., 1999], is closely asso- ciated with the large-scale climate regime shift in the Pacific Ocean around 1976. While the value at the time- point t in the low-passed data is actually averaged from the Figure 7. Composite of the HCA for several La Niña and El Niño years (see text), averaged over the extreme month of each episode. The contour interval is 0.4C, and areas with absolute values over 0.2C are shaded. Figure 8. Time-longitude plot of the reconstructed heat content anomalies along the equator. The reconstructed HCA is from (a) the first NLPCA mode, (b) the first PCA mode, and (c) the first six PCA modes. The contour interval is 0.6C, and areas with absolute values over 0.2C are shaded. 29 - 6 TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN original data over 61 months, thereby precluding fine temporal resolution, it nevertheless seems that the HCA (which involves subsurface temperature changes to 400 m depth) lags the sea surface condition changes around 1976, suggesting that it may take a few years for the surface regime shift to penetrate into the subsurface waters. The second jump in the early 1990s (Figure 10a) is mainly caused by q jumping from p to p, rather than by a physical regime shift like the first one. However, a clear contrast between the 1980s and the 1990s has been found in many observations such as sea level pressure, SST, low- level zonal wind, and subsurface ocean heat content anoma- lies in the Pacific [Kleeman et al., 1996; Latif et al., 1997; Ji et al., 1996]. [19] Decadal dependence of ENSO predictability has been found in many ENSO forecast models. While all models tended to have very good forecast skills in the 1980s, they suffered low skills in the 1990s, even with an improved initialization strategy [Chen et al., 1997]. It has been suggested that the decadal dependence of predictabil- ity may be due to the decadal changes in the mean state leading to the decadal variability of ENSO [e.g., Wang, 1995; Zhang et al., 1997]. Several possible mechanisms for changing the mean state have been suggested by some recent work, including the remote response in the tropical atmosphere to the midlatitude decadal oscillations, anthro- pogenic global warming, and the interaction between trop- ical and extratropical oceans by subduction processes [Kleeman and Power, 1999]. [20] The frequency distribution of the decadal mode (Figure 10b) presents a completely different shape than that of the interannual mode shown in Figure 5b. The FDC of the interannual mode is roughly Gaussian, whereas that of the decadal mode shows several spikes distributed over the full range of phase angles. As we lack sufficient samples to compute the FDC of the decadal mode, the relative short data record leads to the spiky frequency distribution. As such, the spatial patterns associated with these spikes may not be particularly meaningful. Instead, we examine the spatial patterns associated with four phases of the decadal mode, namely those corresponding to maximum p, max- imum q, minimum p, and minimum q (Figure 3b), with their locations in the PC space shown in Figure 9. [21] The spatial anomalies of the NLPCA1 mode corre- sponding to these four phases are shown in Figure 11, where Figures 11b and 11d are roughly the negative version of each other. Their basic pattern, similar to the linear PCA mode 1 (Figure 1b), is characterized by an anomaly in the midlati- tudes about 40N and one of the same sign in the western tropical Pacific, and by a weak anomaly of the opposite sign in the eastern Pacific. The anomaly in the midlatitudes appears to connect to the anomaly in the western tropical Pacific by a clockwise circulation. Hence this ‘‘subtropical gyre’’ pattern depicts a possible link of the decadal oscil- lation from the middle latitudes to the tropical Pacific. Such a pathway of decadal signals from midlatitudes to the tropics has also been proposed by other researchers through data Figure 9. The first NLPCA.cir mode for low-passed HCA data plotted as (overlapping) asterisks in the PC1-PC2-PC3 3-D space. The linear (PCA) mode is shown as a dashed line, and the data points are shown as dots. The circle denotes the point corresponding to min(q), the diamond corresponds to max( p), the pentagram corresponds to max(q), and the hexagram corresponds to min( p). The plotted PCs have been scaled up by a factor of 10. Figure 10. (a) NLPC1, q, and (b) NLPC1 FDC. The data have been smoothed by a 61-month running mean prior to performing NLPCA.cir. Note that q is periodically bounded within (p, p). TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN 29 - 7 analysis and modeling [e.g., Kleeman et al., 1999; Deser et al., 1999]. [22] In contrast to Figures 11b and 11d, the other pair of patterns in Figures 11a and 11c do not resemble each other strongly. The pattern in Figure 11c is characterized by an El Niño-like dipole structure along the equator, with positive anomalies in the east and negative anomalies in the west, suggesting that the pattern depicts the decadal variability of the ENSO mode. Our interpretation is that when the pattern in Figure 11c is on, warm phases of ENSO are reinforced, while cold phases are weakened. The prevalence of warm ENSO conditions in the period from 1991 to 1995 offers one example for this type of interaction between interannual and interdecadal variations. There are also notable midlati- tude anomalies in this decadal phase (Figure 11c). [23] The phase in Figure 11a reveals rather weak anoma- lies, though in the tropics, the anomalies are La Niña-like. The phase would enhance cool episodes and weaken warm episodes. But the fact that the phase in Figure 11a is much weaker than that in Figure 11c implies that the decadal variability for La Niña episodes is much less dramatic than for El Niño episodes. [24] This finding is consistent with the study by A. Wu and W. W. Hsieh (Nonlinear interdecadal changes of the El Nino-Southern Oscillation, submitted to Climate Dynamics, 2002) using nonlinear canonical correlation analysis (NLCCA) of wind stress and SST to examine the mid- 1970s climate regime shift. During 1981–1999, the location of the equatorial easterly anomalies during cool phases of ENSO was found to be unchanged from that observed in the 1961–1975 period, but during warm phases of ENSO, the westerly anomalies were shifted eastward by up to 25. From the position of the wind anomalies, the delayed oscillator theory would lengthen the duration of the warm episodes, but leave the cool episodes unchanged. Hence the NLCCA study also found much larger decadal changes in El Niño episodes than in La Niña episodes. [25] To further explore the spatial structure of the NLPCA1 in the time domain, we plot the Hovmöller diagrams for the reconstructed anomalies from the NLPCA1 along 40N and along 10S, the regions of the strongest decadal variability (Figures 12 and 13). For comparison, the reconstructed anomalies from the linear PCA mode 1 are also given. [26] As shown in Figures 12a and 13a, decadal changes can be clearly seen in the NLPCA1. Along 40N (Figures 12a and 12c), the Pacific basin exhibited a positive anomaly during the middle 1960s to 1981 with a magnitude of +0.6C–+1.0C around 1973–1974 centered in the Kur- oshio-extension region. The whole Pacific basin shifted to a large negative anomaly by 1981, which persisted about 10 years until 1990, when a new positive anomaly with a magnitude of +0.4C–+0.6C emerged (Figure 12c). This positive anomaly, which is not as wide as the earlier one in the 1960s to 1970s, has its center shifted 10–15 toward the east compared with the earlier one. Clearly the NLPCA1 (Figure 12a) models the regime shifts in Figure 12c much better than the PCA1 (Figure 12b), which missed the regime shift of the 1990s completely. [27] Along 10S (Figures 13a and 13c), from the mid- 1960s to the late 1970s, a strong positive anomaly in the western Pacific coincided with a weak negative anomaly in the eastern Pacific. Around early 1981, almost the whole Pacific along 10S shifted to a negative anomaly. This Figure 11. Spatial patterns corresponding to the four phases labeled in Figure 9 for the NLPCA mode 1. The contour interval is 0.1C, and areas with absolute values over 0.1C are shaded. 29 - 8 TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN negative anomaly persisted around 10 years in the eastern Pacific until about 1990, when the eastern Pacific shifted to a positive anomaly. In the western Pacific, the negative anomaly persisted until the late 1990s. 5. The 1980s and 1990s Decadal-Scale Changes [28] Over the last two decades, the upper ocean heat content experienced two prominent changes, resulting in generally warm conditions in the 1970s, cool conditions in the 1980s and mixed conditions in the 1990s, as seen in last section and in other works [e.g., Lysne and Deser, 2002]. The large-scale changes in the upper ocean thermal field around 1980 and 1990 can be seen as phase transitions of the decadal mode. Figures 14a and 14b show the differences in the average HCA between the 1970s and 1980s, and between the 1990s and 1980s, respectively. The spatial pattern in Figure 14a strongly resembles one of the phases of the decadal mode (Figure 11d), with a spatial correlation of 0.96, while the pattern in Figure 14b moderately resem- bles Figure 11c, with a correlation of 0.72. [29] There are several hypotheses to explain the mecha- nism of the decadal changes in the upper thermal field in the Pacific Ocean. The most popular one is the decadal changes in the wind stress curl affecting the gyre-scale patterns of the ocean circulation via the Sverdrup balance [Deser et al., 1999; Lysne and Deser, 2002]. The decadal signals in the wind stress curl are first forced into the surface ocean by Ekman pumping, and then transported to the thermocline by Rossby wave adjustment with the time of about 2–5 years [Deser et al., 1999]. [30] The occurrence of the decadal changes in SST (Figure 15) could be almost simultaneous to the changes in the wind around 1976 and 1988. That the decadal change in the HCA occurred 2–5 years after the wind change is probably due to the adjustment time scale of the subsurface ocean to surface changes. [31] As the 1980s decadal changes in the HCA lagged the surface changes by a longer time compared to the 1990s decadal change in the HCA, this suggests that the adjust- ment timescale of the subsurface to surface changes is considerably longer in the 1980 change than in the 1990 change. Possibly the physical processes involved in the two decadal changes were not completely the same. For exam- ple, for the 1990 decadal change, the main anomalies in the subsurface (Figure 14b) and the surface (Figure 15b) Figure 12. Time-longitude plot of the reconstructed heat content anomalies along 40N. The reconstructed HCA is based on (a) the first NLPCA mode, (b) the first PCA mode, and (c) the first six PCA modes (with 93% of the variance of the HCA). The contour interval is 0.2C, and areas with absolute values over 0.1C are shaded. TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN 29 - 9 occurred roughly in the same or neighboring regions, suggesting that the adjustment processes of the subsurface involved considerable vertical mixing and advection. But for the 1980 decadal change, the main anomalous change in the western equatorial subsurface ocean (Figure 14a) is very different from changes in the surface (Figure 15a or 15c), suggesting that the subsurface adjustment involved consid- erable horizontal transmission of the surface signal. Adjust- ment in the horizontal direction could involve the Rossby wave adjustment timescale, resulting in the longer response Figure 13. As for Figure 12, but along 10S. Figure 14. Differences in mean upper ocean heat content by (a) subtracting the mean of the 1970s from the mean of the 1980s and (b) subtracting the mean of the 1980s from the mean of 1990s. The contour interval is 0.2C. Shading denotes the regions where the two-tailed t-test for difference in means exceeds the 95% confidence level. 29 - 10 TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN time of the subsurface to surface in the 1980 change than in the 1990 change. In addition, a much slower process of subduction along the subtropical oceanic gyre may also be involved in the 1980s subsurface decadal change as sug- gested by Figure 14a. 6. Summary and Conclusion [32] We applied the nonlinear principal component anal- ysis technique to the observed upper ocean heat content anomalies in the Pacific basin from 1961 to 2000, and extracted the leading interannual and decadal modes. For the leading nonlinear interannual mode, the spatial anoma- lies are strongest in the equatorial Pacific, with an ENSO east-west seesaw pattern. As the nonlinear mode is not limited to a standing wave spatial anomaly pattern, it reveals considerable asymmetry between strong La Niña and strong El Niño. During strong El Niño, the negative anomaly in the equatorial western Pacific is much stronger than the positive anomaly found in this region during strong La Niña. This nonlinear interannual mode also manifests eastward phase propagation along the equator (Figure 8), in contrast to the standing wave found in the linear mode 1. [33] Four phases of the nonlinear decadal mode were examined. Two of them are roughly mirror images of each other, both showing a subtropical gyre pattern with the large anomaly in the midlatitudes circulating clockwise around the subtropical gyre towards the western tropical Pacific, a possible link from the middle latitudes to the tropical Pacific in the decadal mode. Two other phases of the decadal mode are related to decadal changes in the La Niña and El Niño characteristics. Since the one associated with La Niña has much weaker anomalies than the one associated with El Niño, it follows that the decadal changes in the character- istics of La Niña episodes are much weaker than the changes for El Niño episodes. [34] Over the last 2 decades, the nonlinear decadal mode experienced two phase shifts in 1981 and 1990, respectively, leading to the remarkable decadal changes in the upper ocean heat content in the 1980s and 1990s. From the equatorial to midlatitude Pacific, positive HCA during the mid-1960s to the late 1970s reversed to negative HCA around 1981. The regime shift around 1990 was also well represented by the nonlinear decadal mode; the negative anomalies in the midlatitudes and in the equatorial region in the 1980s reversed to positive anomalies around 1990 in the central midlatitude region and in the eastern equatorial Pacific. Prior to the two decadal changes in HCA, wind stress (curl) also changed in 1976 and 1988. While the SST changes were almost simultaneous with the wind changes, the HCA changes were delayed 2–5 years, corresponding to the Rossby wave adjustment timescale of the subsurface waters to surface changes. [35] The HCA change around 1980 was quite different from the one around 1990 in that the former occurred after the wind change with a much longer time delay than the Figure 15. Differences in mean SST by (a) subtracting the mean of 1967–1976 from the mean of 1977–1986, (b) subtracting the mean of 1979–1988 from the mean of 1989–1997, (c) subtracting the 1970s from the 1980s, and (d) subtracting the 1980s from the 1990s. The years of surface wind changes were around 1976 and 1988, while the HCA changes were around 1980 and 1990; hence Figures 15c and 15d are provided to temporally match Figure 14. The contour interval is 0.2C. Shading denotes the regions where the two-tailed t-test for difference in means exceeds 95% confidence. TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN 29 - 11 latter. The former (Figure 14a) showed the anomaly link from the midlatitudes to the western tropical Pacific via the subtropical gyre, while the latter (Figure 14b) did not. The former was also more different from the corresponding SST anomalies (Figure 15a) than the latter was from SST (Figure 15b), suggesting that the signals involved more horizontal transmission in the former than in the latter, where the surface signals appeared to be transmitted more vertically to the subsurface. The leading linear PCA mode was able to detect the former change but not the latter, which was clearly detected by the leading nonlinear PCA mode. [36] Acknowledgments. The heat content data were kindly supplied by the Joint Environmental Data Analysis Center, headed by Warren White, at the Scripps Institution of Oceanography. This work was supported by research and strategic grants to W. Hsieh from the Natural Sciences and Engineering Research Council of Canada. References Chen, D., S. E. Zebiak, and M. A. Cane, Initialization and predictability of a coupled ENSO forecast model, Mon. Weather Rev., 125, 773–788, 1997. Cybenko, G., Approximation by superpositions of a sigmoidal function, Math. Control Signals Syst., 2, 303–314, 1989. Deser, C., M. A. Alexander, and M. S. Timlin, Evidence for a wind-driven intensification of the Kuroshio Current Extension from the 1970s to the 1980s, J. Clim., 12, 1697–1706, 1999. Giese, B., and J. Carton, Interannual and decadal variability in the tropical and midlatitude Pacific ocean, J. Clim., 12, 3402–3418, 1999. Hamilton, K., and W. W. Hsieh, Representation of the quasi-biennial oscil- lation in the tropical stratospheric wind by nonlinear principal component analysis, J. Geophys. Res., 107, 4232, doi:10.1029/2001JD001250, 2002. Hsieh, W. W., Nonlinear principal component analysis by neural networks, Tellus, Ser. A, 53, 599–615, 2001. Ji, M., A. Leetmaa, and V. E. Kousky, Coupled model prediction of ENSO during the 1980s and 1990s at the National Centers for Environmental Prediction, J. Clim., 9, 3105–3120, 1996. Kirby, M. J., and R. Miranda, Circular nodes in neural networks, Neural Comp., 8, 390–402, 1996. Kleeman, R., and S. B. Power, Modulation of ENSO variability on decadal and longer timescales, in El Nino and the Southern Oscillation, Multi- scale Variability and its Impacts on Natural Ecosystems and Society, edited by H. F. Diaz and V. Markgraf, pp. 413–441, Cambridge Univ. Press, New York, 1999. Kleeman, R., R. Colman, N. R. Smith, and S. B. Power, A recent change in the mean state of the Pacific Ocean: Observational evidence, atmospheric and oceanic responses, J. Geophys. Res., 101, 20,483–20,499, 1996. Kleeman, R., J. P. McCreary, and B. A. Klinger, A mechanism for the decadal variation of ENSO, Geophys. Res. Lett., 26, 1743, 1999. Kramer, M. A., Nonlinear principal component analysis using autoassocia- tive neural networks, AIChE J., 37, 233–243, 1991. Latif, M., R. Kleeman, and C. Eckert, Greenhouse warming, decadal varia- bility, or El Niño? An attempt to understand the anomalous 1990s, J. Clim., 10, 2221–2239, 1997. Lysne, J., and C. Deser, Wind-Driven thermocline variability in the Pacific: A model/data comparison, J. Clim., 15, 829–845, 2002. Miller, A. J., W. B. White, and D. R. Cayan, North Pacific thermocline variations on ENSO timescales, J. Phys. Oceanogr., 27, 2023–2039, 1997. Miller, A. J., D. R. Cayan, and W. B. White, Awestward-intensified decadal change in the North Pacific thermocline and gyre-scale circulation, J. Clim., 11, 3112–3127, 1998. Monahan, A. H., Nonlinear principal component analysis by neural net- works: Theory and application to the Lorenz system, J. Clim., 13, 821– 835, 2000. Monahan, A. H., Nonlinear principal component analysis, tropical Indo- Pacific sea surface temperature and sea level pressure, J. Clim., 14, 219– 233, 2001. Monahan, A. H., L. Pandolfo, and J. C. Fyfe, The preferred structure of variability of the Northern Hemisphere atmospheric circulation, Geophys. Res. Lett., 28, 1019–1022, 2001. Schneider, N., A. J. Miller, M. A. Alexander, and C. Deser, Subduction of decadal North Pacific temperature anomalies, observations and dynamics, J. Phys. Oceanogr., 29, 1056–1070, 1999. Tang, Y., and W. W. Hsieh, ENSO simulation and prediction in a hybrid coupled model with data assimilation, J. Met. Soc. Jpn., in press, 2003. Trenberth, K., and J. M. Hurrell, Decadal atmosphere-ocean variations in the Pacific, Clim. Dyn., 9, 303–319, 1994. Wallace, J., E. Rasmusson, T. Mitchell, V. Kousky, E. Sarachik, and H. von Storch, On the structure and evolution of ENSO-related climate variabil- ity in the tropical Pacific: Lessons from TOGA, J. Geophys. Res., 103, 14,241–14,259, 1998. Wang, B., Interdecadal changes in El Niño onset in the last four decades, J. Clim., 8, 267–285, 1995. White, W. B., Design of a global observing system for gyre-scale upper ocean temperature variability, Prog. Oceanogr., 36, 169–217, 1995. Zhang, R.-H., L. M. Rothstein, and A. J. Busalacchi, Interannual and dec- adal variability of the subsurface thermal structure in the Pacific Ocean: 1961–90, Clim. Dyn., 15, 703–717, 1999. Zhang, Y., J. M. Wallace, and D. S. Battisti, ENSO-like interdecadal varia- bility, J. Clim., 10, 1004–1020, 1997.  W. W. Hsieh, Department of Earth and Ocean Sciences, University of British Columbia, Vancouver, B.C., V6T 1Z4, Canada. (whsieh@eos.ubc. ca) Y. Tang, Courant Institute of Mathematical Sciences, New York University, 251 Mercer Street, New York, NY 10012, USA. ( ytang@ cims.nyu.edu) 29 - 12 TANG AND HSIEH: NONLINEAR MODES IN SUBSURFACE PACIFIC OCEAN"@en ; edm:hasType "Article"@en ; edm:isShownAt "10.14288/1.0041795"@en ; dcterms:language "eng"@en ; ns0:peerReviewStatus "Reviewed"@en ; edm:provider "Vancouver : University of British Columbia Library"@en ; dcterms:publisher "American Geophysical Union"@en ; ns0:publisherDOI "10.1029/2001JC001236"@en ; dcterms:rights "Attribution-NonCommercial-NoDerivatives 4.0 International"@en ; ns0:rightsURI "http://creativecommons.org/licenses/by-nc-nd/4.0/"@en ; ns0:scholarLevel "Faculty"@en ; dcterms:title "Nonlinear modes of decadal and interannual variability of the subsurface thermal structure in the Pacific Ocean"@en ; dcterms:type "Text"@en ; ns0:identifierURI "http://hdl.handle.net/2429/32827"@en .