Abstract

The Pacific decadal oscillation (PDO), defined as the leading empirical orthogonal function of North Pacific sea surface temperature anomalies, is a widely used index for decadal variability. It is shown that the PDO can be recovered from a reconstruction of North Pacific sea surface temperature anomalies based on a first-order autoregressive model and forcing by variability of the Aleutian low, El Niño–Southern Oscillation (ENSO), and oceanic zonal advection anomalies in the Kuroshio–Oyashio Extension. The latter results from oceanic Rossby waves that are forced by North Pacific Ekman pumping. The SST response patterns to these processes are not orthogonal, and they determine the spatial characteristics of the PDO. The importance of the different forcing processes is frequency dependent. At interannual time scales, forcing from ENSO and the Aleutian low determines the response in equal parts. At decadal time scales, zonal advection in the Kuroshio–Oyashio Extension, ENSO, and anomalies of the Aleutian low each account for similar amounts of the PDO variance. These results support the hypothesis that the PDO is not a dynamical mode, but arises from the superposition of sea surface temperature fluctuations with different dynamical origins.

1. Introduction

The Pacific decadal oscillation (PDO) index (Mantua et al. 1997) is widely used to characterize North Pacific decadal variability and anomalies of Northern Hemisphere climate and the North Pacific ecosystem. Yet, the processes underlying the temporal and spatial characteristics of the PDO have not been clarified. Here, we present a statistical analysis of North Pacific sea surface temperature (SST), sea level pressure (SLP), and wind stress to identify the drivers of the PDO and their relative contributions to its temporal and spatial patterns.

The PDO is defined1 as the leading empirical orthogonal function (EOF) of monthly anomalies (deviations from the climatological annual cycle) of sea surface temperature in the Pacific poleward of 20°N (Davis 1976; Mantua et al. 1997). Its temporal evolution (Fig. 1) is marked by variability on interannual and decadal time scales, with several pronounced shifts—the 1976/77 transition being the most famous and extensively studied (Trenberth 1990; Graham 1994; Trenberth and Hurrell 1994; Miller et al. 1994). The spatial pattern of the PDO (Fig. 2) resembles a horseshoe, with temperature anomalies in the central North Pacific surrounded by anomalies of opposite sign in the Alaska gyre, off California, and toward the Tropics.

Fig. 1.

Time series of monthly values of (from the top) the PDO, Niño-3.4, NPI, PDEL, and PAVG. All time series are normalized to unit standard deviation. The thick line overlaid on the top curve shows the PDO based on Jul–Jun annual averages of SSTA. The correlation of annual averages of the monthly values of the PDO, and of the PDO estimated from annual averages of SSTA, is 0.99, and the latter is used in this study. The dotted lines overlaid on the lower two curves denote PDEL and PAVG obtained from North Pacific wind stress curl that had the regression with low-pass time series of Niño-3.4 and NPI removed.

Fig. 1.

Time series of monthly values of (from the top) the PDO, Niño-3.4, NPI, PDEL, and PAVG. All time series are normalized to unit standard deviation. The thick line overlaid on the top curve shows the PDO based on Jul–Jun annual averages of SSTA. The correlation of annual averages of the monthly values of the PDO, and of the PDO estimated from annual averages of SSTA, is 0.99, and the latter is used in this study. The dotted lines overlaid on the lower two curves denote PDEL and PAVG obtained from North Pacific wind stress curl that had the regression with low-pass time series of Niño-3.4 and NPI removed.

Fig. 2.

(top) PDO pattern from observations, (middle) derived as the leading EOF of the reconstructed SSTA, and (bottom) derived from the reconstruction obtained with a spatially constant value of α̃. The magnitude depicts the SSTA (K) for a unit deviation of the principal component; the contour level is 0.2, with negative contours dotted.

Fig. 2.

(top) PDO pattern from observations, (middle) derived as the leading EOF of the reconstructed SSTA, and (bottom) derived from the reconstruction obtained with a spatially constant value of α̃. The magnitude depicts the SSTA (K) for a unit deviation of the principal component; the contour level is 0.2, with negative contours dotted.

The PDO and its associated sea level pressure pattern have attracted considerable attention because their time evolutions are correlated with anomalies of the Pacific ecosystem (Mantua et al. 1997); North American precipitation, streamflow, and surface temperature anomalies (Mantua and Hare 2002; Dettinger et al. 1998; Cayan et al. 1998); surface temperature anomalies in northeastern Asia (Minobe 2000); fluctuations of the Asian monsoon (Krishnan and Sugi 2003); and a modulation of El Niño–Southern Oscillation (ENSO) teleconnections (Gershunov and Barnett 1998). This covariation can result from the PDO forcing the atmosphere or from a common forcing of the PDO and covarying climate anomalies. Model experiments suggest the latter explanation for the decadal modulation of ENSO teleconnections (Pierce 2002). However, this question has not been clarified from observations, and requires the identification of the forcing mechanism of the PDO.

The PDO pattern is qualitatively consistent with the atmospheric forcing associated with fluctuations of position and strength of the Aleutian low. In the central North Pacific, a deepened Aleutian low decreases SST by advection of cool and dry air from the north, by increases of westerly winds and ocean-to-atmosphere turbulent heat fluxes, and by strengthened equatorward advection of temperature by Ekman currents. In the eastern regions, a deepened Aleutian low enhances poleward winds and leads to warm anomalies of surface temperature.

However, this qualitative picture hides the complex processes governing changes of the Aleutian low and North Pacific SST anomalies. The North Pacific atmosphere affects the surface heat budget via anomalies of surface fluxes and Ekman advection (Alexander 1992), and the excitation and propagation of Rossby waves that affect the upper ocean in the Kuroshio–Oyashio Extension (KOE) region (Miller et al. 1998; Deser et al. 1999; Xie et al. 2000; Miller and Schneider 2000; Seager et al. 2001; Schneider et al. 2002; Qiu 2003). The anomalies of the atmosphere reflect its intrinsic variability (Davis 1976; Pierce et al. 2001; Barnett et al. 1999b; An and Wang 2005), interannual variability associated with ENSO (Newman et al. 2003) that is teleconnected to the North Pacific (Alexander et al. 2002, 2004), and remote forcing from the western Pacific and eastern Indian Ocean (Deser et al. 2004) and possibly elsewhere. Tropical ENSO variability also influences North Pacific SST anomalies via poleward-propagating coastal Kelvin waves along the Pacific coast of the Americas (Chelton and Davis 1982; Clarke and Lebedev 1999).

Here, we quantify the forcing mechanism of North Pacific SST and the PDO, and investigate the hypothesis that the PDO time evolution and spatial pattern result from a superposition of these different forcing mechanisms. To this end, we reconstruct North Pacific SST anomalies (SSTA) using a first-order autoregressive (AR-1) model, forced indices that track ENSO, sea level pressure variations of the Aleutian low, and oceanic changes in the Kuroshio–Oyashio Extension. The leading EOF of this reconstruction is nearly identical to the observed PDO, and allows for an investigation of its properties.

In the following, we introduce the data and forcing indices, and present the method, skill, and interpretation of the SSTA reconstruction. This reconstruction is then used to determine and discuss the PDO: its kinematics; evolution, including the relative importance of tropical versus North Pacific forcing spectrum; and predictability.

2. Data

Throughout this study, we use monthly anomalies from the average seasonal cycle of SST, sea level pressure, and wind stress from the National Centers for Environmental Prediction–National Center for Atmospheric Research (NCEP–NCAR) reanalysis of the period from 1950 to 2003 (Kalnay et al. 1996). The PDO used here is based on the leading EOF of July–June annual averages of SSTA, rather than on monthly anomalies (Mantua et al. 1997). The July–June annual average focuses on the winter season when both teleconnected forcing from the Tropics and intrinsic atmospheric variability are most pronounced (Newman et al. 2003). Annual averages of the leading principal component based on monthly anomalies are nearly identical with the PDO used here, with correlation coefficients of 0.994 (Fig. 1, top). To extend the study back in time, we also employ the 1900–97 time series of the observed PDO time index, the North Pacific index of sea level pressure in the Aleutian low (NPIT; Trenberth and Hurrell 1994), the cold tongue index (CTI) of equatorial Pacific sea surface temperature anomalies (all of which are available online at http://tao.atmos.washington.edu/data_sets), and the optimal tropical index (OTI; Deser et al. 2004; information available online at http://www.cgd.ucar.edu/~cdeser/Data/).

3. Forcing indices

North Pacific SST anomalies are associated with ENSO, with atmospheric variations that are intrinsic to the North Pacific, and, in the Kuroshio–Oyashio region, with oceanic Rossby waves that are excited by North Pacific Ekman pumping. The time evolution of each of these processes will be represented by an index time series, to be identified in this section by a combination of physical reasoning and statistical optimization.

a. ENSO and North Pacific index

The North Pacific SST is composed of at least two temporally independent anomaly patterns: one associated with ENSO, and one resulting from variations of atmospheric pressure intrinsic to the North Pacific (Davis 1976; Alexander 1992; Nakamura et al. 1997; Zhang et al. 1997; Nakamura and Yamagata 1999; Tomita et al. 2001; Barlow et al. 2001; Nakamura and Kazmin 2003; An and Wang 2005). The impact of ENSO, represented by the Niño-3.4 temperature anomaly in the eastern equatorial Pacific (5°S–5°N, 170°–120°W) has been documented by assuming that annual-averaged values of the PDO amplitude are governed by an AR-1 process, forced by the Niño-3.4 time series (Newman et al. 2003). Correlations of North Pacific anomalies of sea level pressure and SST and with changes in the warm pool region and South Pacific convergence zone, indicate that processes in these areas also affect the North Pacific (Deser et al. 2004; Garreaud and Battisti 1999). In addition, intrinsic variability in the extratropical atmospheric circulation, particularly changes in the position and strength of the Aleutian low (Trenberth and Hurrell 1994) and associated storm-track activity, forces North Pacific SST anomalies (e.g., Pierce et al. 2001; Pierce 2001).

To identify forcing indices of the PDO, we search for regions where sea level pressure anomalies have skill in reconstructing the observed time series of the PDO amplitude. We extend the approach of Newman et al. (2003) and assume that the July–June averages of the PDO amplitude D are governed by an AR-1 process, forced by the concurrent July–June averages of anomalies of sea level pressure (FSLP),

 
formula

where t is time in years. The damping coefficients α′ and regression coefficients γ′ are determined by a least squares fit (see appendix A) of Eq. (1) to the observed values D and FSLP, one location at a time. Using this fit, the entire history of the observed PDO amplitude D is reconstructed from the observed initial condition and SLP anomaly time series.

The correlation between the observed and reconstructed PDO amplitude, plotted at the location of the forcing time series (Fig. 3, upper panel), shows high values in the North Pacific (>0.8), eastern South Pacific, western Pacific warm pool, eastern Indian Ocean, and tropical Atlantic; and indicates that SLP at these locations has skill in the reconstruction of the PDO. The tropical signals are associated with ENSO. Removing the regression with Niño-3.4 from the SLP anomaly time series and repeating the fitting exercise removes all of the centers of high PDO amplitude reconstruction skill in the Tropics, with only the North Pacific center remaining, albeit with reduced correlation (Fig. 3, lower panel). This suggests that sea level pressure in the central North Pacific and tropical forcing, as described by Niño-3.4, account for the majority of the variations of the PDO time amplitude.

Fig. 3.

Correlation times 10 of observed and reconstructed PDO time evolution, with a contour interval of 1. At every location, the PDO is fit to and reconstructed from an AR-1 process of Eq. (1), forced by the local sea level pressure time series. High correlation indicates regions in which sea level pressure anomalies are most influential to the observed PDO. Results for (top) full sea level pressure and (bottom) sea level pressure with the Niño-3.4 temperature index removed are shown. The removal is achieved by subtracting the (least squares fit) linear regression of anomalous sea level pressure and Niño-3.4. Dithering in the North Pacific marks the area of spatial average used to estimate the NPI.

Fig. 3.

Correlation times 10 of observed and reconstructed PDO time evolution, with a contour interval of 1. At every location, the PDO is fit to and reconstructed from an AR-1 process of Eq. (1), forced by the local sea level pressure time series. High correlation indicates regions in which sea level pressure anomalies are most influential to the observed PDO. Results for (top) full sea level pressure and (bottom) sea level pressure with the Niño-3.4 temperature index removed are shown. The removal is achieved by subtracting the (least squares fit) linear regression of anomalous sea level pressure and Niño-3.4. Dithering in the North Pacific marks the area of spatial average used to estimate the NPI.

The average of sea level pressure in the area where reconstruction skill without Niño-3.4 exceeds 0.55 (dithered region in Fig. 3, lower panel) serves as an index for the North Pacific atmospheric forcing, and will be referred to as North Pacific index (NPI), because it is approximately collocated (and highly correlated, at 0.85) with the index by the same name introduced by Trenberth and Hurrell (1994), here called NPIT. These indices describe anomalies of the Aleutian low, the primary mode of variability of annual-averaged anomalies of North Pacific sea level pressure. The leading EOF of July–June averages of surface pressure explains more than half of the sea level pressure variance, and has a correlation with the Trenberth and Hurrell index in excess of 0.9. In contrast to the primarily interannual variations of Niño-3.4, the monthly values of the NPI show variability on all time scales (Fig. 1).

b. Ocean circulation anomalies in the Kuroshio–Oyashio Extension

In the boundary currents of the northwestern Pacific, fluctuations of the wind-driven circulation affect sea surface temperature (Xie et al. 2000). These circulation anomalies result from the accumulation of North Pacific Ekman pumping FEk along Rossby wave trajectories

 
formula

where p is the amplitude of a first baroclinic mode Rossby wave and represents an index for the depth of the thermocline, and c is the Rossby wave speed (Sturges and Hong 1995; Schneider and Miller 2001). The anomalies of p affect the position and sharpness of the Kuroshio–Oyashio Extension and cause anomalies of surface temperatures, either directly through changes of thermocline depth and meridional shifts of the Kuroshio–Oyashio frontal system (Seager et al. 2001; Schneider et al. 2002), or by anomalies of zonal advection (Qiu 2003). The former process is most effective in the winter season, when deep mixed layers couple thermocline depth and SST (Xie et al. 2000; Schneider et al. 2002), and is simulated by assuming that anomalies of SST and thermocline depth p are proportional. The latter process, anomalous zonal advection, invokes the meridional gradient of pressure, and is affected by the meridional dependence of the forcing and of the Rossby wave speed c (Qiu 2003). Anomalous advection acts year round, and should, therefore, be more important for annual averages.

To derive indices for the effect of this gyre adjustment on SST, anomalies of p are estimated from Eq. (2) for the North Pacific, with the Ekman pumping calculated from anomalies of the NCEP–NCAR reanalysis wind stress (Schneider and Miller 2001). The Rossby wave speeds are based on satellite observations of sea level (Qiu 2003). The estimates p are averaged in the Kuroshio–Oyashio Extension between 140° and 170°E, and a thermocline depth index PAVG is then obtained by the average between 35° and 40°N; the zonal advection index PDEL is evaluated by differencing the zonal-averaged p between 38° and 40°N. The latitudes correspond to the location of the Kuroshio–Oyashio Extension region and were refined by trial and error, a selection procedure that is considered in the estimation of the significance (see below).

The Kuroshio–Oyashio Extension thermocline depth index PAVG and zonal advection index PDEL are dominated by variance at decadal frequencies (Fig. 1), as expected from the integration of atmospheric forcing along Rossby wave trajectories (Frankignoul et al. 1997), and can only account for a significant fraction of SSTA at long time scales. The KOE indices are largely independent of wind stress curl anomalies that are associated with Niño-3.4 and NPI. Removal of the regression of the wind stress curl with these indices before the solution of Eq. (2) has very little impact on PDEL and merely removes a trend from PAVG (dotted lines in Fig. 1).

4. Reconstruction of SST

The PDO is defined as the leading mode of SSTA, so we begin by reconstructing SSTA, and will later return to examine how well the leading EOF of the reconstructed SSTA compares to the actual PDO.

a. Approach

The budget for surface temperature anomalies T is assumed to be governed by an AR-1 model (Hasselmann 1976)

 
formula

where t is time; α̃ is the damping rate; and the different forcing terms, denoted by index i, have the time evolution i with a spatial footprint γi; and a repeated index i implies the Einstein summation convention. The first term on the right-hand side represents the damping of SST anomalies by the air–sea heat fluxes, with a decay rate α̃ determined by the feedback of the air–sea heat flux and by the mixed-layer depth. While these processes and thus α̃ have a seasonal cycle, we focus on the time scales that are longer than a year, so that α̃ is taken as a constant representative value for the entire year. Considering heat content from the surface to winter mixed-layer depths would remove the complications that are associated with the seasonal cycle of the mixed layer and the reemergence of temperature anomalies (Deser et al. 2003). However, the PDO is defined with respect to SST, and we need to determine the relationship between seasonal variations resulting from reemergence and the representative annual value α̃ of the feedback. In appendix B, an extension of the Hasselmann (1976) model to include reemergence is used to show that the α̃ derived from annual averages of SST lies between summer and winter damping rates, but close to winter values.

Equation (3) is discretized to describe the year-to-year evolution of July–June averages of temperature anomalies. With α̃ constant in time, Eq. (3) is integrated over 1 yr, and then averaged from July to the subsequent June. This yields for the annual-averaged value of temperature anomalies at every location x, y

 
formula
 
formula

where δ is the 1-yr time step of the discretization, and the overline represents July–June 12-month averages. The annual-averaged value of SSTA results from its value of the previous year, discounted at a rate of α̃, and the annual averages of the forcing fields, which have been filtered by a one-sided exponential with decay α̃ resulting from the damping of temperature anomalies. For ease of reference, Fi is defined in Eq. (5) as the annual average of the filtered forcing time series with filter coefficients varying in x, y, and α = exp(−α̃δ) is the discretized autoregressive parameter.

The discretized temperature in Eq. (4) preserves the lead–lag relationships of Eq. (3) between the forcing and temperature, unlike the ad hoc approach of Eq. (1), which relates the value of the PDO to contemporaneous forcing fields. A slight complication arises from the spatial dependence of α̃ that renders the filtered forcing function Fi as a function of time and space, and produces nonseparable forcing for the evolution of annually averaged temperature anomalies. However, we will find that α̃ has only a weak spatial dependence so that the spatial variance of forcing functions Fi are small and not important.

The reconstruction of North Pacific SST anomalies reduces to an exercise of fitting the observed evolution of to the right-hand side of Eq. (4) by adjusting the damping α̃(x, y) and spatial footprints of the forcing γi(x, y) at every location. Equation (4) is nonlinear in α̃ as a result of the filtering of the forcing time series, and the coefficients are determined iteratively. Starting from a first guess of α̃, the time integrals of i and June–July averages are formed to yield Fi. We then determine α and γi from Eq. (5) by multilinear regression using a least squares approach (see appendix A), and repeat the procedure with the new estimate of α̃ obtained from α. Typically, this procedure converges in less than 10 iterations.

The skill of the fit is determined by comparison of the observed and fitted temperatures, and is obtained in two related ways: as a 1-yr “forecast” of year n from the observed previous year’s SSTA and the forcing time series2 (the procedure followed by Newman et al. 2003 when comparing the observed and fitted time evolution of the PDO), or from a 50-yr “hindcast,” initializing Eq. (5) with the observed 1950/51 July average and determining the evolution until 2002/03 using only the forcing time series. The two approaches yield essentially the same results, albeit with the hindcast—the more demanding test—exhibiting overall lower skill values than that of the forecast, and we will present the hindcast skills only.

The significance of the reconstruction and of each regression pattern is estimated by comparison to the performance of the null hypothesis: that a selection of random red-noise series can do as well as the physically meaningful series chosen here. Figure 4 shows the 95% highest correlation level from 1000 reconstructions of the annual-averaged values of the PDO from the red-noise time series that have the same 1-yr lag autocorrelation as the annual averages of the index series Fi. To account for artificial skill, we repeat with the noise time series the “prospecting” used in choosing the NPI, PAVG, and PDEL by assuming that the NPI was selected out of North Pacific sea level pressure from 20 independent regions or modes, and that the PDEL and PAVG were selected each from a total of five possible latitude bands. These estimates of the degrees of freedom are deliberately high, so as to yield a conservative estimate of significance.

Fig. 4.

(top) The significance levels (95%) of the correlation of 50-yr hindcast of observed SSTA and best-fit AR-1 model forced by red-noise time series. (bottom four panels) The 95% significance levels for the contributions of Niño-3.4, NPI, PDEL, and PAVG, respectively. The significance levels are based on 1000 reconstructions of SSTA, and include the effect of artificial skill from the selection used in deriving some of the forcing indices. Contour interval is 0.1.

Fig. 4.

(top) The significance levels (95%) of the correlation of 50-yr hindcast of observed SSTA and best-fit AR-1 model forced by red-noise time series. (bottom four panels) The 95% significance levels for the contributions of Niño-3.4, NPI, PDEL, and PAVG, respectively. The significance levels are based on 1000 reconstructions of SSTA, and include the effect of artificial skill from the selection used in deriving some of the forcing indices. Contour interval is 0.1.

Specifically, each reconstruction of the PDO from the noise time series is obtained just as in the analysis of the observations. The observed PDO is fit to an AR-1 process that is forced by one “NINO3.4” noise realization. The residual is then reconstructed using 20 independent “NPI” noise time series, and the realization with the best fit is chosen. We then choose the best of five noise realizations each for PDEL and PAVG that capture the SSTA variations of a representative time series in the western North Pacific. For this set of “best” noise forcings, we obtain the PDO reconstruction skill, and then determine at every point in the North Pacific the combined skill of the noise series in the reconstruction of observed SSTA and the individual skill of each of the forcing indices.

For the PDO reconstruction, this yields 95% (90%) significance levels of the hindcast correlation of 0.68 (0.66). The 95% significance levels of the SSTA reconstruction skill (Fig. 4) reach values of 0.6–0.65 in the central and eastern North Pacific, and mimic the horseshoe pattern of the PDO. The 95% significance levels of the individual forcing indices vary by about 0.4 for Niño-3.4, and between 0.45 and 0.55 for NPI, PDEL, and PAVG.

b. Skill

The reconstruction of SSTA based on the forcing indices of ENSO, NPI, PDEL, and PAVG is very successful in the North Pacific. Correlations of the hindcast and observations of SSTA (Fig. 5) are larger than 0.75 in the areas of large loading of the PDO in the central and eastern North Pacific (Fig. 2, top) and mimic closely the variance of SSTA that is explained by the PDO itself. The skill in these regions is significant because it exceeds the 95% significance levels that are derived above (Fig. 4). The reconstruction passes the significance test except south of 25°N, and along the zero line of the PDO pattern.

Fig. 5.

(top) Correlation times 10 of 50-yr AR-1 model SSTA hindcast and observed Jul–Jun averages of SSTA. (bottom four panels) Correlation of observed and hindcast SSTA using only forcing by Niño-3.4, NPI, PDEL, and PAVG, respectively. Contour interval is 2.

Fig. 5.

(top) Correlation times 10 of 50-yr AR-1 model SSTA hindcast and observed Jul–Jun averages of SSTA. (bottom four panels) Correlation of observed and hindcast SSTA using only forcing by Niño-3.4, NPI, PDEL, and PAVG, respectively. Contour interval is 2.

Central to our purpose here is the comparison of the leading EOF of observed SSTA and its reconstruction from the hindcast. The spatial patterns of original and reconstructed fields are nearly identical (Fig. 2), with a pattern correlation in excess of 0.9 (Table 1). Small differences occur in the amplitudes of the spatial pattern in the central North Pacific, where the ridge of high loading extends further east in the reconstruction. The principal component also shows a good fit (Fig. 6), with all transitions captured. The correlation of observed and reconstructed PDO is 0.93, which is highly significant (Table 1). The only differences occur in slightly reduced amplitudes in the late 1950s and 1990. This suggests that the processes underlying the PDO are captured, and the reconstruction can be used to explore the properties of the PDO.

Table 1.

Correlation coefficients of observed PDO (top) patterns and (bottom) principal component with leading EOF of the SSTA reconstruction using spatially varying α (column 1) or constant α0 (column 2). Column 3 shows the results obtained from setting the off-diagonal elements of the covariance matrix of the forcing indices to zero, and column 4 shows results from considering only the forcing by Niño-3.4* and NPI*.

Correlation coefficients of observed PDO (top) patterns and (bottom) principal component with leading EOF of the SSTA reconstruction using spatially varying α (column 1) or constant α0 (column 2). Column 3 shows the results obtained from setting the off-diagonal elements of the covariance matrix of the forcing indices to zero, and column 4 shows results from considering only the forcing by Niño-3.4* and NPI*.
Correlation coefficients of observed PDO (top) patterns and (bottom) principal component with leading EOF of the SSTA reconstruction using spatially varying α (column 1) or constant α0 (column 2). Column 3 shows the results obtained from setting the off-diagonal elements of the covariance matrix of the forcing indices to zero, and column 4 shows results from considering only the forcing by Niño-3.4* and NPI*.
Fig. 6.

(top) Time series of Jul–Jun annual averages of PDO (black line). Overlaid in gray is the PDO based on SSTA reconstruction. (bottom four) The Fi(t) for Niño-3.4, NPI, PDEL, and PAVG, respectively, are shown. Thin solid lines denote the spatial averages of the filtered index time series, and tiny vertical bars their spatial standard deviation. Dashed lines show the rotated time series Fi(t)*, and the thick gray lines denote their projection onto the observed pattern of the PDO [see Eq. (7)]. All time series have been normalized to unit standard deviation.

Fig. 6.

(top) Time series of Jul–Jun annual averages of PDO (black line). Overlaid in gray is the PDO based on SSTA reconstruction. (bottom four) The Fi(t) for Niño-3.4, NPI, PDEL, and PAVG, respectively, are shown. Thin solid lines denote the spatial averages of the filtered index time series, and tiny vertical bars their spatial standard deviation. Dashed lines show the rotated time series Fi(t)*, and the thick gray lines denote their projection onto the observed pattern of the PDO [see Eq. (7)]. All time series have been normalized to unit standard deviation.

The hindcast skill of the individual components (Fig. 5) shows that Niño-3.4, NPI, and PDEL account for significant fractions of North Pacific SSTA variability. The skill of Niño-3.4 is largest in the eastern portion of the basin and along the coast of North America, while NPI accounts for the largest variations along 40°N and in the central and the western portion of the Pacific. The zonal advection index has skill in the Kuroshio–Oyashio Extension, as expected. The contribution of PAVG, while revealing a physically plausible pattern (as will be discussed below), nowhere exceeds the skill that is expected from an arbitrary red-noise process, and PAVG will therefore not be considered in the reconstruction of the PDO.

c. Damping time scale

The lag-1 autocorrelation of annual averages of SSTA α is largely a constant in space with typical values of 0.2–0.3, and corresponds to a damping time scale α̃−1 of 7.5–10 months (Fig. 7). Off the coasts of Japan and China the damping times are larger, but in these areas skill is small, and not significant (Figs. 4, 5). The damping time reflects the feedback ε to surface temperature of the air–sea heat flux, normalized by an effective mixed-layer depth H, and the density ρ0 and heat capacity cp of seawater. With ε = 15 W m−2 K−1 (Deser et al. 2003), ρ0 = 1000 kg m−3 and cp = 4 · 103 J kg−1 K−1, the effective annually averaged mixed-layer depth is H = ε(ρ0cpα̃)−1, and varies between 74 and 98 m. This is shallower than the deep wintertime mixed-layer depths that are required to explain the autocorrelation of the upper-ocean heat content (Deser et al. 2003), but is consistent with our focus on annually averaged sea surface temperature anomalies and typical mixed-layer depths of the North Pacific (see appendix B).

Fig. 7.

Estimate of damping time α̃−1 of Eq. (3) in months. Contour interval is 10 months.

Fig. 7.

Estimate of damping time α̃−1 of Eq. (3) in months. Contour interval is 10 months.

d. Forcing time series

Given the damping time scales, the filtered forcing time series Fi can be obtained and are shown in Fig. 6, normalized to unit standard deviation. The filtered Niño-3.4 series shows the major El Niño and La Niña events, slightly delayed, as expected, from the integration. The 1976/77 change is marked by several years of depressed Niño-3.4, followed after the shift by a recovery to near-normal forcing. The NPI time series displays more year-to-year variability and also shows a trend to more negative values, punctuated by dramatic variations of either sign. The transition to higher values after 1986 corresponds to an intensified Pacific storm track that is reported by Nakamura et al. (2002, their Fig. 5).

The spatial dependence of each Fi is shown in Fig. 6 by its spatial average (thin black line) and its standard deviation (vertical bars). As expected from the small variability of α̃, the scatter of Fi is small, so small, in fact, that the standard deviation is only perceptible for NPI, which has the largest variance at high frequencies and is therefore most sensitive. The pressure (KOE) indices are hardly changed from their unfiltered i, because their variance is dominated by interannual and decadal time scales, and are not affected by the filtering operation. Thus, our calculations and the forcing footprints γi are insensitive to the exact values of α̃ in the filtering of the monthly forcing indices i to the annual averages Fi, and much the same results are obtained by filtering the forcing functions with a constant value of α̃−1 of 3/4 yr.

e. Orthogonalization matrix

The forcing time series Fi are not independent, and we have to account for their shared variance. While the independent components of each Fi are sufficient to determine their unique regression patterns their correlated components imply that any particular forcing has an indirect impact on SSTA through the other time series. To quantify these relationships, we determine the rotation matrix 𝗕ij that relates Fi to an orthogonal base F*j,

 
formula

with rotated regression pattern γ*j = γi𝗕ij. The matrix 𝗕ij is obtained from a Gram–Schmidt orthogonalization that leaves the Niño-3.4 time series unchanged, removes from NPI the regression with Niño-3.4, removes from PDEL the multilinear regression with Niño-3.4 and NPI, and so on. This order reflects the hypothesis that anomalies of the Aleutian low pressure, the NPI, result in part from teleconnections from the Tropics (i.e., Niño-3.4) and that anomalies in the Kuroshio–Oyashio Extension are in part a response to forcing of ENSO and the NPI.

The rotation matrix 𝗕ij (Table 2) is estimated with both Fi and F*j normalized to unit variance, so that the square of any element 𝗕ij represents the variance of Fi captured by the F*j. Only the following two nontrival relationships emerge: Niño-3.4 is negatively correlated with NPI and accounts for 30% of its variance, and PDEL is weakly related to PAVG and accounts for 6% of its variance; all of the remaining elements are small. Thus, the rotated F*j, shown in Fig. 6 as dashed lines, only differ from Fi for NPI and PAVG.

Table 2.

Rotation matrix 𝗕ij that converts the orthogonal forcing time series F*j, denoted by top row to the forcing time series Fi, denoted by the leftmost column. The original and rotated time series are all normalized to unit standard deviations, so that the variance explained by F*j is given by the square of the values of the table. The values are the spatial mean of the rotation matrix at every location; however, the spatial variations are very small, and the standard deviation of the nonzero elements are of the order of 0.001 or smaller.

Rotation matrix 𝗕ij that converts the orthogonal forcing time series F*j, denoted by top row to the forcing time series Fi, denoted by the leftmost column. The original and rotated time series are all normalized to unit standard deviations, so that the variance explained by F*j is given by the square of the values of the table. The values are the spatial mean of the rotation matrix at every location; however, the spatial variations are very small, and the standard deviation of the nonzero elements are of the order of 0.001 or smaller.
Rotation matrix 𝗕ij that converts the orthogonal forcing time series F*j, denoted by top row to the forcing time series Fi, denoted by the leftmost column. The original and rotated time series are all normalized to unit standard deviations, so that the variance explained by F*j is given by the square of the values of the table. The values are the spatial mean of the rotation matrix at every location; however, the spatial variations are very small, and the standard deviation of the nonzero elements are of the order of 0.001 or smaller.

Note that the Gram–Schmidt orthogonalization could start with NPI and attribute the covariance of NPI and Niño-3.4 to changes of the NPI. Physically, this suggests that the extratropics modulate ENSO via teleconnections to the Tropics (Barnett et al. 1999a; Pierce et al. 2000) or provide stochastic forcing for ENSO (Vimont et al. 2001, 2003). However, it seems unlikely that ENSO variability is completely determined by the NPI, and so we do not treat this option in detail, although it is statistically appealing.

f. Forcing footprints

The forcing terms γi(x, y) Fi(x, y, t) in Eq. (5) represent the annually averaged SSTA that is the result of the forcing i after 1 yr. Results are shown for Fi(x, y, t) normalized to unit standard deviation at every location, so that the maps of the regression coefficients γi(x, y) show the pattern of SSTA caused by a typical amplitude of the filtered forcing time series Fi. Overall, the patterns are consistent with underlying physical processes and the results from prior studies.

The NPI regression pattern dominates the regression (Fig. 8, bottom), with large amplitudes over the western North Pacific and the Kuroshio–Oyashio frontal regions, surrounded by anomalies of the opposite sign in the eastern Pacific and in the Gulf of Alaska. A high sea level pressure perturbation of the Aleutian low is associated with increased temperatures in the central North Pacific along 40°N, and cooler temperatures in the eastern Pacific, the Gulf of Alaska, and the subpolar region (Davis 1976). These changes are associated with anomalies of the Pacific storm track (Nakamura et al. 2002; An and Wang 2005), and are consistent with a decrease of westerlies that reduces turbulent heat losses and causes poleward Ekman advection in the central North Pacific. Cooling tendencies in the east are consistent with advection of cold and dry air from the high latitudes (Tanimoto et al. 2003).

Fig. 8.

Forcing footprints γi(x, y) of (top) Niño-3.4 and (bottom) NPI. The maps represent the response of SST (K) after the forcing of unit magnitude has been applied for 1 yr. Contour interval is 0.1, and negative contours are dotted.

Fig. 8.

Forcing footprints γi(x, y) of (top) Niño-3.4 and (bottom) NPI. The maps represent the response of SST (K) after the forcing of unit magnitude has been applied for 1 yr. Contour interval is 0.1, and negative contours are dotted.

The pattern associated with Niño-3.4 (Fig. 8, top) has the most loading in the eastern Pacific and along the coast of North America, with opposite-signed anomalies farther offshore. This pattern captures coastal wave processes and atmospheric teleconnections that are independent of the NPI, and, as we will show below, together with the correlated contribution of the NPI, forms the well-known SST response pattern to tropical variations.

The zonal advection index in the northwestern Pacific (PDEL) is associated with temperature anomalies in the Kuroshio–Oyashio Extension (Fig. 9, top), and overlaps with the regions of influence of the Niño-3.4 and NPI. A positive index, corresponding to an increase of eastward flow, leads to warming centered west of the date line, spreading along 40°N. The warming is shifted to the east-northeast of the region of the ocean pressure index (indicated by dotted lines in Fig. 9), consistent with downstream advection by the mean flow, which has a east-northeasterly direction in the Kuroshio–Oyashio Extension (Qiu 2003).

Fig. 9.

Regression pattern with Kuroshio–Oyashio Extension circulation indices of (top) zonal advection and (bottom) thermocline depth. Units are in kelvins and correspond to the response of SST after the forcing of unit magnitude has been applied for 1 yr. Contour interval is 0.1, negative contours are dotted, and the zero contour is omitted. Thin dashed lines indicate the areas used in the calculation of the indices.

Fig. 9.

Regression pattern with Kuroshio–Oyashio Extension circulation indices of (top) zonal advection and (bottom) thermocline depth. Units are in kelvins and correspond to the response of SST after the forcing of unit magnitude has been applied for 1 yr. Contour interval is 0.1, negative contours are dotted, and the zero contour is omitted. Thin dashed lines indicate the areas used in the calculation of the indices.

The thermocline depth index (PAVG) in the Kuroshio–Oyashio Extension shows little effect on annual-averaged SSTA (Fig. 9, bottom). Yet, the dipole pattern of warm anomalies at 42°N and cool anomalies at 35°N in response to a positive anomaly of PAVG is consistent with a northward shift of the Kuroshio–Oyashio Extension axis. The weak signal is consistent with the annual cycle of mixed-layer depth, which links thermocline depth and surface temperature only during the wintertime, when mixed layers are deep (Schneider and Miller 2001). For annual averages, this process accounts for only a small fraction of the variance.

g. Rotated forcing footprints

According to Eq. (6), the rotated forcing footprint γ* of Niño-3.4* augments the unrotated ENSO footprint with the unrotated pattern of NPI, PDEL, and PAVG, with weights 𝗕ij given by the first column of Table 2. The resulting Niño-3.4* pattern (Fig. 10, top) extends the unrotated Niño-3.4 pattern (Fig. 8, top) to the west, primarily by the addition of the unrotated NPI patterns (Fig. 8, bottom). It is the canonical response of the North Pacific SST to ENSO and corresponds to the composite ENSO response of the North Pacific to the ensemble-averaged response obtained from models, and to the pattern removed from North Pacific SSTA by regression with an index for the cold tongue SSTA (Zhang et al. 1997).

Fig. 10.

Rotated footprints for γi(x, y)* of (top) Niño-3.4* and (bottom) NPI*. The maps represent the response of SST (K) after the forcing of unit magnitude has been applied for 1 yr. Contour interval is 0.1, and negative contours are dotted.

Fig. 10.

Rotated footprints for γi(x, y)* of (top) Niño-3.4* and (bottom) NPI*. The maps represent the response of SST (K) after the forcing of unit magnitude has been applied for 1 yr. Contour interval is 0.1, and negative contours are dotted.

The rotated γ* of Niño-3.4* captures the annually averaged footprint resulting from the atmospheric bridge (Alexander et al. 2002, 2004) between ENSO (Niño-3.4) and the Aleutian low (NPI) and the atmosphere of the eastern North Pacific, and resulting from the oceanic connection by coastal wave propagation off of North America equatorward of 40°N (Chelton and Davis 1982; Clarke and Lebedev 1999; Lluch-Cota et al. 2001). Overall, El Niño is associated with low temperatures in the central North Pacific in boreal fall and winter, and along the Kuroshio–Oyashio Extension in boreal summer (Alexander et al. 2004), and with warm anomalies in the eastern Pacific and in the Gulf of Alaska. Alexander et al. (2002) attribute the anomalies of SST to changes of the surface heat flux, entrainment, and Ekman advection, together of the order of 20 W m−2 in the cooling region east of the date line. Similar numbers are obtained from ship observations (Tanimoto et al. 2003). To compare this heat flux with the magnitude of the regression pattern of 0.2–0.3 K, we interpret the time series of Niño-3.4 (Fig. 1) as a heat flux time series, with anomalies during ENSO of 20 W m−2 that are typical during ENSO events and equivalent to the two standard deviations. Applying this heat flux to a layer of 80-m depth, and estimating the corresponding F time series via the definition implied by Eq. (5), yields a standard deviation of annual-averaged values of SSTA of 0.4 K, which are roughly consistent with Fig. 10 (top).

The rotated NPI* pattern is almost identical to its unrotated counterpart, because the contributions from PDEL and PAVG are small (second column of Table 2). The γ* associated with NPI* indicates that anomalies at 40°N (Fig. 10, bottom) are consistent with decadal wintertime temperature anomalies independent of ENSO, resulting from shifts of the North Pacific subarctic fronts (Nakamura et al. 1997; Nakamura and Yamagata 1999; Nakamura and Kazmin 2003).

The rotated pattern of PDEL (not shown) is enhanced slightly between 40° and 45°N west of the date line because of the contribution from PAVG (Table 2, third column). Applying instead the variance that is shared between the two KOE indices to PAVG enhances its northern loading. However, this increase is not sufficient to raise the variance that is associated with PAVG beyond the significance level.

5. PDO properties

a. Kinematics

Using the autoregressive model for SSTA and the regression patterns, we now investigate the evolution equation for the PDO. The equation governing the time evolution of the PDO is derived by projecting Eq. (5) on the leading empirical orthogonal function of annual-averaged values of temperature anomalies

 
formula

where ak is the principal component and Φk the orthonormal spatial loading pattern of mode k, with the PDO being k = 0, and 〈 〉 denoting the spatial projection. In general, the spatial dependence of α couples the evolution of the PDO and the other principal components. The forcing of the PDO amplitude, the last term on the right-hand side of Eq. (7), can be written for each i as the product of the regression coefficient 〈γiΦ0〉 and the weighted spatial average of the forcing time series 〈γiΦ0−1γiΦ0Fi〉. Because the spatial dependence of Fi is weak, the weighting has little effect, and the forcing time series for the PDO are indistinguishable from the spatial means of Fi(x, y, t) (cf. black and gray lines in Fig. 6).

Using the estimates of α and γi, and the observed Φ0, the coefficients of Eq. (7) are obtained and presented in Table 3. The coupling of the PDO to higher modes is weak, and an autoregressive model for the PDO is appropriate, as assumed by Newman et al. (2003) and in Eq. (1). The damping coefficient for the PDO is 0.26, corresponding to a damping time of 8.9 months. This coefficient is about half the size of that obtained by Newman et al. (2003). In their analysis, the autocorrelation is determined after the regression with Niño-3.4 is removed from the PDO time series. Thus, their residual includes the effect of the low-frequency forcing by the ocean adjustment and therefore overestimates α.

Table 3.

Coefficients of Eq. (7) governing the PDO evolution estimated from the fit of the observed SSTA evolution to Eq. (4). The first row shows the projection onto observed patterns of the PDO; column 1 depicts the autoregressive parameter of the PDO amplitude, and rows 2–4 show the coupling to higher observed EOFs. The last four columns show the coefficients multiplying the forcing indices corresponding to Niño-3.4, NPI, PDEL, and PAVG. The second row shows results for the forcing fields orthogonalized with matrix 𝗕ij. All forcing time series were normalized to unit variance. Rows 3 and 4 show the results from fitting the observed time evolution of the PDO to Eq. (7), while neglecting the coupling to the other EOFs.

Coefficients of Eq. (7) governing the PDO evolution estimated from the fit of the observed SSTA evolution to Eq. (4). The first row shows the projection onto observed patterns of the PDO; column 1 depicts the autoregressive parameter of the PDO amplitude, and rows 2–4 show the coupling to higher observed EOFs. The last four columns show the coefficients multiplying the forcing indices corresponding to Niño-3.4, NPI, PDEL, and PAVG. The second row shows results for the forcing fields orthogonalized with matrix 𝗕ij. All forcing time series were normalized to unit variance. Rows 3 and 4 show the results from fitting the observed time evolution of the PDO to Eq. (7), while neglecting the coupling to the other EOFs.
Coefficients of Eq. (7) governing the PDO evolution estimated from the fit of the observed SSTA evolution to Eq. (4). The first row shows the projection onto observed patterns of the PDO; column 1 depicts the autoregressive parameter of the PDO amplitude, and rows 2–4 show the coupling to higher observed EOFs. The last four columns show the coefficients multiplying the forcing indices corresponding to Niño-3.4, NPI, PDEL, and PAVG. The second row shows results for the forcing fields orthogonalized with matrix 𝗕ij. All forcing time series were normalized to unit variance. Rows 3 and 4 show the results from fitting the observed time evolution of the PDO to Eq. (7), while neglecting the coupling to the other EOFs.

The weights of the forcing fields that are projected onto the PDO are determined with Fi normalized to unit standard deviation, and are shown as gray lines in Fig. 6. The forcing functions are indistinguishable from the individual Fi and their spatial means, showing the weak variations imparted by the spatial dependence of α̃. While the forcing by NPI clearly dominates, a significant fraction of the NPI response results from teleconnections with Niño-3.4. The orthogonalization to Niño-3.4* and NPI* yields about equal weights for these forcings, with opposite signs. Thus, ENSO primarily affects the PDO via changes of the Aleutian low (and NPI), and leads to the known relationships between ENSO and the PDO—warm conditions in the eastern equatorial Pacific are associated with a cool central North Pacific and a warm eastern North Pacific. Positive values of the NPI* indicate a reduced Aleutian low, with warmer temperatures in the central North Pacific (and, thus, negative anomalies of the PDO). The zonal advection index PDEL in the Kuroshio–Oyashio Extension has about one-third of the weight compared to the NPI*, and indicates that increased advection leads to warm conditions and negative values of the PDO.

As a consistency test, we fit the observed evolution of the PDO amplitude to Eq. (7), while neglecting the coupling to higher EOFs of SSTA. This yields a slightly reduced damping coefficient, and very similar coefficients of the forcings (Table 3, bottom two rows).

In the following, we will restrict the discussion to the rotated forcing fields and regression pattern, noting that a significant fraction of the ENSO forcing affects North Pacific SSTA and the PDO via the NPI.

b. PDO evolution

The evolution of the PDO results from the sum of the temporal accumulation of the forcing by Niño-3.4, and changes of the Aleutian low and the Kuroshio–Oyashio Extension, and their relative contributions can be estimated (Fig. 11). The shift in 1976/77 resulted primarily from Niño-3.4*, that is, changes in the Tropics (Trenberth 1990), with a small addition from changes of the NPI*. In contrast, the shift from positive to negative conditions in 1998 was a result of the NPI*, with a smaller contribution from Niño-3.4*.

Fig. 11.

(top) The PDO time series and reconstruction (gray) based on Eq. (7). (middle three) The contributions to the PDO from Niño-3.4*, NPI*, and P*DEL, respectively. (bottom) The index for thermocline depth estimate P*AVG does not explain a significant fraction of the SSTA variability anywhere and is therefore deemed insignificant for the reconstruction of the PDO. Dotted vertical lines mark the winters of 1976/77 and 1998/99.

Fig. 11.

(top) The PDO time series and reconstruction (gray) based on Eq. (7). (middle three) The contributions to the PDO from Niño-3.4*, NPI*, and P*DEL, respectively. (bottom) The index for thermocline depth estimate P*AVG does not explain a significant fraction of the SSTA variability anywhere and is therefore deemed insignificant for the reconstruction of the PDO. Dotted vertical lines mark the winters of 1976/77 and 1998/99.

Overall, the reconstruction recovers about 83% of the PDO variance. The identification of the variance that is accounted for by individual time series is complicated by the fact that, even though the forcing time series are orthogonal, their integrated contributions to the PDO are not. Thus, we estimate the variance that is explained by a particular index by the difference of the PDO variance obtained with and without this index. This shows that the unrotated Niño-3.4 and NPI account for 21% and 73%, respectively, while the rotated Niño-3.4* and NPI* account for approximately equal fractions of the PDO variance, 40% and 44%, respectively. The contributions of the zonal flow index PDEL and P*DEL are much smaller, and correspond to 6% and 7%, respectively of the year-to-year PDO variance. Thus, the PDO is, to a large part, a reflection of the NPI, with Niño-3.4 affecting the PDO both through the NPI and directly.

c. PDO spectrum

The frequency dependence of the forcing processes is best seen by a reconstruction of the frequency spectrum of the PDO (Fig. 12). To this end, results are applied to the monthly time series, and the spectral power P2 of the PDO is estimated as

 
formula

where ω is the frequency, α̃0 = −δ−1 loge(〈αΦ0Φ0〉), and P2i are the power spectra of the forcing indices. The spectrum of the reconstruction has a very good agreement with the observed spectrum (Fig. 12), with a slight overestimation at interannual frequencies, which can also be gleaned from the comparison of the reconstructed and observed PDO time series (Fig. 11). At annual and shorter time scales, the reconstructed spectrum slightly, but systematically, underestimates the observed spectrum. This is consistent with reemergence, observed in the PDO time series (Newman et al. 2003), and its impact on the spectrum (appendix B).

Fig. 12.

Power spectrum as a function of frequency in cycles per year of the observed and reconstructed PDO, and its contributions (ω2 + α20)−1γ2iP2i resulting from NPI*, Niño34*, and P*DEL (as indicated). The solid black line denotes the spectrum of the observed PDO, the light gray line with black center that of the reconstruction. Spectra have been smoothed by three successive applications of a five-point running mean. Fluctuations of the spectra for frequencies above 0.6 cpy are not significant.

Fig. 12.

Power spectrum as a function of frequency in cycles per year of the observed and reconstructed PDO, and its contributions (ω2 + α20)−1γ2iP2i resulting from NPI*, Niño34*, and P*DEL (as indicated). The solid black line denotes the spectrum of the observed PDO, the light gray line with black center that of the reconstruction. Spectra have been smoothed by three successive applications of a five-point running mean. Fluctuations of the spectra for frequencies above 0.6 cpy are not significant.

The contribution (ω2 + α20)−1γ2iP2i of each of the forcing components i shows that the NPI* contribution dominates the response at annual and higher frequencies, because the spectrum of the NPI* itself is primarily white, with only a slight reddening tendency expected from thermodynamic air–sea coupling (Barsugli and Battisti 1998; Blade 1997), and a broad, gentle peak at interannual time scales, likely the result of teleconnections to the Tropics. The spectrum of Niño-3.4* has most power at interannual frequencies, with a strong low-frequency component. Its contribution to the PDO increases rapidly, approximately as ω−4, toward interannual frequencies, and then decreases slightly toward decadal time scales. The spectrum of the zonal advection index P*DEL is red with all of its power at lower-than-interannual frequencies, and P*DEL adds to the power of the PDO only at decadal time scales. There, however, its contribution to the PDO is on par with the NPI* and Niño-3.4*. Thus, the frequency spectrum of the PDO (Fig. 12) reflects the effects of the NPI* at periods of a year and higher. At interannual frequencies, Niño-3.4* and NPI* explain the variability, while at lower frequencies NPI*, Niño-3.4*, and P*DEL are of equal importance.

d. Relation of PDO and regression patterns

Because the evolution of SSTA over the North Pacific is determined by three different forcing indices and patterns, so is the PDO. The regression patterns have large overlap, and are not orthogonal, and an empirical orthogonal function expansion cannot separate these into distinct modes. In fact, we hypothesize that the PDO and the relative roles of its individual forcings are determined primarily by the spatial covariance of the forcing footprints, even if the time evolution of the forcing indices are independent.

To test this hypothesis, we note that the solution to Eq. (4) is separable if α is a constant. This approximation is consistent with the lack of coupling of the PDO with higher modes. A reconstruction of SSTA for α replaced at all locations with α0 = 〈αΦ0Φ0〉 results in almost identical estimates of the PDO (Fig. 2), with a pattern correlation of 0.98 and a temporal correlation of 0.93 with the observed PDO (Table 1).

The solution for the annual-averaged temperature anomaly at every location and time step t = is obtained by integration of Eq. (5),

 
formula

and is the sum of separable products of the forcing footprints γi and the time integrals fi of the AR-1 solution to each forcing index. The initial values fi are obtained from a multilinear fit of SSTA to the γi patterns at the initial time n = 0.

The covariance matrix 𝗥 of temperature anomalies, whose leading eigenvector is the PDO, is then

 
formula

where , γ*, and 𝗳* are written in matrix notation, with prime denoting the transpose. The matrix 𝗥 is constructed from the covariance of the regression fields γ*, weighted by the covariance (correlation) of 𝗳. The latter is dominated by the diagonal elements. To test whether the spatial pattern of the PDO is determined by the covariance of γ even if the forcing indices are independent, a test covariance matrix 𝗥̈ is made using Eq. (10) by setting off-diagonal elements of 𝗳′𝗳 to zero. Again, the PDO pattern from this test matrix is almost identical to the observed pattern, with spatial and temporal correlations with observations of 0.98 and 0.93, respectively (Table 1), suggesting that the pattern of the PDO and the weights of the separate forcings result from the nonorthogonality of the forcing footprints. In fact, the pattern of the PDO can be recovered by considering only Niño-3.4* and NPI* in the estimation of 𝗥 (Table 1, rightmost column) with a reduction of the skill of the principal component reconstruction commensurate with the lack of the contribution of PDEL.

e. Robustness of PDO reconstruction

Deser et al. (2004) report that the relationship between winter values of Trenberth and Hurrell’s (1994) NPIT and Niño-3.4, or the CTI (defined as SSTA in 6°S–6°N, 180°–90°W, with the global mean SSTA subtracted) was weaker in the first half of the twentieth century than in the second half. However, the OTI—the leading EOF of a number of indices for western Pacific and eastern Indian Ocean SSTA, anomalies of cloudiness and rainfall close to the South Pacific convergence zone, anomalies of stratus clouds in the eastern tropical Pacific, and sea level pressure in the tropical Pacific and Indian Oceans—shows consistent variations with the NPI throughout the twentieth century.

For the PDO, the question arises if its reconstruction is sensitive to the choice of the tropical index, and if the relationship between PDO and its forcing is robust in the earlier part of the twentieth century. To this end, 1900–97 time series of the observed PDO, NPIT, and CTI or OTI are fitted to Eq. (4). Indices for the KOE cannot be extended into the earlier half of the century because the wind stress estimates are not available, and forcing of Eq. (2) by a reconstruction of the Pacific wind stress from regressions with Niño-3.4 and NPI fails to recover the KOE indices (see section 3b).

The relationship of the PDO to the tropical indices alone is sensitive to the choice of the latter. Fitting the observed evolution of the PDO to the CTI only (Newman et al. 2003) leads a hindcast correlation of 0.52, while the fit of PDO to OTI captures more of the variance with a hindcast correlation of 0.65. However, when Eq. (4) is forced by NPIT and either CTI or OTI, the hindcast skills for the time period 1900–97 are 0.74 and 0.68, respectively. Thus, the variance of the NPIT that can be attributed to tropical forcing is sensitive to the tropical index, but the tropical indices do not add any information to the reconstruction of the PDO once the NPIT (or NPI) is provided. The evolution of the PDO is largely determined by the NPI (Davis 1976), while the changes of the NPI are attributed, in part, to teleconnections from the Tropics (Deser et al. 2004). However, please note that this does not imply that changes of the Aleutian low physically account for all changes of the PDO. Rather, prescription of the NPI captures part of the tropical to extratropical teleconnections through the correlation of NPI and tropical indices.

The skill of the fit is degraded in the earlier part of the twentieth century, compared to the latter, with correlations for the period before (after) 1950 being 0.67 (0.85) for CTI and 0.61 (0.88) for OTI. Some degradation is expected because all indices are constrained by the fewer observations in the earlier period.

f. PDO predictability

Using the framework that is presented above, a forecast of the PDO requires predictions of NPI, ENSO, and the gyre processes in the Kuroshio–Oyashio Extension. ENSO forecasts have skill at lead times of 1 and possibly up to 2 yr (Barnston et al. 1999; Chen et al. 2004). The ENSO forecasts fix part of the NPI; the component of NPI* resulting from the intrinsic atmospheric variability of the Aleutian low is expected to have very limited predictability, likely less than a year. The predictability of remote forcing from other areas, such as perturbations in the western Pacific or eastern Indian Ocean (Deser et al. 2004), remains to be determined. In the Kuroshio–Oyashio Extension, predictive skill with a lead time of 1–2 yr has been documented for wintertime SSTA by virtue of the Rossby wave propagation (Schneider and Miller 2001). The predictability of PDEL itself should be better, however, it only accounts for a small fraction of the year-to-year variability of the PDO. Thus, it appears that a forecast of annually averaged values of the PDO is possible at a lead time of up to a few years. This analysis does not indicate much hope for a skillful forecast of the PDO at multiyear to decadal time scales.

6. Conclusions

Causes of decadal variability of SSTA in the North Pacific, and its leading empirical orthogonal function, have long been the subject of intense debate. Here, we show that the winter-centered, annually averaged values of SSTA and the PDO can be reconstructed with high fidelity from an autoregressive model forced by indices that track sea level pressure of the Aleutian low, ENSO, and ocean circulation anomalies in the Kuroshio–Oyashio Extension region. This implies that North Pacific SSTA and the PDO are a response to changes of the North Pacific atmosphere resulting from its intrinsic variability, remote forcing by ENSO and other processes, and ocean wave processes associated with ENSO and the adjustment of the North Pacific Ocean by Rossby waves. The relative importance of these forcing processes for the PDO is frequency dependent: at periods shorter than 1 yr, intrinsic North Pacific variability dominates; at interannual frequencies, changes of the Aleutian low and of ENSO are essential; and at decadal periods, the ocean circulation anomalies in the Kuroshio–Oyashio Extension gain in importance to par with the other processes.

Under this hypothesis, the PDO is not governed by a single physical process that defines a climate mode, akin to ENSO, but results from at least three different processes. Of course, correlation is not casuality, and statistical relations like this one can only be disproved, not proved, but we have presented physical arguments for the hypothesized influences and have shown that randomly chosen red-noise processes do not approach the skill of the three indices used here. This hypothesis also exemplifies the well-known limitation of empirical orthogonal functions that represent the variance optimally, but are not designed to separate multiple forcing fields with nonorthogonal spatial characteristics, as in the case of the PDO.

These results suggest a model for the origin of the PDO and quantify that both tropical and North Pacific processes are at play. The famous 1976/77 shift of the PDO is of tropical origin, as suggested in many previous studies (e.g., Nitta and Yamada 1989; Trenberth and Hurrell 1994; Graham 1994). The symmetry of the SSTA pattern between the northern and southern extratropics is a reflection of the common tropical forcing (Evans et al. 2001).

This study further clarifies the causes of the “impact” of the PDO on climate of the adjacent continental areas and teleconnection patterns from the Tropics. As found from modeling experiments by Pierce (2002), the PDO does not excite the climate modulation, but the PDO and these climate anomalies share the same forcing from the Tropics and from intrinsic variability of the Aleutian low. The impact of oceanic anomalies in the Kuroshio–Oyashio Extension complicates this picture at decadal time scales. The atmospheric response to these anomalies is likely small at best, and we do not expect a large trace of these in the climate over North America (Kushnir et al. 2002; Yulaeva et al. 2001). Furthermore, the relative roles of forcing are likely to be site dependent, and our results suggest, therefore, that the stratification of climate anomalies or teleconnection patterns be based on the underlying indices of ENSO and the NPI, rather than ENSO and the PDO.

A similar conclusion applies to attempts to reconstruct Northern Hemisphere climate from paleodata, such as tree rings (Biondi et al. 2001; Gedalof et al. 2002). The atmospheric conditions controlling tree growth and tree ring properties are a function of the tropical Pacific and NPI, rather than of the PDO. Thus, conditions in the tropical Pacific and the NPI can be reconstructed from an array of tree rings. A reconstruction of the PDO alone will likely miss the impact of ocean circulation anomalies in the Kuroshio–Oyashio Extensions, and has to be viewed with caution for decadal time scales that correspond to the adjustment time of the extratropical ocean circulation.

Fig. B1. Autocorrelation of SST based on the idealized reemergence model, Eq. (B1) for H1 = 100 m and H2 = 40 m. The x axis denotes the month, with winter occupying months 1 through 6, and summer the remainder of the year. The y axis denotes the lag in months. Contour level is 0.2.

Fig. B1. Autocorrelation of SST based on the idealized reemergence model, Eq. (B1) for H1 = 100 m and H2 = 40 m. The x axis denotes the month, with winter occupying months 1 through 6, and summer the remainder of the year. The y axis denotes the lag in months. Contour level is 0.2.

Fig. B2. Spectrum of the simulated SST based on Eq. (B1) for the reference run H1 = H2 = 100 m (without reemergence, solid thick line), and with H2 = 40 m (dotted spectrum). The damping time scale for the winter mixed layer is marked by the solid vertical line, while the damping time scale for the summer mixed layer is shown as dotted vertical lines. The dashed lines correspond to the damping rate α̃ derived from winter-centered annual-averaged data, for the summer mixed-layer depth H2 = 40 m, and for the model without reemergence (H2 = H1).

Fig. B2. Spectrum of the simulated SST based on Eq. (B1) for the reference run H1 = H2 = 100 m (without reemergence, solid thick line), and with H2 = 40 m (dotted spectrum). The damping time scale for the winter mixed layer is marked by the solid vertical line, while the damping time scale for the summer mixed layer is shown as dotted vertical lines. The dashed lines correspond to the damping rate α̃ derived from winter-centered annual-averaged data, for the summer mixed-layer depth H2 = 40 m, and for the model without reemergence (H2 = H1).

Acknowledgments

The authors thank Drs. Arthur J. Miller, Jim Potemra, Bo Qiu, and Shang-Ping Xie for helpful discussions, and Drs. Clara Deser, Nathan Mantua, Hisashi Nakamura, and one anonymous reviewer for comments on an earlier draft. This research was supported by the National Science Foundation Grant OCE00-82543, by the Office of Science (BER), U.S. Department of Energy Grants DE-FG03-01ER63255 and DE-FG02-04ER63862, and by the National Oceanic and Atmospheric Administration Grant NA17RJ1231. The views expressed herein are those of the authors and do not necessarily reflect the views of NSF, DOE, or NOAA, or any of their subagencies.

REFERENCES

REFERENCES
Alexander
,
M. A.
,
1992
:
Midlatitude atmosphere–ocean interaction during El Niño. Part I: The North Pacific.
J. Climate
,
5
,
944
958
.
Alexander
,
M. A.
,
I.
Blade
,
M.
Newman
,
J. R.
Lanzante
,
N. C.
Lau
, and
J. D.
Scott
,
2002
:
The atmospheric bridge: The influence of ENSO teleconnections on air–sea interaction over the global oceans.
J. Climate
,
15
,
2205
2231
.
Alexander
,
M. A.
,
N-C.
Lau
, and
J. D.
Scott
,
2004
:
Broadening the atmospheric bridge paradigm: ENSO teleconnections to the North Pacific in summer and to the tropical west Pacific–Indian Oceans over the seasonal cycle. Earth Climate: The Ocean-Atmosphere Interaction, Geophys. Monogr., No. 147, Amer. Geophys. Union, 85–104
.
An
,
S-I.
, and
B.
Wang
,
2005
:
The forced and intrinsic low-frequency modes of the North Pacific.
J. Climate
,
18
,
876
885
.
Barlow
,
M.
,
S.
Nigam
, and
E. H.
Berbery
,
2001
:
ENSO, Pacific decadal variability, and U.S. summertime precipitation, drought, and streamflow.
J. Climate
,
14
,
2105
2128
.
Barnett
,
T. P.
,
D. W.
Pierce
,
M.
Latif
,
D.
Dommenget
, and
R.
Saravanan
,
1999a
:
Interdecadal interactions between the tropics and midlatitudes in the Pacific basin.
Geophys. Res. Lett.
,
26
,
615
618
.
Barnett
,
T. P.
,
D. W.
Pierce
,
R.
Saravanan
,
N.
Schneider
,
D.
Dommenget
, and
M.
Latif
,
1999b
:
Origins of the midlatitude Pacific decadal variability.
Geophys. Res. Lett.
,
26
,
1453
1456
.
Barnston
,
A. G.
,
Y.
He
, and
M. H.
Glantz
,
1999
:
Predictive skill of statistical and dynamical climate models in SST forecasts during the 1997–98 El Niño episode and the 1998 La Niña onset.
Bull. Amer. Meteor. Soc.
,
80
,
217
244
.
Barsugli
,
J. J.
, and
D. S.
Battisti
,
1998
:
The basic effects of atmosphere–ocean thermal coupling on midlatitude variability.
J. Atmos. Sci.
,
55
,
477
493
.
Biondi
,
F.
,
A.
Gershunov
, and
D. R.
Cayan
,
2001
:
North Pacific decadal climate variability since 1661.
J. Climate
,
14
,
5
10
.
Blade
,
I.
,
1997
:
The influence of midlatitude ocean–atmosphere coupling on the low-frequency variability of a GCM. Part I: No tropical SST forcing.
J. Climate
,
10
,
2087
2106
.
Cayan
,
D. R.
,
M. D.
Dettinger
,
H. F.
Diaz
, and
N. E.
Graham
,
1998
:
Decadal variability of precipitation over western North America.
J. Climate
,
11
,
3148
3166
.
Chelton
,
D. B.
, and
R. E.
Davis
,
1982
:
Monthly mean sea-level variability along the west coast of North America.
J. Phys. Oceanogr.
,
12
,
757
784
.
Chen
,
D.
,
M. A.
Cane
,
A.
Kaplan
,
S. E.
Zebiak
, and
D. J.
Huang
,
2004
:
Predictability of El Nino over the past 148 years.
Nature
,
428
,
733
736
.
Clarke
,
A. J.
, and
A.
Lebedev
,
1999
:
Remotely driven decadal and longer changes in the coastal Pacific waters of the Americas.
J. Phys. Oceanogr.
,
29
,
828
835
.
Davis
,
R. E.
,
1976
:
Predictability of sea surface temperature and sea level pressure anomalies over the North Pacific Ocean.
J. Phys. Oceanogr.
,
6
,
249
266
.
Deser
,
C.
,
M. A.
Alexander
, and
M. S.
Timlin
,
1999
:
Evidence for a wind-driven intensification of the Kuroshio Current Extension from the 1970s to the 1980s.
J. Climate
,
12
,
1697
1706
.
Deser
,
C.
,
M. A.
Alexander
, and
M. S.
Timlin
,
2003
:
Understanding the persistence of sea surface temperature anomalies in midlatitudes.
J. Climate
,
16
,
57
72
.
Deser
,
C.
,
A. S.
Phillips
, and
J. W.
Hurrell
,
2004
:
Pacific interdecadal climate variability: Linkages between the tropics and North Pacific during boreal winter since 1900.
J. Climate
,
17
,
3109
3124
.
Dettinger
,
M. D.
,
D. R.
Cayan
,
H. F.
Diaz
, and
D. M.
Meko
,
1998
:
North–south precipitation patterns in western North America on interannual-to-decadal timescales.
J. Climate
,
11
,
3095
3111
.
Dommenget
,
D.
, and
M.
Latif
,
2002
:
Analysis of observed and simulated SST spectra in the midlatitudes.
Climate Dyn.
,
19
,
277
288
.
Evans
,
M. N.
,
M. A.
Cane
,
D. P.
Schrag
,
A.
Kaplan
,
B. K.
Linsley
,
R.
Villalba
, and
G. M.
Wellington
,
2001
:
Support for tropically-driven Pacific decadal variability based on paleoproxy evidence.
Geophys. Res. Lett.
,
28
,
3689
3692
.
Frankignoul
,
C.
,
P.
Müller
, and
E.
Zorita
,
1997
:
A simple model of the decadal response of the ocean to stochastic wind forcing.
J. Phys. Oceanogr.
,
27
,
1533
1546
.
Garreaud
,
R. D.
, and
D. S.
Battisti
,
1999
:
Interannual (ENSO) and interdecadal (ENSO-like) variability in the Southern Hemisphere tropospheric circulation.
J. Climate
,
12
,
2113
2123
.
Gedalof
,
Z.
,
N. J.
Mantua
, and
D. L.
Peterson
,
2002
:
A multi-century perspective of variability in the Pacific Decadal Oscillation: New insights from tree rings and coral.
Geophys. Res. Lett.
,
29
.
2204, doi:10.1029/2002GL015824
.
Gershunov
,
A.
, and
T. P.
Barnett
,
1998
:
Interdecadal modulation of ENSO teleconnections.
Bull. Amer. Meteor. Soc.
,
79
,
2715
2725
.
Graham
,
N. E.
,
1994
:
Decadal-scale climate variability in the tropical and North Pacific during the 1970s and 1980s—Observations and model results.
Climate Dyn.
,
10
,
135
162
.
Hasselmann
,
K.
,
1976
:
Stochastic climate models. Part 1. Theory.
Tellus
,
28
,
473
485
.
Kalnay
,
E.
, and
Coauthors
,
1996
:
The NCEP/NCAR 40-Year Reanalysis Project.
Bull. Amer. Meteor. Soc.
,
77
,
437
471
.
Krishnan
,
R.
, and
M.
Sugi
,
2003
:
Pacific decadal oscillation and variability of the Indian summer monsoon rainfall.
Climate Dyn.
,
21
,
233
242
.
Kushnir
,
Y.
,
W. A.
Robinson
,
I.
Blade
,
N. M. J.
Hall
,
S.
Peng
, and
R.
Sutton
,
2002
:
Atmospheric GCM response to extratropical SST anomalies: Synthesis and evaluation.
J. Climate
,
15
,
2233
2256
.
Liebelt
,
P. B.
,
1967
:
An Introduction to Optimal Estimation.
Addison-Wesley, 273 pp
.
Lluch-Cota
,
D. B.
,
W. S.
Wooster
, and
S. R.
Hare
,
2001
:
Sea surface temperature variability in coastal areas of the northeastern Pacific related to the El Niño–Southern Oscillation and the Pacific Decadal Oscillation.
Geophys. Res. Lett.
,
28
,
2029
2032
.
Mantua
,
N. J.
, and
S. R.
Hare
,
2002
:
The Pacific decadal oscillation.
J. Oceanogr.
,
58
,
35
44
.
Mantua
,
N. J.
,
S. R.
Hare
,
Y.
Zhang
,
J. M.
Wallace
, and
R. C.
Francis
,
1997
:
A Pacific interdecadal climate oscillation with impacts on salmon production.
Bull. Amer. Meteor. Soc.
,
78
,
1069
1079
.
Miller
,
A. J.
, and
N.
Schneider
,
2000
:
Interdecadal climate regime dynamics in the North Pacific Ocean: Theories, observations and ecosystem impacts. Progress in Oceanography, Vol. 27, Pergamon, 257–260
.
Miller
,
A. J.
,
D. R.
Cayan
,
T. P.
Barnett
,
N. E.
Graham
, and
J. M.
Oberhuber
,
1994
:
Interdecadal variability of the Pacific Ocean: Model response to observed heat flux and wind stress anomalies.
Climate Dyn.
,
9
,
287
302
.
Miller
,
A. J.
,
D. R.
Cayan
, and
W. B.
White
,
1998
:
A westward intensified decadal change in the North Pacific thermocline and gyre-scale circulation.
J. Climate
,
11
,
3112
3127
.
Minobe
,
S.
,
2000
:
Spatio-temporal structure of the pentadecadal variability over the North Pacific. Progress in Oceanography, Vol. 47, Pergamon, 381–408
.
Nakamura
,
H.
, and
T.
Yamagata
,
1999
:
Recent decadal SST variability in the northwestern Pacific and associated atmospheric anomalies.
Beyond El Niño: Decadal and Interdecadal Climate Variability, Springer Verlag, 49–72
.
Nakamura
,
H.
, and
A. S.
Kazmin
,
2003
:
Decadal changes in the North Pacific oceanic frontal zones as revealed in ship and satellite observations.
J. Geophys. Res.
,
108
.
3078, doi:10.1029/1999JC000085
.
Nakamura
,
H.
,
G.
Lin
, and
T.
Yamagata
,
1997
:
Decadal climate variability in the North Pacific during recent decades.
Bull. Amer. Meteor. Soc.
,
78
,
2215
2225
.
Nakamura
,
H.
,
T.
Izumi
, and
T.
Sampe
,
2002
:
Interannual and decadal modulation recently observed in the North Pacific storm track activity and east Asian winter monsoon.
J. Climate
,
15
,
1855
1874
.
Newman
,
M.
,
G. P.
Compo
, and
M. A.
Alexander
,
2003
:
ENSO-forced variability of the Pacific decadal oscillation.
J. Climate
,
16
,
3853
3857
.
Nitta
,
T.
, and
S.
Yamada
,
1989
:
Recent warming of tropical sea surface temperature and its relationship to the Northern Hemisphere circulation.
J. Meteor. Soc. Japan
,
67
,
375
383
.
Pierce
,
D. W.
,
2001
:
Distinguishing coupled ocean–atmosphere interactions from background noise in the North Pacific. Progress in Oceanography, Vol. 49, Pergamon, 331–352
.
Pierce
,
D. W.
,
2002
:
The role of sea surface temperatures in interactions between ENSO and the North Pacific Oscillation.
J. Climate
,
15
,
1295
1308
.
Pierce
,
D. W.
,
T.
Barnett
, and
M.
Latif
,
2000
:
Connections between the Pacific Ocean tropics and midlatitudes on decadal time scales.
J. Climate
,
13
,
1173
1194
.
Pierce
,
D. W.
,
T. P.
Barnett
,
N.
Schneider
,
R.
Saravanan
,
D.
Dommenget
, and
M.
Latif
,
2001
:
The role of ocean dynamics in producing decadal climate variability in the North Pacific.
Climate Dyn.
,
18
,
51
70
.
Qiu
,
B.
,
2003
:
Kuroshio Extension variability and forcing of the Pacific decadal oscillations: Responses and potential feedback.
J. Phys. Oceanogr.
,
33
,
2465
2482
.
Schneider
,
N.
, and
A. J.
Miller
,
2001
:
Predicting western North Pacific ocean climate.
J. Climate
,
14
,
3997
4002
.
Schneider
,
N.
,
A. J.
Miller
, and
D. W.
Pierce
,
2002
:
Anatomy of North Pacific decadal variability.
J. Climate
,
15
,
586
605
.
Seager
,
R.
,
Y.
Kushnir
,
N.
Naik
,
M. A.
Cane
, and
J. A.
Miller
,
2001
:
Wind-driven shifts in the latitude of the Kuroshio–Oyashio Extension and generation of SST anomalies on decadal timescales.
J. Climate
,
14
,
4249
4265
.
Sturges
,
W.
, and
B. G.
Hong
,
1995
:
Wind forcing of the Atlantic thermocline along 32°N at low frequencies.
J. Phys. Oceanogr.
,
25
,
1706
1715
.
Tanimoto
,
Y.
,
H.
Nakamura
,
T.
Kagimoto
, and
S.
Yamane
,
2003
:
An active role of extratropical sea surface temperature anomalies in determining anomalous turbulent heat flux.
J. Geophys. Res.
,
108
.
3304, doi:10.1029/2002JC001750
.
Tomita
,
T.
,
B.
Wang
,
T.
Yasunari
, and
H.
Nakamura
,
2001
:
Global patterns of decadal-scale variability observed in sea surface temperature and lower-tropospheric circulation fields.
J. Geophys. Res.
,
106
,
26805
26815
.
Trenberth
,
K. E.
,
1990
:
Recent observed interdecadal climate changes in the Northern Hemisphere.
Bull. Amer. Meteor. Soc.
,
71
,
988
993
.
Trenberth
,
K. E.
, and
J. W.
Hurrell
,
1994
:
Decadal atmosphere–ocean variations in the Pacific.
Climate Dyn.
,
9
,
303
319
.
Vimont
,
D. J.
,
D. S.
Battisti
, and
A. C.
Hirst
,
2001
:
Footprinting: A seasonal connection between the tropics and mid-latitudes.
Geophys. Res. Lett.
,
28
,
3923
3926
.
Vimont
,
D. J.
,
D. S.
Battisti
, and
A. C.
Hirst
,
2003
:
Seasonal footprinting in the Pacific: Implications for ENSO.
J. Climate
,
16
,
2668
2675
.
Xie
,
S. P.
,
T.
Kunitani
,
A.
Kubokawa
,
M.
Nonaka
, and
S.
Hosoda
,
2000
:
Interdecadal thermocline variability in the North Pacific for 1958–1997: A GCM simulation.
J. Phys. Oceanogr.
,
30
,
2798
2813
.
Yulaeva
,
E.
,
N.
Schneider
,
D. W.
Pierce
, and
T.
Barnett
,
2001
:
Modeling of North Pacific climate variability forced by oceanic heat fluxes anomalies.
J. Climate
,
14
,
4027
4046
.
Zhang
,
Y.
,
J. M.
Wallace
, and
D. S.
Battisti
,
1997
:
ENSO-like interdecadal variability: 1900–93.
J. Climate
,
10
,
1004
1020
.

APPENDIX A

Least Squares Solution

To find the coefficients of Eqs. (1) and (5) at a point, we rewrite the problem in terms of the unknown parameters z = (α, γ1, γ2, . . .) and the matrix 𝗔 with columns that are the time series of each index, the temperature or PDO amplitude average of year tδ, and the forcing time series. The vector d is the time series of observed temperature or PDO amplitude at year t, and the residual is r,

 
formula

The coefficients z are determined by minimizing the cost function J,

 
formula

The norms ||· · ·||2 are chosen as the inverse covariance matrices

 
formula

where R is the covariance of the uncertainty in the residuals r, and P is the covariance of the model parameters z.

Here J is minimized at the point :

 
formula

(Liebelt 1967). The inverse of the curvature at the minimum is , and it determines the uncertainty of the solution—how much could the solution vary while still doing a good job on the objective function.

The solution and uncertainty for the overdetermined system of the application here is obtained by inversion in model space

 
formula

with uncertainty covariance

 
formula

In this study, we employ a large signal-to-noise ratio, that is, P−1 = 0 and R = I, so that the cost function involves the misfit (residual) variance only. This implies that the parameter variance far exceeds that of the residuals, and that the elements of z, or the columns of A, have equal weight. Our solutions are insensitive to any large, but finite, choice of signal-to-noise ratio.

APPENDIX B

The Impact of Reemergence

The focus on annual averages filters seasonal processes. The evolution of SST is particularly affected by the seasonal cycle of the mixed-layer depth. Temperature anomalies formed in deep winter and spring mixed layers persist under the seasonal thermocline during summer being insulated from the influence of the air–sea fluxes. During the erosion of the seasonal thermocline in the subsequent cold season, the deep anomalies are entrained into the mixed layers, and reemerge at the surface. This reemergence is found in the extratropical SSTA (e.g., Deser et al. 2003) and the PDO time series (Newman et al. 2003), and represents an important departure from the AR-1 hypothesis for SST anomalies. Surface temperature anomalies not only depend on the anomalies of the prior month, but also on the anomalies of the previous cold season. This influences the power of SSTA at biannual to shorter-than-annual time scales (Dommenget and Latif 2002).

To investigate the impact of reemergence on our analysis of winter-centered annual averages of SST, we explore a simple extension of the Hasselmann (1976) model that includes reemergence. Consider a mixed layer with depth H1 in winter and H2 in summer months, with a summer seasonal thermocline in the layer between H1 and H2 (H1 > H2). The rate of change of the surface temperature T is governed by

 
formula

with H = H1 from fall to spring, and H = H2 from spring to fall equinox. The air–sea heat flux consists of a forcing term Q, and feedback component of strength ε, normalized by the density of water ρ0 and its heat capacity cp. The last term on the right-hand side describes the erosion of the seasonal thermocline and is nonzero only during the entrainment season, which commences at the fall equinox tfall and deepens the mixed layer to winter values during the period δt. This switch is indicated by the function Λ that is one only if its argument is true. The temperature of the seasonal thermocline Tst equals the surface values from the end of the fall entrainment period to spring equinox when the seasonal thermocline is nonexistent and absorbed in the deep winter mixed layer. At the spring equinox at tspring the mixed layer abruptly shallows to its summer depth, and the seasonal thermocline becomes insulated from the surface fluxes. Thus, Tst remains at its spring value throughout the summer until the seasonal thermocline is reentrained into the mixed layer in the fall.

This model was implemented with parameters that are typical for the North Pacific: ε = 15 W m2 K−1, ρ0 = 1000 kg m−3, cp = 4000 J kg−1 K−1, H1 = 100 m, tfall = 0 months, tspring = 6 months, and Q being Gaussian white-noise Q with a standard deviation of 20 W m−2. As the simplest case, we focus on a short entrainment season, and choose δt to be zero. In the fall, the mixed layer instantly deepens to its winter depth. This idealized case recovers the qualitative influence of reemergence. More gradual spring shallowing and fall entrainment seasons would simply smooth the autocorrelation function shown below. The model was run for 1024 yr for the same realization of the random Q forcing for a reference case, H2 = H1, identical to the Hasselmann (1976) model, and for summertime mixed-layer depths of H2 = 40 m and of H2 = 20 m.

The autocorrelation function of simulated T has the characteristic structure of reemergence. While for the reference case the autocorrelation is independent of the calendar month, shallow summer mixed layers result in an autocorrelation for winter months that have secondary peaks in subsequent winters (Fig. B1). These secondary peaks are exacerbated for shallower summer mixed layers because the heat anomaly in the seasonal thermocline is larger relative to the heat anomaly stored in the summer mixed layer.

The spectrum of T for the reference case (Fig. B2, thick line) shows the expected frequency to the “−2” slope for frequencies below ε(2πρ0cpH1)−1 (denoted by a vertical thick line). Shallower summer thermoclines enhance the power at time scales of a few years and shorter (dotted spectrum in Fig. B2), but do not affect the variance of SST at low frequencies. Furthermore, spectra of winter-centered annual averages of SST faithfully reproduce the low-frequency spectrum of T (not shown).

To explore the relationship between the damping rate α̃ derived from annual-averaged values and the damping rates ε(ρ0cpH)−1 for winter H = H1 and summer H = H2 (shown in Fig. B2 as the vertical, dotted line), Eq. (4) is fit to winter-centered annual averages of the simulated T and = Q. The estimate of α̃ is identical to the winter values for the reference run (dashed thick vertical line in Fig. B2). For cases II and III (not shown), the value of α lies between the summer and winter damping rates, and is closer to the former. The exact dependence of the effective annual damping rate on the seasonal evolution of the mixed layer is nonlinear and complicated [see Eq. (B1)]. For our purposes it suffices that the damping rate α̃−1 is shorter than, but close to winter values.

Footnotes

* International Pacific Research Center Contribution Number 335 and School of Ocean and Earth Science and Technology Contribution Number 6594

Corresponding author address: N. Schneider, IPRC, University of Hawaii at Manoa, 1680 East West Road, Honolulu, HI 96822. Email: nschneid@hawaii.edu

1

In the original definition of the PDO, the global monthly mean SST anomaly is removed prior to the calculation of the EOF. The leading EOF of this field and of the raw monthly SST anomalies is barely different, and we consider the raw SST anomalies. This is dictated by our desire to develop a regression model for the North Pacific SST anomalies, without attempting a similar model for the global mean of the SST.

2

Note that this is not strictly a prediction because the values of the forcing time series are assumed to be known up to time n.