The influence of Pacific-North American teleconnection on the North Pacific SST anomalies in Wintertime under the global warming

The impact of the enhancing Pacific-North American teleconnection (PNA) on the intensity of the first mode of SST anomalies (SSTa) in wintertime North Pacific [known as the Pacific decadal oscillation (PDO)-like SST pattern] under global warming is studied using reanalysis datasets and 12 selected CMIP5 models. The robust observational result of the PNA one month ahead of the PDO shows that the PNA has a great effect on the PDO-like SSTa. The intensities of PNA and PDO are defined with elimination of the ENSO signals to evaluate the direct impact of PNA on the PDO in the North Pacific under global warming. By comparing RCP8.5 and historical scenario of selected multi-models (PNA enhanced models), future projection illustrates that the PDO intensity will intensify 2.2 times its internal variability. Approximately 67% of the increase in PDO variations is contributed by the PNA-induced SSTa in the future warming scenario of multi-models ensemble mean. Models with stronger intensification of PNA variability tends to have larger magnitude of SSTa response in the North Pacific, which has more contributions to the enhancement of PDO intensity change in a warmer climate. This may shed some light on the projection of PDO variability and the relative role of PNA forcing under the global warming.


Introduction
The dominant pattern of observed sea surface temperature (SST) variability in the wintertime North Pacific is referred to as the Pacific decadal oscillation (PDO), which has been first introduced by Mantua et al. (1997) as the leading empirical orthogonal function (EOF) of North Pacific (20° N-70° N) SST anomalies (SSTa). The PDO is characterized with the time scales in the observed records of 15-25 years and 50-70 years (Minobe 1997). The positive PDO pattern is featured by a "horseshoe" shape, with below-normal SSTa in the central North Pacific and the opposite SSTa along the west coast of North American and the eastern tropical Pacific. The PDO has great climatic influences on the Pacific Rim. Positive PDO phases are associated with less precipitation and warmer air temperature anomalies in the northwestern and central America, and more precipitation and cooler air temperature anomalies in the southwestern America (Mantua and Hare 2002). The negative phase of the PDO is exactly opposite influences to the positive phase. Such influences have drawn attentions on its changes under the future global warming.
Subsequent research, has found that the PDO is not a single physical mode of climate variability but instead largely represents the combination of different physical processes, including ocean Rossby waves in the ocean gyre circulation, teleconnections from the tropical Pacific, and internal atmospheric variability forcing, which contribute to the PDO variability on a variety of time scales and regions in the North Pacific (Iwasaka and Wallace 1995;Qiu and Chen 2005;Alexander and Scott 2008;Newman et al. 2016). A body of research exists focusing on the decadal time scale of the PDO, and it is mainly provided by the westward ocean Rossby wave propagation along the Kuroshio-Oyashio Extension (KOE) region, which is the most important contributor to the decadal component of PDO, although the relative importance of propagation is different among models (Zhang and Delworth 2015). The meridional displacement of the subtropical and subpolar gyre boundary generates the KOE SST anomalies. The positive SSTa can feed back to the overlying atmosphere, resulting above-normal sea level pressure anomalies that can further enhance SSTa via Ekman transport. This positive air-sea feedback over the North Pacific makes great contributions to the decadal time scale and it is in agreement with observations (e.g., Seager et al. 2001;Nakamura and Kazmin 2003;Frankignoul and Sennéchael 2007;Zhang et al. 2011).
The PDO is also associated with variability outside the North Pacific (Deser et al. 2004). Teleconnections between the tropics and extratropics have an integrating effect on the PDO variability (Gu and Philander 1997;Schneider et al. 1999). For the sensitivity of North Pacific SSTa to highfrequency atmospheric variability (Pierce 2001;Newman et al. 2003), the Pacific-North America (PNA) teleconnection which is induced by the El Niño-Southern Oscillation (ENSO) makes contributions to the North Pacific SSTa via ocean-atmosphere interactions Barnett 1994, 1996). During the ENSO events, tropical Pacific SSTa induces tropical precipitation shifts forcing PNA wave trains (e.g., Trenberth et al. 1998;Liu and Alexander 2007). This ENSO induced teleconnection causes changes in SST and ocean currents by altering air temperature, wind, humidity, and clouds in the North Pacific (e.g., Lau and Nath 2001;Alexander et al. 2002).
In addition, the PDO variability is response to internal atmospheric variability forcing, which is associated with the white noise forcing that generates much of the SST variability over the North Pacific (Frankignoul and Reynolds 1983). This dominant atmospheric circulation pattern is closely linked to interannual variability in surface fluxes and SSTa which modulates the PDO variability (Cayan 1992;Iwasaka and Wallace 1995). Note that, in an atmospheric general circulation model (AGCM) coupled to a mixed layer ocean model with no currents (no ENSO variability or ocean dynamics), the dominant atmospheric pattern is associated with fluctuations in the strength of the Aleutian low pressure system primarily resulting from internal atmospheric dynamics, including large-scale dominant atmospheric teleconnection pattern such as the PNA (Pierce 2001;Alexander et al. 2010).
The PNA is a prominent mode of planetary-scale atmospheric variability over the North Pacific basin in boreal winter, which has powerful influences on the hydroclimate of the North Pacific rim. These dynamical changes are expected to impact the surface climate of North America especially in winter (e.g., Wallace and Gutzler 1981;Leathers et al. 1991). Furthermore, the PNA's signature on the Aleutian Low is closely related to PDO (e.g., Latif and Barnett 1996;Deser et al. 2004). When El Niño events peak during boreal winter, the North Pacific center of PNA deepens and results a positive PDO phase by causing variations of evaporation, wind stress and air temperature in the North Pacific (Strong and Magnusdottir 2009). In addition, the PNA plays an active role in connecting the climate variations between the tropical and extratropical Pacific. Such as the ENSO variability influences surface climate over the North Pacific through atmospheric PNA-like changes (Alexander et al. 2002). Furthermore, PNA also have a close relationship of the atmospheric variability between midlatitude and polar climate. For example, Arctic Oscillation is significantly related to the geopotential height (GPH) anomalies around the signature on the North Pacific center of PNA (Wang et al. 2005;Leung et al. 2017Leung et al. , 2022aJian et al. 2020).
Based on the multi-model outputs of CMIP5, the capability of climate models to simulate the spatiotemporal features of the PNA in the twentieth century has been evaluated. The assessment reveals that 14 out of 30 models can well capture the geographical distribution and amplitude of the PNA pattern. For the quadrupole structure of the PNA, terminologies have been named after PNAa (the center near Hawaii), PNAb (the North Pacific center), PNAc (the center over Alberta, Canada) and PNAd (the center in the Gulf Coast region) by the propagation direction of PNA wave train. It is manifested that the change of PNA intensity is mainly caused by the PNAb. And there have been 12 out of 14 selected models show the enhancement of PNAb under global warming (Chen et al. 2018). In this case, how the SSTa of the North Pacific will response to the enhanced PNAb in the 12 CMIP5 models (see Sect. 2.1 for details) is our interest.
The purpose of this study is to examine the changes in the PDO's spatial structure under global warming and the consequent changes in relationship with the enhancing PNA intensity in future projection. We also briefly discuss the influences of PNA on the PDO and the possible dynamical processes that contribute to the PDO intensity change under global warming. The paper is organized as follows. Section 2 describes the observational datasets, multi-model outputs and analysis methods. The 21st-century projected changes of PDO and the relative role of the enhancing PNA to drive significant PDO intensity change under global warming are shown in Sect. 3. Summary and discussion is given in Sect. 4.

Observations and multi-model outputs
To investigate the observed PDO and the oceanic anomalies forced by PNA, we employ SST data taken from the Met Office Hadley Centre Sea Ice and Sea Surface Temperature version 1 (HadISST1) dataset, with a horizontal resolution of 1° × 1° (Rayner et al. 2003). This dataset has been checked against other reanalysis datasets (such as COBE, ERSST) and found relatively minor differences both in time series and pattern of PDO from the years 1901-2014 (Newman et al. 2016). The PNA index and the spatial pattern in observation are investigated based on GPH data from the Twentieth-Century Reanalysis dataset version 2c (20CRv2c), with a horizontal resolution of 2° × 2° and 24 vertical pressure levels (Compo et al. 2011).
To explore the North Pacific SSTa response to the intensified PNAb under global warming, we select 12 CMIP5 models (Table 1) which robustly demonstrate the intensification of PNAb (Chen et al. 2018), organized by the Program for Climate Model Diagnosis and Intercomparison for IPCC AR5 (Taylor et al. 2012). The multi-model outputs of CMIP5, including monthly GPH and SST are used for analysis. The pre-industrial (PI) control simulation represents the unforced natural variability of ocean-climate system with CO 2 concentration fixed at 280-290 ppm; the historical simulation indicates the observed oceanic and atmospheric composition changes in the twentieth century; the future climate simulation under the Representative Concentration Pathways (RCP) 8.5 shows the worst global warming scenario, with the radiation forcing reaching approximately 8.5 W m −2 in 2100 (equivalent to CO 2 concentration > 1370 ppm). To ensure the outputs are comparable, we choose 95 years of data in both historical  and RCP8.5 (2006RCP8.5 ( -2100 scenario to represent the past and future climate. Only one member (''r1i1p1'') run for each model is used in this study.
In Sect. 3, we use SST and GPH from the historical and RCP8.5 scenario of the 12 selected models to study the changes of PDO under global warming. Then we capture the internal variability of PDO in PI run to determine if the changes are significant. These outputs are also employed to further investigate the projected future SSTa response to the enhanced PNA forcing and its contribution to the changes of PDO.

PDO, PNA and its forced SSTa patterns along with their intensity
Previous studies have defined the PDO as the leading empirical orthogonal function (EOF) of the monthly averaged SSTa of the North Pacific (e.g., Mantua et al. 1997;Newman et al. 2016). In this study, the observed PDO pattern is extracted by regressing SSTa onto the normalized leading principal component (PC) of the North Pacific (20° N-70° N), which represents the typical amplitude of PDO (i.e., one standard deviation of the index). Note that the SST field is weighted by the square root of the cosine of latitude to account for the decrease of grid area poleward. In order to measure the intensity of PDO, each PDO pattern has been adjusted to the positive phase, and the central North Pacific domain for the SSTa is set as a box (11° in latitude and 21° in longitude) on the below-normal SSTa center of PDO (as shown in Fig. 2a).
The absolute area-averaged value of regression coefficients in the box is defined as the PDO intensity (− PDO i ).
As to the PNA pattern, monthly December-January-February (DJF) 500-hPa GPHa are used for the EOF and regression analyses. Monthly anomalies are first obtained by subtracting the climatological monthly means and removing the linear trend, and then they are smoothed by 3-month running mean to extract the DJF 500-hPa GPHa. The PNA pattern is derived by the regression of the GPHa field onto the normalized PC1 of 500-hPa GPHa. The PNA intensity is calculated as a linear combination of four centers intensities of the PNA quadrupole (i.e. PNAa-PNAd) in its positive phase. The intensity of the PNA quadrupole is computed as the absolute area-averaged value of regression coefficients in the central area of the centers. Only PNAb intensity (− PNAb i ) is counted in the present study for it is manifested that the change of PNA intensity is mainly caused by the PNAb (Chen et al. 2018). It is also well known that the North Pacific atmosphere-ocean system fluctuates with periods of 2-7 years in response to El Niño-Southern Oscillation (ENSO) with the narrow cold tongue located in the eastern equatorial Pacific between 28° N and 28° S (e.g., Bjerknes 1969;Horel and Wallace 1981;Rasmusson and Wallace 1983;Halpert and Ropelewski1992;Deser et al. 2004). And the associated PNA teleconnection causes changes in the wind-driven mixing, heat fluxes, and Ekman transport in the upper North Pacific Ocean, which all acts to significantly influence the PDO pattern (Alexander et al. 2002;Alexander and Scott 2008). Consider this circumstance, we use partial least squares (PLS) regression method to eliminate the ENSO effects and focus on the SSTa forced by PNA in the central North Pacific. The PNA-forced SSTa pattern is calculated by PLS regressing SSTa onto the normalized PNA index and removing Niño 3.4 index. As the same with the PDO intensity, the PNA-forced SSTa intensity is computed by absolute area-averaged value of the identical domain (i.e., magenta boxes in Fig. 2e).
Note that the PDO, PNA and PNA-forced SSTa pattern are calculated for the historical and RCP8.5 scenario separately, the difference in individual intensity indicates the response to greenhouse warming. Meanwhile the multimodel ensemble mean (MEM) statistics are calculated as the mean of the individual statistic being analyzed. For example, the MEM pattern of PDO is computed as the mean of individual patterns in positive phase which are derived by the regression of SSTa field onto the normalized PDO index from individual models.

Estimation of unforced internal variability in CMIP5 models
To test the significance of the projected PDO variability changes under greenhouse warming, we calculate the unforced internal variability of both modes based on the PI run of each model. First, we construct 41 segments with time increment of 95-year, for which the starting year of each one was taken 10-year apart from the 500-year long (after spin-up) outputs of the PI run. These 41 segments are considered independent of each other under the assumption that atmospheric states among different segments are uncorrelated. Then, a total of 41 PDO patterns are derived from 41 segments. Indeed, further inspection finds no significant correlations of PDO index between two consecutive segments. Finally, we calculate the intensity of each PDO pattern, and the standard deviation of intensities across 41 patterns is used to measure the unforced internal variability of PDO. Also the unforced internal variabilities of PNA and PNAforced SSTa patterns are calculated with the same method.

The response of PDO to PNA in observation and model simulation
To measure the observed North Pacific SST variability generated by PNA, we first use the seasonal cycle of cross correlation between PDO (HadISST1 dataset) and PNA (20CRv2c dataset) monthly expansion coefficient (EC) time series (both are filtered with 3-month running mean, Fig. 1b). The observed correlations show that atmospheric variations generally lead SST variations (e.g., Davis 1976;Deser and Timlin 1997). The prevernal PDO has a high correlation with prior wintertime PNA variability. The maximum correlation value appears between January PNA and February SST, indicating that the PNA has a strong influence on the North Pacific SST variability in the later month during the dominating season. To confirm how dose prevernal PDO (JMF) with a 1-month lag change by PNA forcing in wintertime (DJF), we calculate the regression of the Pacific SSTa (30° S-70° N, 110° E-60° W) in JFM onto the observed PNA index in DJF (as shown in Fig. 2c). Unless otherwise noted, the season is fixed on DJF for atmospheric field and JFM for oceanic field through the article. The regression pattern is similar to the observed PDO pattern (Fig. 2a), except for the weaker amplitude in SST variation. Note that the tropical Pacific SSTa excites atmospheric Rossby wave trains during ENSO events which force extratropical large-scale atmospheric teleconnections (e.g., Trenberth et al. 1998;Liu and Alexander 2007). And the PDO variability is affected by ENSO as well. Figure 1a shows that the observed PDO and ENSO index (here, defined as the Niño 3.4 index) are simultaneously correlated throughout the year. The high correlations at the lags from − 6 to 6 suggest significant interactions between ENSO and PDO. These results indicate that the ENSO signal entangles in both PNA and PDO signals.
To illustrate the PNA-induced SSTa in the North Pacific, the PLS regression is employed to eliminate the ENSO effect (Fig. 2e). Note that the PNA can be considered as the atmospheric internal variability, especially after removing the tropical signals, which still induces PDO-like SSTa pattern. To further inspect the relationship between the PDO and the PNA-forced SSTa, we compare the time evolution of the intensities of the North Pacific center of PDO (defined as the absolute regression coefficients averaged over the box in Fig. 2a) and the PNAforced SSTa (the same as the former, but in Fig. 2e) as shown in Fig. 1c. The time evolutions of the intensities are constructed over 33 80-year sliding windows (i.e., 1901-1980, 1902-1981, …, 1933-2012). The windows should be selected to ensure that the length of segments are longer than the maximum time scale of PDO cycle and that there are enough samples for processing. The red (blue) line represents the time evolution of the North Pacific center intensity of the PDO (PNA-forced SSTa).
Here the time evolution of PNAb intensity (green line in Fig. 1c) is also constructed by using the same methods (see Sect. 2.2 for details). The results show a robust connection between PDO and PNA-forced SSTa (Corr. = 0.76, p = 0.01). Also PNAb and PDO reveal the significant correlations 0.69 with p-value 0.02. Percentages at the lower right corner of the graphs indicate the contribution of PNA-forced SSTa to PDO intensity, which is defined as: is the PNA-forced SSTa (PDO) intensity in the ith segments and r is the correlation coefficient between these two variables. For example, if P reaches 100%, which means the variable x is exactly the same as y (there is no energy loss from x to y) and x is the only factor to cause the change of variable y. It shows that the wintertime PNA affects the North Pacific SST within one month, and it impacts on the PDO variability via the PNA-forced SSTa.
Similar phenomenon also appears in the patterns (Fig. 2b, d, f) of multi-model ensemble mean (MEM, see Table 1 for the details of the selected models), which demonstrate the model capability to simulate the observed atmospheric and oceanic variability in the tropical and North Pacific. When compared with the observe patterns, the amplitudes in Alaska Current area are overestimated for all the MEM patterns of Fig. 2. Besides that, the SST variations in KOE are also over-amplification in multi-model simulation (Fig. 2b), which may affect the position of the North Pacific center of PDO. In this case, we avoid the KOE area and set the maximum value of SSTa as the North Pacific center of PDO in each selected CMIP5 model (as shown in Fig. 3).

Projected future changes of PDO and the contribution of PNA under global warming
Based on the statistics of PDO patterns in CMIP5 models, we find that the location of the North Pacific center is different from model to model.  . 3) show an enhancement of the North Pacific center of the PDO under the global warming consistently. Besides, the overestimated KOE variability becomes weaker than historical scenario except NorESM1-ME model. To measure the magnitude of the North Pacific SST response to the intensification of PNA under the warming scenario, the historical PLS regression pattern of each model along with the pattern difference between the RCP8.5 and Historical scenarios is drawn in Fig. 4. The multi-model collectively shows the PDO-like SSTa response under the PNA forcing. However, the inter-model diversity promotes us to design tailored box for each model to measure the intensities of PNAb and its forced SSTa centers. Most of the selected models demonstrate good performances and show enhanced PNA-forced SSTa centers.
We have confirmed in our previous study that the PNAb in these 12 models are enhanced in RCP8.5 scenario as well (also shown by the dashed lines of Fig. 4) Fig. 1c location of PNAb in each model is denoted in the green box, and the PNAb intensity is defined as the area-averaged regression coefficients in that box.
The PNA-forced SSTa of the models only contribute to the negative (positive) value of SSTa in the North Pacific (Alaska Current). We use the same location (magenta boxes)    Note that most of the intensity changes are more than twice the internal variability, but the only model (HadGEM2-ES) with a weakened PDO shows that the PNAb and the PNAforced SSTa intensity changes are close to their internal variability. In general, through the numerical statistics of the North Pacific SSTa variability, the PDO will intensify in a warming climate and the enhancement of PNA-forced SSTa is the main cause. Although the enhancement of the North Pacific center of PDO in future projection (RCP8.5 scenario) is robust in the multiple models, the difference between individual models is primarily in magnitude. A naturally raised issue is, to what extent the inter-model differences in the North Pacific center of PDO are related to the PNA variability and associated SST responses.
As shown in Fig. 6, we find a close inter-model relationship among the intensity changes of the PDO, the PNAb and the PNA-forced SSTa in multi-model projections. For these three variables, we examined the correlations with each other, which are 0.78, 0.75, and 0.8 (all significant at the 0.01 significance level). The MEM of the PDO intensity change shows an increase of 0.13 °C with the 95% confidence interval of 0.06 °C to 0.20 °C (Fig. 6a), and the signal-to-noise ratio (i.e., change value divided by the estimation of unforced internal variability) is 2.22. Such an intimate linkage also holds for the MEM result, which shows an increase of 0.10 °C with the 95% confidence interval of 0.05 °C to 0.15 °C (Fig. 6c), and the signal-to-noise ratio is 3.59. Meanwhile the enhanced PNA in warming projection also shows an increase of 11.74 gpm with the 95% confidence range of ± 2.71 gpm, and the signal-to-noise ratio is 2.8. Apparently, the changes among the future projections of the selected models will enhance consistently. High correlations with the models indicate little difference in the intermodel diversity. It gives us a robust conclusion on the cause of PDO enhancement under global warming. The deepened trough of the PNAb will force a PDO-like SSTa pattern, and thus enhance the PDO.

Possible mechanism of the changes in North Pacific SSTa
To further inspect the evolution of the PDO intensity enhancement and the contribution of the PNA forcing under global warming, we examine the relationship between the North Pacific center of PDO and the PNA-forced SSTa bars denote the 95% confidence interval based on a two-tailed Student's t test. Note that the calculation is based on the absolute intensity of each variable, so that the positive (negative) value indicates the intensification (weakening) with their projected intensity changes. The time evolution of the PDO intensity is constructed over 107 95-year sliding windows (i.e., 1900-1994, 1901-1995, …, 2006-2100) with the combination of historical  and RCP8.5 (2006RCP8.5 ( -2100 run. Here the PDO intensity in each sliding window is derived from the absolute area-averaged regressions coefficients of the specific 95-year JFM SSTa against the PDO index in the selected boxes (the same boxes as in Fig. 3). Likewise, the time evolution of PNA and the PNAforced SSTa intensity are calculated based on the absolute area-averaged regressions of DJF 500-hPa GPHa and JFM SSTa against the normalized PNA index in the respective green and magenta boxes (Fig. 4) of each model. Note that similar results are obtained when using larger boxes. Here the statistical significance of correlation coefficient between the time series of PNA and PDO based on 95-year sliding windows is tested with a two-tailed Student's t test using the effective sample size, which takes into account the serial autocorrelation (Bretherton et al. 1999). The effective sample size (N e ) is defined as: N e = N × (1 − r x × r y )∕(1 + r x × r y ) , where N is the length of time series and r x (r y ) is the autocorrelation coefficient of the time series for variable x (y) at lag one.
As shown in Fig. 7, the MEM result clearly illustrates that the long-term evolution of the North Pacific center of PDO is in line with that of the PNAb and the PNA-forced SSTa. All the intensities are increased gradually, particularly from 2030 and onward. It is obvious that there is a close association between the multidecadal evolution of the PDO intensity and PNA-forced SSTa, with high correlation coefficient of 0.96 at 0.01 significance levels. Also the PDO intensity has a close relationship with the PNAb intensity (correlation coefficient 0.93). Note that the time coverage of the RCP8.5 data of CMCC-CESM is up to the year 2095 while that of the other models are all up to the year 2100, therefore the time-length of MEM is adjusted to the shortest one. The definitions of the percentages at the lower right corner of the graphs are the same as in Fig. 1c (see Sect. 3.1 for details). Figure 7 shows that all models except HadGEM2-ES can capture the close relationship between the time evolutions of the North Pacific center intensity of PDO and the PNA along with its forced SSTa intensity. It is noted that PDO intensity is more correlated to the PNA-forced SSTa than that is to PNA. The percentage of contribution in each model ranges from 46.3 to 72.8%. Therefore, the PNA and its forced SSTa will indeed intensify under global warming, and the PNA contributes about 67% of the intensification of PDO intensity change.
Overall, the above results suggest that, on average, models with stronger intensification of PNA variability tend to have larger magnitude of PNA-forced SSTa response, which contributes more to the enhancement of PDO in a warmer climate. For periods with a stronger PNAb, enhanced wind speeds and reduced air temperature and humidity cool the underlying North Pacific Ocean via surface sensible and latent heat fluxes, while northward advection of warm moist air heats the ocean near North America. This flux-driven SSTa pattern forced by PNA in the North Pacific closely resembles the observed PDO. Anomalous Ekman transports tend to amplify the fluxdriven pattern and causing the larger amplitude of SSTa in the North Pacific (Miller et al. 1994;Alexander and Scott 2008). The MEM projected change of − PNAb i is 11.74 ± 2.71 m. c As in (a), but for the PNA forced SSTa (projected change is 0.10 ± 0.05 °C) versus PDO intensity

Summary and discussion
This study explores the relative roles of PNA in the longterm changes of PDO under global warming. The observational reanalysis and selected CMIP5 outputs have been analyzed. The former is for understanding the relationship between PNA and PDO variations. The latter simulates similar temporal and spatial features to that of the observation, and is employed to study the changes of PNA and its effects on PDO in a warming environment. The main findings are summarized as follows: (a) The wintertime PNA variability leads PDO variations by one month with a high correlation coefficient above 0.6. It suggests that the PNA during the dominating season (DJF) has a strong influence on the North Pacific SST variability in JFM for the thermal inertia of the ocean. The CMIP5 models simulate similar spatial features of PNA and PDO patterns to the observation. Note that the simulation of PNA shows a stronger North America ridge and a weaker North Pacific trough than observation, while the simulated PDO pattern shows higher amplitude in Alaska current area and KOE than that in the observation. (b) The 12 selected CMIP5 models which PNAb intensities are increasing in future projection perform well in simulating atmospheric and oceanic variability. By comparing RCP8.5 and historical scenario of each model, we find that 11 out of 12 models exhibit significant changes of PDO intensity and 10 out of 12 models show robust changes of PDO-like SSTa forced by PNA consistently in North Pacific under global warming. Even without the tropical ENSO signals, the PNA still has influence on the North Pacific SSTa, which can be taken as a kind of subtropical atmospheric internal variability. Specifically, the North Pacific trough deepens about 11.47 ± 2.71 m (the 95% confidence inter-val, the same with below) in the twenty-first century, which is 2.8 times the unforced internal variability and drives 0.1 ± 0.05 °C enhancement in the amplitude of PDO-like SSTa variation. Note that the PDO intensity increases 0.13 ± 0.07 °C in the future projection, which is 2.2 times internal variability (Figs. 5,6). It is indicated that the PNA-forced SSTa is one of the main factor causing the enhancement of PDO. In addition, the models show that the changes of PDO, PNAb and PNA-forced SSTa are closely related to each other. (c) The evolution of the PDO intensity is ever-increasing in the twenty-first century. In the North Pacific, the PNAb and the PNA-forced SSTa also show enhancement in the future projection. The intensification of PDO-like SSTa forced by PNA is closely related to the intensification of PDO from early 21-st century onward. The contribution of the PNA-forced SSTa to the PDO intensification is 67.3% in the future projection, which is greater than that in the observation (50.7%). This may indicate a closer bond between the PNA and the PDO under global warming.
Overall, the above results suggest that, on average, models with stronger intensification of PNA variability tend to have larger magnitude of its associated SSTa response, which contributions to the enhancement of PDO in a warmer climate. A possible mechanism of this phenomenon is that the enhanced PNA causes larger wind anomalies which leads to the evaporations and Ekman transports increase in the North Pacific center, and thus reduces air temperature and humidity cool the underlying Ocean via surface sensible and latent heat fluxes, which causes the larger amplitude of North Pacific SSTa. This may shed some light on the projection of PDO and the active role of the PNA within the climate change under the global warming.  Fig. 7 Time evolution of the absolute PDO intensity (red line corresponding to right vertical axis; °C), PNAb intensity (green line corresponding to left vertical axis; m) and PNA-forced SSTa intensity (blue line corresponding to right vertical axis; °C) in winter over 95-year sliding windows from 1900 to 2100 for the MEM and individual models. Correlation coefficient between red line and green (blue) line is shown in green (blue) color, along with the corresponding p value. The percentage in black color at the lower right corner shows the contribution of the PNA-forced SSTa to the change of PDO intensity (the same as in Fig. 1c). The last year of each sliding window is indicated on the horizontal axis ◂