Based on the Simple Ocean Data Assimilation (SODA) product and 37 models from phase 5 of the Coupled Model Intercomparison Project (CMIP5) database, the North Pacific Gyre Oscillation (NPGO) and its decadal generation mechanisms are evaluated by studying the second leading modes of North Pacific sea surface height (SSH) and sea level pressure (SLP) as well as their dynamical connections. It is found that 17 out of 37 models can well simulate the spatial pattern and decadal time scales (10–30 yr) of the NPGO mode, which resembles the observation-based SODA results. Dynamical connections between the oceanic mode (NPGO) and the atmospheric mode [North Pacific Oscillation (NPO)] are strongly evident in both SODA and the 17 models. In particular, about 30%–40% of the variance of the NPGO variability, which generally exhibits a preferred time scale, can be explained by the NPO variability, which has no preferred time scale in most models.
Two mechanisms of the decadal NPGO variability that had been proposed by previous studies are evaluated in SODA and the 17 models: 1) stochastic atmospheric forcing and oceanic spatial resonance and 2) low-frequency atmospheric teleconnections excited by the equatorial Pacific. Evaluation reveals that these two mechanisms are valid in SODA and two models (CNRM-CM5 and CNRM-CM5.2), whereas two models (CMCC-CM and CMCC-CMS) prefer the first mechanism and another two models (CMCC-CESM and IPSL-CM5B-LR) prefer the second mechanism. The other 11 models have no evident relations with the proposed two mechanisms, suggesting the need for a fundamental understanding of the decadal NPGO variability in the future.
Low-frequency variability of the Pacific climate system has been received much attention since the 1990s, owing to its strong influence on both regional and global climate (e.g., Latif and Barnett 1996; Biondi et al. 2001; Schneider et al. 2002). Understanding the dynamics of North Pacific decadal variability is vital for improving large-scale and regional climate predictability, as well as for being able to discern its potential global influence. Previous studies have shown how North Pacific variability is associated with the Pacific decadal oscillation (PDO), which represents the first leading mode of sea surface temperature (SST) in the North Pacific, and has wide-ranging connections with the Pacific climate as well as strong impacts on the marine ecosystem (e.g., Mantua et al. 1997; Latif 1998; Miller and Schneider 2000; Schneider et al. 2002; Kwon and Deser 2007; Newman et al. 2016). Like the PDO, the North Pacific Gyre Oscillation (NPGO), representing the second leading mode of sea surface height (SSH) anomalies in the northeast Pacific, is a basin-scale feature and captures prominent low-frequency changes in the Pacific physical and biological fields including temperature, salinity, sea level, nutrients, chlorophyll-a, and so on (Di Lorenzo et al. 2008, 2009). Moreover, it is an important part of the North Pacific gyre system and linked to the sequential evolution of several prominent patterns of basin-scale climate variability, such as the PDO and ENSO (Di Lorenzo et al. 2013).
Since the time when the NPGO was first identified by Di Lorenzo et al. (2008) for explaining key aspects of North Pacific climate variations, climatologists and oceanographers have proposed many hypotheses about what physical processes may be involved in the complicated system of Pacific dynamics. These new perspectives mainly concern what physical controls are most important in driving the NPGO (e.g., Chhak et al. 2009; Di Lorenzo et al. 2010, 2013, 2015; Yi et al. 2015) and to what extent the NPGO influences variations of the ocean, atmosphere, and marine ecosystem on a regional scale (e.g., Ceballos et al. 2009; Sydeman and Thompson 2010; Cloern et al. 2010; Sydeman et al. 2013).
Concerning the generation of the low-frequency variability of NPGO, a dynamical framework connecting the NPGO with the central Pacific ENSO (CP-ENSO) and the North Pacific Oscillation (NPO) has been proposed as a driving mechanism of the NPGO (Di Lorenzo et al. 2013, their Fig. 1). This mechanism involves low-frequency variability of the NPO being forced by the CP-ENSO (Ashok and Yamagata 2009), which then consequently drives the NPGO (Di Lorenzo et al. 2010; Furtado et al. 2011, 2012). In contrast, by analyzing reanalysis products Yi et al. (2015) suggested that the NPGO may originate from atmosphere stochastic noise, which has no preferred time scale, driving the oceanic integrated response that leads to a spatial resonance with a preferred oceanic time scale, thereby providing a key process in the generation of the low-frequency variability of NPGO. The dynamic framework proposed by Di Lorenzo et al. (2013), akin to an “atmospheric bridge,” was diagnosed by a set of idealized experiments with an atmosphere general circulation model coupled with a mixed layer (AGCM-ML), which suppresses the oceanic dynamic processes. The other mechanism proposed by Yi et al. (2015), however, emphasizes the adjustment of “oceanic integration,” in which ocean dynamics plays a critical role in anomalously subducting and peaking decadal signals of the NPGO in the upper 300 m of the North Pacific Ocean.
A complete understanding of the origins of NPGO low-frequency variability remains elusive. In particular, one of issues yet to be fully addressed involves the extent to which the NPGO and its driving mechanism behind the low-frequency variability are represented in the current generation of global climate models. Furtado et al. (2011) suggested that the simulated temporal and spatial statistics of the NPGO display distinct discrepancies among coupled general circulation models (CGCMs) participating in phase 3 of the Coupled Model Intercomparison Project (CMIP3). Moreover, CMIP3 models failed to reveal the atmospheric teleconnections emanating from the tropical Pacific Ocean to the North Pacific, namely the connection between the CP-ENSO and the NPO, which is critical to the operation of mechanism proposed by Di Lorenzo et al. (2013). At present, phase 5 of CMIP (CMIP5) now provides a new suite of historical simulations of various CGCMs with an overall improvement of model performance with regard to dynamic fields [e.g., SSH (Landerer et al. 2014); wind stress (Lee et al. 2013)], climate modes [e.g., ENSO (Kim and Yu 2012); the Pacific–North American teleconnection (PNA) and the NPO (Chen et al. 2018)], and climate drift (Sen Gupta et al. 2013), among others. Therefore, it is of great interest to evaluate how well the current generation of CGCMs participating in CMIP5 portrays the observed spatiotemporal features of the NPGO in the twentieth century, as well as to examine the range of model performances concerning the generation mechanism of decadal variability of the NPGO. For clarity, we define the decadal time scale here as 10 to 30 years; the multidecadal time scale refers to 30 years and longer.
The rest of the paper is organized as follows. Section 2 describes the datasets and methods used in this study. Section 3 shows the evaluation of the spatiotemporal features of the NPGO in CMIP5 models. The association between North Pacific atmospheric and oceanic variability is examined in section 4, followed by the physical mechanisms of decadal variability of NPGO in sections 5 and 6. The paper is concluded with a summary and further discussion in section 7.
2. Datasets and methods
a. Reanalysis data and multimodel outputs
The NPGO was examined based on the Simple Ocean Data Assimilation version 2.2.4 (SODA 2.2.4) reanalysis covering 1871–2008, which provides three-dimensional oceanic fields with horizontal resolution of 0.5° × 0.5° and 40 vertical levels (Carton and Giese 2008; Giese and Ray 2011; Ray and Giese 2012). Surface boundary conditions for SODA are taken from the Twentieth Century Reanalysis version 2 (20CRv2). In this study, monthly data of both SODA and 20CRv2 are used from 1900 to 2004 and then are remapped to the resolution of 1° × 1° and 2° × 2° for oceanic and atmospheric variables, respectively.
For model simulations, we investigated 37 climate models participating in CMIP5 [Table 1, organized by the Program for Climate Model Diagnosis and Intercomparison (PCMDI) for the Intergovernmental Panel on Climate Change’s Fifth Assessment Report (IPCC AR5; Taylor et al. 2012)]. These models are selected based on data availability for the diagnostics required [e.g., SSH, SST, sea level pressure (SLP), ocean temperature, wind stress, etc.]. We used multimodel outputs from the historical simulation, which incorporates the anthropogenic and natural forcings from the observed atmospheric composition changes in the twentieth century. In each model, outputs covering 1900–2004 were taken from only one ensemble member run. To facilitate model intercomparison, all fields are interpolated onto the same grids (i.e., 1° × 1° for oceanic variables and 2° × 2° for atmospheric variables). All analyses are based on monthly mean data.
b. Statistical methods
The NPGO is defined as the second-leading empirical orthogonal function (EOF-2) of SSH anomalies (SSHa) in the northeast Pacific (25°–62°N, 180°–250°E), and the normalized corresponding principal component (PC) time series is used as the NPGO index (Di Lorenzo et al. 2008). A positive NPGO phase is a dipole pattern with negative SSHa north of 40°N and the opposite south of 40°N. It should be noted that since the NPGO mode physically represents changes in the strength of the subpolar and subtropical gyres in the North Pacific, SSH is chosen to define the dynamic pattern in this study, which is different from the method of Furtado et al. (2011), who used SST as an indicator.
Here monthly anomalies are obtained by subtracting the climatological monthly means and removing the linear trend. All fields are weighted by the square root of the cosine of latitude to account for the decrease of grid area poleward. Extraction of the NPGO-like and NPO-like mode is achieved by applying SVD analysis to monthly SLP anomalies (SLPa) and SSHa that are smoothed by a 3-month running mean. Then the December–February (DJF) values of their fields are retained. The reason for using the winter data only is that North Pacific atmosphere–ocean connections are most prominent in the boreal winter months.
To investigate the temporal variability of the NPGO-like and NPO-like mode in SVD, we applied a multitaper spectrum analysis (Mann and Lees 1996) on the expansion coefficient (EC) time series. Following Furtado et al. (2011), we constructed the monthly EC time series by projecting the NPGO-like (NPO-like) pattern onto the unfiltered monthly SSHa (SLPa) fields. Spectral power is normalized to unity for an individual spectrum, which thereby denotes the contributing percentage of this frequency to total variance. A significant period is determined if the peak exceeds the 95% confidence level with respect to a red noise null hypothesis. The significance of correlation coefficients (regression coefficients) is assessed by a two-sided Student’s t test (F test) using the effective temporal degrees of freedom by considering serial autocorrelation.
To evaluate how closely the modeled pattern resembles the referenced one, we plotted a Taylor diagram (Taylor 2001), which is commonly utilized in assessing the relative merits of models. Herein, three statistics including spatial correlation, standard deviation, and root-mean-square (RMS) difference are calculated to quantify pattern similarity. The RMS can be used to evaluate whether two patterns have similar structure and amplitude of variation. In the Taylor diagram, simulated models that agree well with the referenced counterpart will lie nearest the reference point. Note that for the sake of direct comparison, RMS differences and standard deviations for each model are normalized by the referenced standard deviation.
The multimodel ensemble-mean statistics are calculated as the mean of the individual statistic being analyzed. For example, the ensemble-mean NPGO pattern is computed as follows: 1) the NPGO pattern is derived by regressing the monthly SSHa field onto the normalized NPGO index; 2) the sign of the NPGO mode for each model is kept in a positive phase; and 3) the individual NPGO pattern is averaged collectively.
3. Spatiotemporal features of the NPGO: Reanalysis data versus model simulations
First, we evaluate the capability of the 37 models to portray the spatial features of the NPGO in the twentieth century (1900–2004) by employing the Taylor diagram. As shown in Fig. 1, the normalized RMS errors vary greatly from 0.41 (CNRM-CM5) to 1.77 (GISS-E2-R), with 27 models showing values lower than 1. Note that around 90% of the models display a normalized standard deviation in amplitude of the NPGO pattern within the range of 1 ± 35% compared with the reference data. However, the spatial correlations range from 0.06 (GFDL-ESM2M) to 0.92 (CNRM-CM5), with only 51% of the models displaying spatial correlations higher than 0.7. Therefore, it is clear that the ability to simulate the spatial structure or the amplitude of the NPGO pattern varies greatly among the models. In the evaluation among 37 models, CNRM-CM5 agrees best with SODA data, with the highest correlation (0.92) and the lowest RMS error (0.41).
The multimodel ensemble-mean (MEM) of the NPGO pattern from the 37 models (the light purple dot in Fig. 1) very well captures the spatial structure of the NPGO pattern, with the spatial correlation reaching 0.89, but underestimates the amplitude deviation (standard deviation is 0.57). Specifically, as shown in Fig. 2, the MEM generally resembles the SODA values, with a negative SSH anomalous pole near the Alaska region (denoted as NPGOn) and a positive SSH anomalous pole over the subtropics (denoted as NPGOs). However, the simulated amplitude and location of the dipole differ somewhat from the referenced counterparts. The north (south) pole magnitude is much smaller than the SODA pattern of −0.03 m (0.025 m), by 47% (44%). Additionally, the MEM pattern exhibits a robust westward shift of its northern center and a slight eastward shift of the southern center.
Next, we compared the temporal variability of the NPGO between SODA and historical simulations by applying spectral analysis to the NPGO indices. Figure 3 shows that the SODA power spectrum exhibits significant power on decadal time scales (10–30 yr). Compared with SODA, the MEM power spectrum shows lower power in the period range of 7–20 years, whereas individual model generally exhibit higher power in a wider range for low frequencies. More specifically, 7 out of the 37 models (BCC-CSM1.1-M, GFDL-CM2.1, GFDL-ESM2M, HadGEM2-ES, INMCM4, MPI-ESM-LR and MRI-ESM1) do not show significant peaks on the observed decadal time scales. However, among the remaining 30 models, 25 models have peak power about 4.1%–150% higher than the reference, and 22 models show power at longer (multidecadal) time scales. On interannual time scales, over half of the models have significant peaks in the period range of 5–10 years, with mostly higher power than the (nonsignificant) reference.
Overall, results from the temporal and spatial statistics of the NPGO suggest that most models can reasonably capture the spatial pattern and decadal variability of NPGO. To further investigate the mechanism driving the decadal variability of NPGO, we selected 17 models (Table 2) that have significant decadal-period energy, and spatial correlations higher than 1 of MEM (0.66). These selected models are highlighted with red fonts in Fig. 3.
Figure 4 displays how well the selected 17 individual models simulate the NPGO pattern. It can clearly be seen that the basinwide pattern in SODA (Fig. 4a) is characterized by a south–north dipole of SSHa in the central and northeast Pacific, which is well reproduced by all the selected models (Figs. 4b–s), although most models exhibit notable differences. First, the variance contributions of the NPGO in individual models, ranging from 9.9%–16%, are all overestimated compared with the reanalysis counterpart (7.1%). Second, the selected model simulations generally underestimate the NPGO amplitude, as 16 out of 17 (10 out of 17) models have smaller central magnitudes of NPGOn (NPGOs) than the SODA pattern of −0.030 m (0.025 m) by 7.1%–51% (2.4%–34%). In addition, the center locations shift and the pole shapes distort, with 14 out of 17 and 10 out of 17 models tending to shift westward for NPGOn and NPGOs, respectively. The shifts are larger in the east–west direction than in the south–north direction.
Next, we inspected the vertical structure of the NPGO along meridional and zonal sections. Three sections, as indicated in Figs. 4a and 4b, are taken along 155°W (north to south), 48°N (east to west), and 29°N (east to west). Figure 5–7 show the meridional–zonal section of temperature anomalies associated with the NPGO. The climatological mean depth of mixed layer and specific isopycnal surfaces (26.2, 25.2, and 24) are given as a reference to generally localize the maximum signal of NPGO. Overall, the NPGO is a robust climate mode extending to 350 m deep and confined to the upper ocean, specifically above the 26.2 isopycnal surface. The depth of the strongest signal associated with the NPGO is in the interior ocean at around 100 m, which is close to the climatological mixed layer depth (Figs. 5a,b, 6a,b, and 7a,b), implying that the NPGO mode is closely related to active air–sea exchanges within the mixed layer.
The meridional (north–south) section in SODA (Fig. 5a) is characterized with stronger and shallower negative temperature anomalies in the northern Gulf of Alaska (GOA) and relatively weaker and deeper positive anomalies in the southern California Current System (CCS), which are well captured by all selected models (Figs. 5c–s). The MEM pattern of the 17 models (Fig. 5b) has a spatial correlation coefficient of 0.91 and an RMS difference of 0.4 compared with SODA. Correlation coefficients of almost all the models are larger than 0.65, except for CESM1-CAM5 (0.04). The maximum SST anomalies centered at the GOA and the CCS of Fig. 5 are evaluated in zonal sections, as shown in Figs. 6 and 7 respectively.
The north zonal [first east–west (E-W1)] section in SODA (Fig. 6a) exhibits consistent negative anomalies above the depth of 100 m stretching from the central Pacific to the California coast, while the south zonal [second east–west (E-W2)] section (Fig. 7a) displays positive anomalies at a deeper layer reaching 350 m at the central subtropics, and also negative anomalies east of 145°W. The structures of zonal section are generally found in most simulations, with spatial correlations at 0.86 (Fig. 6b) and 0.85 (Fig. 7b) for MEM models at two latitudes, respectively. However, some models fail to simulate the depth and location of each center, such as CESM1-CAM5 and NorESM1-M. The spatial correlations range from −0.01 to 0.87 (0.03 to 0.87), and only 7 out of 17 (4 out of 17) models display spatial correlations higher than 0.7 for the northern (southern) part.
In summary, SODA products and most of the 17 selected CMIP5 simulations suggest that the vertical structure of temperature variation associated with the NPGO signal is characterized by both a shallower anomalous signature in the northern GOA that extends to about 150 m deep and a deeper opposite-signed anomaly in the southern CCS that extends to about 350 m in the upper ocean.
4. Association between NPGO and NPO in the models
Previous studies have suggested that the NPGO is principally linked to the NPO, which is one of the prominent modes of planetary-scale atmospheric variability over the North Pacific basin (Di Lorenzo et al. 2008; Chhak et al. 2009; Ceballos et al. 2009). This section provides tests for the covariance of these patterns, by examining whether space–time coherencies are simulated, and to what degree ocean–atmosphere dynamical interaction is captured in the climate models.
a. The second leading modes in the North Pacific
The links were assessed by singular value decomposition (SVD) of the cross-covariance matrix between SSHa over the northeast Pacific (25°–62°N, 180°–250°E) and SLPa over the North Pacific basin (15°–70°N, 120°–270°E). Following Yi et al. (2015), analysis is limited by focusing on the second leading SVD mode (SVD-2), since the first leading mode usually relates to the PDO–Aleutian low pattern. Generally speaking, as we can see from Fig. 8, the second SVD spatial patterns in SODA correspond to NPO-like (Fig. 8a) and NPGO-like modes (Fig. 8b). The SVD-2 spatial patterns of MEM compare well visually with SODA, with correlations of 0.85 (0.96) for the SSHa (SLPa) field. Specifically, the spatial correlations of northern SSHa (SLPa) in each model range from 0.53 to 0.88 (0.8 to 0.97) as shown in the second (third) column of Table 2. However, there are two main differences between SODA and MEM: 1) the north region of the NPO-like pattern is visibly enhanced and the south region is weakened in MEM (comparing Figs. 8b and 8d) and 2) the strength of the NPGO-like pattern is weaker and its center is farther away from the North American coast in the models (Figs. 8a,c).
To further assess the temporal variability of the North Pacific atmosphere–ocean SVD modes, we applied spectral analysis to the EC time series. Before projecting the unfiltered monthly SSHa/SLPa field onto the NPGO-like and NPO-like patterns, these SVD spatial patterns are compared with the original definitions of the two patterns. In this study, the NPO is defined as EOF-2 of SLPa in boreal winter season from DJF over the North Pacific basin (20°–70°N, 120°–270°E), and the normalized corresponding PC time series is used as the winter NPO index (Linkin and Nigam 2008). For comparison, the corresponding winter NPGO pattern is calculated as the EOF-2 of SSHa using boreal winter (DJF) data over the northeast Pacific basin (25°–62°N, 180°–250°E), and the normalized corresponding PC time series is the winter NPGO index. As shown in Table 2 (Table 2; COR3 SSH and COR4 SLP), the spatial correlations of SLP field between SVD-2 and EOF-2 are high and commonly larger than that of SSH field in 17 models and SODA.
The power spectra of the EC times series corresponding to SVD-2 of SSHa (EC-SSHa or EC-NPGO) and SLPa (EC-SLPa or EC-NPO) are examined next. As shown in Fig. 9a, the spectral peaks of NPO in SODA and model simulations are not statistically significant at the 95% confidence level, except for seven models (CESM1-CAM5, CESM1-WACCM, CMCC-CM, CMCC-CMS, CNRM-CM5, HadGEM2-CC, and NorESM1-M). Note that only CESM1-CAM5 and NorESM1-M exhibit significant peaks for NPO indices on decadal time scales. It is indicative that the low frequency presents in the atmosphere through air–sea interaction, remote forcing, or random sampling error. On the other hand, the power spectra of the NPGO reflected in Fig. 9b have widely varying ranges of power on interannual to interdecadal time scales (4–70 yr), which is similar to the spectral results in Fig. 3 but with wider ranges. Most models, except for four (CCSM4, ESM1-FASTCHEM, CNRM-CM5, and HadGEM2-CC), display higher power than SODA, especially on the interdecadal time scales from 10 to 50 years.
It should be pointed out that the weak and insignificant white noise of atmospheric forcing may act at low frequencies to redden the forced oceanic response. The power spectra rates of EC-SSHa to EC-SLPa are calculated to assess whether there are relatively nondominant peaks in the oceanic NPGO related to the atmospheric NPO. As shown in Fig. 10, SODA displays a robust peak at 5 years and some enhanced variance at decadal time scales (~18 yr). The peak at 18 years is a weak one, suggesting that the random atmospheric variability may contribute to the SSH variability. Compared with the referenced data, the MEM spectrum ratio exhibits a wide but unobvious peak around 15–50 years, while the individual spectrum ratio generally shows a significant peak at low frequencies. The majority of CMIP5 models (11 out of 17) have obvious peaks at decadal time scale, such as CESM1-CAM5 and CNRM-CM5, revealing that the insignificant decadal spectra of atmospheric forcing have greatly limited contributions to the corresponding decadal time scales of this oceanic mode in this relatively short record.
b. Reconstruction of the NPGO index using an AR-1 model
The dynamical ocean–atmosphere link between the North Pacific SSH and SLP is reflected by the percentage of squared covariance (COV) explained by the SVD-2 (Table 2). In SODA, the NPGO–NPO patterns are important modes that link the North Pacific SSH field and SLP field, which explain 11% of total covariance. Although the model data exhibit various covariance contributions ranging from 2.9% (CMCC-CESM) to 21% (CNRM-CM5), the MEM (with a value of 14%) is higher than that in SODA, suggesting the association of the second leading modes strengthened in the CMIP5 simulations.
The strength of the association is further measured in a simple model following Chhak et al. (2009) and Furtado et al. (2011). The forced ocean variability in this simple model is computed by a first-order autoregressive process (AR-1), which approximately decays as an exponential function in SODA and models (not shown).
In Eq. (1), the reconstruction of the NPGO indices, , is driven by atmospheric variability and uncorrelated white noise . Here, the EC times series corresponding to SVD-2 of SSH and SLP are used to reconstruct the . The coefficients are obtained using least squares fitting method.
As shown in Fig. 11, the reconstructed result for EC-NPGO in SODA indicates that the NPO variability contributes about 40% of the variance of the EC-NPGO (Fig. 11a, red line). For model simulations (Fig. 11b, red lines), the reproduced results are consistent with the SODA association, as the temporal correlations in all models are fluctuating at around 0.6, accounting for about 36% of the variance. These results are also supported by another definition of atmospheric forcing , defined as the region-averaged SLP anomalies difference between two nodes of NPO (south node minus north node; black box in Figs. 8b,d). This forcing index is similar to EC-NPO, with temporal correlation at 0.92 for the MEM. The good reconstruction of NPGO using this forcing index is reflected by SODA with correlation at 0.58 (Fig. 11a, blue line) and supported by the MEM (correlation is 0.60) and all the models (Fig. 11b, blue line). The coefficient from these reconstructions indicates that the damping time scale of NPGO (Furtado et al. 2011) in the MEM of CMIP5 models is about 4.1 months, which is slightly longer than that in SODA (3.8 months). The coefficients and are not shown here.
In summary, results from the SVD analysis and the AR-1 simple model suggest that the North Pacific second leading climate mode is well simulated in 17 CMIP5 models. The oceanic signal/NPGO with significant peak at interannual and decadal time sales is directly influenced by the atmospheric forcing/NPO, which is relatively weak in the low-frequency band and generally regarded as white noise. About 30%–40% of the total variance of the NPGO variability can be explained using a simple linear model forced by the NPO variability.
5. Driving mechanism of the NPGO decadal variability: Stochastic spatial resonance
In this section, we test whether the “stochastic spatial resonance” mechanism proposed by Yi et al. (2015) functions in the selected CMIP5 models. As a mechanism for turning white-noise atmospheric forcing with a preferred spatial structure (NPO) into a preferred time scale in the oceanic response (NPGO), spatial resonance requires three conditions: 1) stochastic forcing, 2) propagation of the forced oceanic signals (e.g., induced by the advection or planetary wave) along the positive/negative antipodal forcing structures of the atmospheric forcing (here, in the meridional direction), and 3) phase reversal matching the time scale of propagation.
a. Life cycle of the decadal oscillation
As mentioned in section 4, the atmospheric forcing indicates relatively stochastic variability, and potentially has a strong impact on the North Pacific. According to Yi et al. (2015), the upper ocean SSHa associated with decadal NPGO mainly arise from the wind-induced Ekman pumping and meridional advection. During the NPGO positive phase, the anomalous cyclonic (anticyclonic) wind over the subpolar (subtropical) region can apparently induce an anomalous upwelling (downwelling) to support negative (positive) anomalies in that area. Meanwhile, the anomalous westerly wind located around 40°N can create equatorward Ekman transport, which also favors cooling the midlatitude ocean. With the contribution from ocean dynamics (e.g., advection, mixing) and surface heat flux, the corresponding SST and heat content (HC) vary with the SSH changing.
Following Yi et al. (2015), the heat budget of the anomalous surface temperature associated with the NPGO is diagnosed in CMIP5 models. Our results (not shown) highlight the dominant effects of temperature advection by the anomalous Ekman pumping and the anomalous meridional Ekman current, along with the inconsistent roles of surface net heat flux in the selected CMIP5 models. The result comes from four CMIP5 models (CMCC-CM, CMCC-CMS, CNRM-CM5, and CNRM-CM5.2) that exhibit good simulations of the decadal variability of NPGO.
The evolution of the NPGO decadal oscillation, which has been discussed in SODA by Yi et al. (2015), is next inspected in CMIP5 models. They suggested that the southwestward propagation of ocean heat content (OHC) anomalies in conjunction with the enhancement of local wind stress curl (WSC) forcing, associated with the meridional dipole pattern of the NPO stochastic forcing, could provide a potential mechanism of spatial resonance (Jin 1997; Saravanan and McWilliams 1998) for decadal NPGO. The similar life cycles of NPGO are examined in the selected four CMIP5 models by regressing NPGO index against multiple indicators, including OHC, SSH, WSC, SST, surface net heat flux, and SLP, which are dynamically connected. The OHC is integrated above the climatological 26.2 isopycnal given that the NPGO variability is mainly confined above this potential density surface (Yi et al. 2015; also evident in Figs. 5–7). Here all variables are bandpass filtered to retain decadal variability at 10- to 30-yr periodicity. Given each model has an individual peak at a decadal time scale, exhibiting the coherent phase evolution from multimodel data, the relative lags at 0%, 25%, 50%, 75%, and 100% of the half cycle are selected based on dominant period per model.
As shown in Fig. 12, compared with SODA (Yi et al. 2015, their Fig. 7), it is obvious that the MEM regression map generally resembles the observational one, displaying a clear southwestward movement of OHC in the northeastern Pacific (Fig. 12, the first column), the matched SST anomalies changing on the surface, and the changes in the overlying wind stress and the related SSH. The net heat flux anomalies of the MEM simulation are inconsistent with the SODA results for the model diversity mentioned above and so must be discussed individually.
To evaluate the model differences in how surface net heat flux affects the decadal evolution of NPGO, we performed lead–lag regressions of heat flux and the influenced SST anomalies against the NPGO index over two key regions (Fig. 15, solid black box), the NPGOn and NPGOs regions, located at the cores of the dipole structure. All the variables are averaged in these regions, and are then regressed against the band-passed NPGO index. Figure 13 shows the lagged regressions for SODA and individual models. The results suggest the influence of heat flux anomalies occurs in three ways. Both SODA and CMCC-CMS exhibit forcing effects in the north but damping effects in the south, although there is a 2–3-yr lead for the heat flux in SODA. In contrast, the heat flux behaviors of CNRM-CM5 and CNRM-CM5.2 involve damping the SST anomalies over the NPGOn region and forcing the surface changes over the NPGOs region. In the CMCC-CM, surface heat flux damps SST variations over most of the basin for decades.
In summary, in SODA and the selected four CMIP5 models, the decadal SST variations associated with the NPGO are generally dominated by similar surface winds anomalies but different surface thermal processes. With the oceanic signals enhanced by the atmospheric forcing when the phase reverses, a decadal evolution of NPGO is established.
b. Propagation of oceanic signals
Given that understanding, we further examined the NPGO transit by investigating its vertical structure to figure out whether there is robust southward propagation in both SODA and historical simulations. Considering the evolution of the OHC anomalies (e.g., Fig. 12, the first column), a representative meridional band is selected in the eastern North Pacific (130°–145°W; Fig. 15, dashed black box). As shown in Fig. 14, we compared the decadal evolution of subsurface temperature anomalies associated with the NPGO in SODA and the selected four models. The relative lags at 0%, 25%, 50%, 75%, and 100% of the half cycle per model are performed to keep phase consistency. It can clearly be seen that the MEM (Fig. 14, the second column) shows a remarkable meridional propagation of oceanic signals from the midlatitudes to low latitudes along isopycnals, which is similar to SODA [Fig. 14, the first column; similar to Fig. 8 in Yi et al. (2015) but in different lags]. The signals subduct around 45°N and then penetrate to about 20°N, mainly confined above the 26.2 isopycnal, although there are deeper isopycnals in the models. For the individual models, the southward propagation tendency is consistent with SODA in some models such as CMCC-CMS.
Furthermore, from SODA and the selected CMIP5 models, the preferred time scale for the resonance mechanism is evaluated using the advection speed of within subsurface layers. The zonal mean flow velocity is calculated within a meridional band (130°–145°W; Figs. 15a,b), which is regarded as the strongest response of ocean–atmosphere interaction of the NPGO. As shown in Fig. 15c, considering the southward velocity and basin scale, roughly calculating from the latitude where velocity is zero to latitude at 20°N, it takes about 14, 14, 12, 15, and 13 years for SODA, CMCC-CM, CMCC-CMS, CNRM-CM5, and CNRM-CM5.2, respectively. The estimated time scales basically match the decadal time scales, but not the dominant peak for each model.
It should be noted that our evaluation of the decadal mechanism of NPGO is crude for two main reasons: 1) all CMIP5 models are 10- to 30-yr bandpass filtered to highlight the physical mechanism of decadal variability but sacrifice the focus on the specific decadal band in each model; and 2) the decadal time scales of the NPGO cycle could not be computed accurately given the uncertainty of the meridional pathway for each model. Nevertheless, it appears that the spatial resonance mechanism for the generation of the NPGO decadal variability with a preferred time scale is active in the selected four CMIP5 models.
6. Driving mechanism of the NPGO decadal variability: Tropical forcing
The previous two sections evaluated the extratropical ocean–atmosphere dynamical relations and the local physical mechanism for driving regional NPGO decadal variability. In this section the influence of atmospheric teleconnections from the tropical Pacific on extratropical NPGO variability is evaluated. This represents another physical mechanism for the origin of decadal variability of the NPGO [proposed by Di Lorenzo et al. (2010)], with far-reaching effects on the extratropics including ocean–atmosphere interaction.
a. The second leading mode of tropical representation on North Pacific
Previous studies have shown that the NPGO can be affected by El Niño Modoki (Ashok et al. 2007) or central Pacific warming (CPW; Di Lorenzo et al. 2010) via atmospheric teleconnection acting on the NPO mode. Following Furtado et al. 2011, the CPW mode is defined as the EOF-2 of SST anomalies in the boreal winter season (DJF) over the Pacific tropical region (−20°S to 20°N, 110°–290°E), and the related PC and EC (calculated like the EC-NPO in section 4a) time series are used as the winter CPW index and EC-CPW index, respectively. As shown in reanalysis data (Fig. 16a), a negative phase of the CPW pattern is characterized as a maximum anomalous cooling centered in the central equatorial Pacific and a smaller anomalous warming located in the eastern and western Pacific. The MEM of the CPW pattern from the 17 models basically resembles that of SODA, with the spatial correlation reaching 0.77, but it fails to simulate the warming of the western Pacific. Note that the spatial correlations vary considerably among individual models, ranging from 0.12 (CanESM2) to 0.78 (CNRM-CM5.2). For some models (9 out of 17 models; e.g., NorESM1-ME), the CPW mode is not tracked as the second but rather as the first leading mode in the tropics (ENSO), which echoes the similar results in a previous study of CMIP3 (Furtado et al. 2011).
Next, the CPW mode is linked to North Pacific atmosphere by displaying the correlation maps between the SLPa in boreal winter (DJF) season and the winter CPW index. The black dots denote temporal correlations exceed the 95% significance level. For the 20CRv2, the SLP pattern (Fig. 17a) displays a northeast–southwest dipole of with negative and positive centers over the Alaska and the Hawaii regions, respectively, resembling the corresponding NPO-likewise mode (Fig. 18a; extracted by correlating the winter SLP time series of the SVD-2 with the Pacific-wide SLPa). The spatial correlation between the SLPa pattern related to the CPW (Fig. 17a) and the NPO-likewise pattern (Fig. 18a) is small but significant (0.33; >95%). The same correlation maps are presented for the 17 models. As shown in Fig. 17b, the MEM displays similar structure to its observational counterpart, with spatial correlation of 0.67. However, there is substantial difference in the SLPa pattern for individual models (Figs. 17c–s), since only five models have correlation exceeding 0.5.
Taken together, most models fail to simulate the links for the CPW-NPO in the North Pacific, which undermines the idea that the modeled CPW variability has important connections with the North Pacific oceanic signal via atmospheric bridge at the subtropics.
b. North Pacific air–sea interaction induced by the tropical forcing
With tropical variability connected with the North Pacific via the subtropical atmospheric bridge, the degree to which the NPGO variance can be explained is assessed in some models. To focus on the CPW forcing, only four models are selected: CMCC-CESM, CNRM-CM5, CNRM-CM5.2, and IPSL-CM5B-LR. They exhibit good simulations of the CPW features compared to SODA, with the spatial correlation coefficients larger than 0.5, both in the basin-scale modes (Fig. 16) and the correlation maps (Fig. 17).
Figures 18a–d exhibit the key region of the atmospheric bridge that links the North Pacific air–sea interaction to the tropical variability for SODA and the MEM of the selected four models. By correlating the Pacific-wide SLPa with the winter SLP time series of the SVD-2, the atmospheric expressions of the NPGO are shown in Figs. 18a and 18b. The atmospheric patterns of the CPW are also displayed in Figs. 18c and 18d using a similar method. As shown in Figs. 18a–d, the results indicate that the region with black rectangle in the proximity of Hawaii (10°–30°N, 191°–231°E) shares the same SLP anomalous variation. Therefore, the SLPa averaged over the black box (Figs. 18a–d) are used to drive the AR-1 model, which was introduced in section 4b. The reconstruction of the NPGO index from SODA (red line, Fig. 18e) is significantly correlated with the original NPGO index (black line, Fig. 18e; EC-NPGO time series, from SVD analysis in section 4b) at 0.46. That is, it can explain 21% of the variance of the NPGO variability. Similarly, four models and their ensemble-mean model display comparable contributions to the NPGO variance shown in Fig. 18f, with the maximum value for CESM-CAM52 (31%) and the minimum for IPSL-CM5B-LR (12%).
As shown in Fig. 18f, the reconstruction of the NPGO index does not exhibit a large enhancement in the percentage of NPGO variance explained for the low-frequency variability in SODA (temporal correlations are 0.46 and 0.40 for the monthly and 9-yr low-pass index, respectively) and the MEM (0.48 vs 0.57). Only one model (CNRM-CM5) shows that the contribution to the NPGO variance is greatly reinforced. These results suggest that the hypothesis [proposed by Di Lorenzo et al. (2010), their Fig. 2] that low-frequency variation of SLP near the Hawaii region plays an important role in the decadal variability of NPGO is well captured in some models.
In summary, most CMIP5 models fail to show the observed CPW–NPO–NPGO dynamical connection. Only four of the models and SODA exhibit indications for the low-frequency part of the NPGO, which originates from couple ocean–atmosphere interactions associated with the CPW dynamics.
7. Summary and discussion
Climate modes such as ENSO, the PDO, and the NPGO strongly impact ocean dynamics and marine ecosystems, and are linked to other extreme weather events on a global scale, particularly over North America and East Asia. In this study, common NPGO metrics are used to compare coupled model behavior to SODA, including the frequency/amplitude of SSHa in the eastern North Pacific and the dynamical links to the atmosphere, as well as the horizontal and vertical structure of temperature anomalies associated with the NPGO. Additionally, we complete those evaluations by exploring the model performance in terms of the generation mechanism of the decadal variability of NPGO.
Through analyses of the spatial and temporal statistics of NPGO using Taylor diagrams and power spectral methods, we provide evaluations for the simulation of the NPGO mode in the 37 CMIP5 models and suggest that there are remarkable differences across the various models, with varying degrees of resemblance to the observational characteristics of the NPGO. When evaluating the power spectra of the normalized NPGO indices, each model shows stronger decadal power than SODA and a wide range of behavior at the decadal time scales (Fig. 3). The results of the SVD analyses reveal dynamically consistent dipole structures of the SLPa and the SSHa. Taken together by comparing the NPGO and NPO patterns (Fig. 8), the CMIP5 CGCMs tend to display consistent pattern features, indicating a plausible explanation in generation of NPGO from the atmospheric forcing (NPO). The oceanic signal (NPGO), with a significant peak at interannual and decadal time scales, is directly forced by the atmospheric mode/NPO, which is relatively weak in low-frequency band and generally regarded as white noise. About 40% of the variance of the NPGO variability can be explained using a linear AR-1 model response to this forcing.
Although the structural similarity of basin-scale patterns and vertical sections are largely captured in these historical simulations (17 out of 37), we selected only a few models for the evaluation of the two proposed mechanisms of the decadal NPGO variability: 1) stochastic atmospheric forcing and oceanic spatial resonance, and 2) low-frequency variability of NPO forced by the CPW.
The physical process of the decadal NPGO in observation and in four selected CMIP5 models displays good similarity and indicates that a stochastic spatial resonance mechanism might be active, as previously demonstrated by Yi et al. (2015) for SODA. The southwestward propagation of oceanic signatures in conjunction with enhancing effects of local WSC associated with the meridional dipole pattern of the NPO forcing are noted to provide the necessary components of the mechanism of spatial resonance.
On the other hand, there is also an indication of the CPW–NPO–NPGO physical mechanism (Di Lorenzo et al. 2010) being active in the CMIP5 models. The dynamical connections between the ocean–atmosphere interactions in North Pacific and the CPW mode are well simulated in the SODA ocean analysis product and also in the four selected models, which contribute about 20% of the NPGO variance.
Di Lorenzo et al. (2015) proposed a Pacific climate null-hypothesis framework, in which the Pacific meridional mode (PMM; Chiang and Vimont 2004) is driven by subtropical stochastic noise associated with the NPO so that the tropical climate system is thereby charged with decadal scale variability. In the tropics, this low-frequency variability is enhanced by local ocean–atmosphere feedbacks and then affects the extratropics through atmospheric teleconnections. Our results suggest that the decadal-scale energy of the NPGO can be enhanced both by the oceanic spatial resonance of local stochastic atmospheric forcing and through the low-frequency teleconnections excited by the equatorial Pacific. Evaluation of 17 models with reasonable simulations of NPGO decadal variability reveals that two models (CNRM-CM5 and CNRM-CM5.2) can be explained both by two mechanisms, and two models (CMCC-CM and CMCC-CMS) tend to display an oceanic spatial resonance mechanism, while two other models (CMCC-CESM and IPSL-CM5B-LR) are more likely forced by the tropical SST through the teleconnection.
It should be noted that the subsurface thermal anomalies associated with the NPGO signal and its southward propagation can be studied further. According to Schneider et al. (1999), this subsurface temperature response might be due to the heaving of thermocline or/and the spiciness variability. However, the southward-propagating signal cannot be found in the lagged regression of salinity anomalies in SODA and the most CMIP5 models (figures are not shown). Therefore, we propose a hypothesis that for some CMIP5 models (e.g., CMCC-CMS and CNRM5.2), the southward propagation result from both the heaving of thermocline and spiciness variability. For SODA and some models, compared with the spiciness variability, other possibilities might play a decisive role, such as the planetary wave in a higher baroclinic mode [e.g., A mode (advective mode); Liu 1999].
The main conclusions of our work focus on the evaluation of the NPGO in CMIP5 models and the dynamical consistency of the climate mode variability with SODA in the context of the generation mechanisms for NPGO decadal energy. Although this study considered multiple models, the runs were confined to the 105-yr historical period, so longer control simulations or ensembles of historical runs could have increased statistical reliability and potentially decreased the noise that obscures the dynamical signal. In addition, targeted coupled ocean–atmosphere modeling experiments are needed to separate the fraction of decadal variability that is generated intrinsically in the extratropics versus in the tropics.
We acknowledge the CGCM modeling groups, the PCMDI, and the WCRP’s Working Group on Coupled Modeling for their endeavor in making available the WCRP CMIP5 datasets. This work is supported by National Natural Science Foundation of China (NSFC) Projects (41490643, 41490640, 41506009, 41521091, U1606402). AJM was supported by the NSF Earth System Modeling Program (OCE1419306), the NOAA Climate Variability and Prediction program (NA14OAR4310276), and the NOAA Modeling, Analysis, Predictions, and Projections program (NA17OAR4310106). Discussions with Prof. Niklas Schneider were greatly appreciated. We thank the three referees for extensive comments and suggestions that greatly improved the manuscript. The authors acknowledge financial support from the China Scholarship Council.