The index “Niño 3.4” (5°N–5°S, 140°E–145°W) defined over this region is widely used to represent the ENSO index, which is available at http://www.esrl.noaa.gov/psd/gcos_wgsp/Timeseries/Nino34/. The “Pacific Decadal Oscillation” (PDO) is a long-lived El Niño-like pattern of Pacific climate variability (Mantua et al., 1997; Christensen, 2013). The PDO index is defined as the projections of monthly mean SST anomalies onto their first EOF vectors in the North Pacific (north of 20ºN) (Mantua et al., 1997). The EOF vectors are derived for the period from 1900 to 1993, and climatological value is defined as monthly mean for the same period. Globally averaged monthly mean SST anomaly is subtracted from monthly mean SST anomaly of the Pacific Ocean north of 20°N before calculation of the first EOF vector for the purpose to remove the effect of global warming. The PDO index is available at https://www.esrl.noaa.gov/psd/gcos_wgsp/Timeseries/PDO/. The Interdecadal Pacific Oscillation (IPO) and the PDO are closely related modes of decadal to interdecadal climate variability in the Pacific Ocean (Mantua et al., 1997; Power et al., 1999; Folland et al., 2002; Christensen, 2013; Henley et al., 2015). The PDO can be considered as the North Pacific node of the Pacific-wide IPO (Power et al., 1999; Folland et al., 2002). The IPO is associated with a distinct “tripole” pattern of SST anomalies, and the IPO Tripole Index (TPI) is calculated as follows based on monthly global SST data (Henley et al., 2015): Firstly, subtract the monthly climatological value from SST in each grid cell to eliminate the seasonal cycle and calculate monthly mean SST anomalies (SSTA_i) in the three TPI regions using a chosen base period (1971–2000 used here). TPI regions are defined as follows: region 1 (25°N–45°N, 140°E–145°W); region 2 (10°S–10°N, 170°E–90°W); region 3 (50°S–15°S, 150°E–160°W). Secondly, the TPI is calculated as follows: T P I = S S T A 2 − S S T A 1 + S S T A 3 2

The TPI is also available at https://www.esrl.noaa.gov/psd/data/timeseries/IPOTPI/.

SST data are taken from the US National Oceanic and Atmospheric Administration Extended Reconstruction of SST (NOAA-ERSST; Smith & Reynolds, 2003), which is a global gridded dataset reconstructed based on historical observations of SST data using statistical methods. This dataset has spatial resolution of 2°×2° and temporal resolution of monthly, and the scope of this paper is global (0°–360°, 90°N-90°S).

GPI has important indicative significance for the generation of TCs and has been widely used in the study of TCs. Its definition is as follows: GPI = | 10 5 η 850 | 3 2 × ( 1 + 0.1 × V s h e a r ) − 2 × ( H 600 50 ) 3 × ( P I 70 ) 3 Where η 850 represents the vorticity of 850 h Pa, V s h e a r represents the vertical wind shear between 850 h Pa and 200 h Pa, H 600 represents the relative humidity of 600 h Pa, and PI is the potential intensity. PI describes the thermodynamically-based maximum TC intensity that the environment will support and thus is closely associated with TC genesis and intensification. In conjunction with the work of Emanuel (1999), it is defined as follows: P I 2 = C k C D × T s T 0 × ( C A P E ∗ − C A P E b ) where C k and C D represent exchange coefficient for enthalpy and drag coefficient respectively, T s and T 0 represent sea surface temperature and average outflow temperature respectively. CAPE is the convective available potential energy, and C A P E * represents the CAPE value when an air mass within the radius of the maximum wind speed reaches saturation for the first time under the current SST and barometric conditions. C A P E b represents the value of CAPE of the air parcel in the surrounding boundary layer atmosphere reduced along the isotherm to the radius of the maximum wind speed. In our study, large-scale environmental flow is mainly represented by 500 hPa geopotential height and 850–300 hPa vertical-mean wind.

Note that the GPI, PI and large-scale environmental flow are all calculated from the National Centers for Environmental Prediction–National Center for Atmospheric Research Global Reanalysis (NCEP–NCAR) (Kalnay et al., 1996) data. This reanalysis dataset is using a state-of-the-art system to perform data assimilation with spatial resolution of 2.5×2.5 and temporal resolution of monthly. In this paper, these data are selected from 1980 to 2015, and the spatial range is the Western North Pacific region (100°E-180°, 0°-60°N).

Best-track data are obtained from the International Best Track Archive for Climate Stewardship (IBTrACS) v04r00 (Knapp et al., 2010; Knapp et al., 2018), which is available at http://www.ncdc.noaa.gov/ibtracs/. There are four best-track data provided by four different agencies in the WNP: the Joint Typhoon Warming Center (JTWC), Japan Meteorological Agency (JMA), China Meteorological Administration (CMA), and Hong Kong Observatory (HKO). Note that only TCs with LMI ≥35 kt and ≤1,000 hPa are used in our analysis.

To further examine the impact of anthropogenic forcing on meridional migration of WNP TCs in future, we use an ensemble of TCs generated by 14 CMIP5 models, which includes the twenty-first-century RCP4.5 and RCP8.5 projection simulations (Camargo, 2013). RCP4.5 (RCP8.5) is a scenario that stabilizes radiative forcing at 4.5 (8.5) W m⁻² at the end of the century. CMIP5 models used in this paper include CanESM2, CCSM4, CSIRO Mk3.6.0, Fgoals-G2, GFDL CM3, GFDL-ESM2M, HADGEM2-ES, and INM-CM4.0. Ipsl-cm5a-lr, MIROC-ESM, MIROC5, MPI-ESM-LR, MRI-CGCM3, and NORESM1-M. RCP4.5 is a stabilization scenario and assumes that climate policies are invoked to achieve the goal of limiting emissions, while RCP8.5 corresponds to the pathway with the highest greenhouse gas emissions and is based on a presumption of no to little change from the present trajectories of greenhouse gas concentrations. RCP4.5 and RCP8.5 scenarios have better model coverage compared to the other scenarios, which is available at https://esgf-node.llnl.gov/search/cmip5/.

The significance of residual analysis is to check whether the trend of TC poleward movement is still significant after removing the trend of natural variability from the trend of TC poleward migration. ENSO (PDO, IPO) variability is excluded from the annual-mean TC latitude at the time of LMI ( φ LMI ) time series by regressing each individual φ LMI series from the various best-track data onto the index of the Niño-3.4 (PDO, TPI) and analyzing the residuals named ε _Niño-3.4 (ε _PDO, ε _TPI). ENSO, PDO, and IPO variabilities are removed from the φ LMI series by regressing the φ LMI series onto the indices of Niño-3.4, PDO and TPI, and analyzing the residuals (ε _ALL3).

The measurement of φ LMI is much less uncertain, since the absolute LMI can be easily determined based on when the TC has reached its maximum intensity during its lifetime (Kossin et al., 2014). Because of the relatively short period of observation after 1980, the estimated trend of φ LMI is highly sensitive to the sample size. For example, looking at the TCs in the IBTrACS best-track data during 1980–2013, the annual-mean φ LMI (TC intensity is measured by maximum wind speed) in the western North Pacific shows a poleward migration trend (0.48 ± 0.46 latitude decade⁻¹), which is basically consistent with results of the previous studies (Kossin et al., 2014; Kossin et al., 2016). However, if the sample length is extended to 2020, the poleward migration trend decreases significantly (0.29 ± 0.38 latitude decade⁻¹) due to the low φ LMI in 2014 and 2015 (Figure 1A; Table 1). Results based on other TC best-track data (i.e., JTWC, JMA, and CMA) also indicate that the poleward migration trend is no longer significant at the 95% confidence level when observations in 2014 and 2015 are considered (Table 1; Figures 1B–F) except HKO. More importantly, further analysis suggests that, rather than anthropogenic forcing, it is the natural variability that plays a critical role in modulating the migration of WNP TC track during the past 41 years.

The variability of the WNP TC activity has received a large amount of attention during the past 2 decades, and many studies have linked TC track variability to ENSO (Camargo et al., 2007a; Liu & Chan, 2008) on interannual timescale. In recent years, the PDO has been found to have important influences on large-scale flow patterns and thus the TC track variability over the WNP (Liu & Chan, 2008; Mei et al., 2015). Here, averages of the ENSO (Niño-3.4), PDO, IPO (IPO tripole index, TPI) (Henley et al., 2015) indices are taken over July-November (JASON), which is considered as the peak TC season (about 80% of annual WNP TCs occur in these months) (Kossin et al., 2016). As expected, there is an inverse relationship between φ LMI and the two indices (i.e., Niño-3.4 and PDO) with the correlation coefficients of −0.62 and −0.59 respectively, which indicate ENSO and PDO both have strong impacts on TC migration (Figure 2A). Moreover, compared with the ENSO and PDO, the IPO plays an even more important role in determining the WNP TC migration with a correlation coefficient of −0.70 between φ LMI and TPI. Almost all the peaks of IPO negative (positive) values correspond to high (low) values of φ LMI and φ LMI increased dramatically during the transition of the IPO from positive to negative phase in the late 1990s (Meehl et al., 2013; Henley et al., 2015) (Figure 2A). Note that the low φ LMI in 2014 and 2015 correspond to positive values of the TPI, which makes the poleward migration trend of φ LMI insignificant when the study period is changed from 1980–2013 to 1980–2015. On the other hand, as the low frequency component of the IPO signal can be obtained by applying a 13-years Chebyshev low-pass filter to TPI (Henley et al., 2015), we calculated the filtered φ LMI using the same filter, and found that the correlation coefficient between the filtered TPI and the filtered φ LMI is −0.52, which is significant at the 99% confidence level (Figure 2B). This result strongly suggests that there is a significant inverse correlation between IPO and the WNP TC migration on interdecadal timescale. In addition, the map of correlation coefficient between φ LMI and SST presents a pattern highly similar to the map of correlation coefficient between TPI and SST (Figures 2C,D). The latter displays a typical IPO pattern (see Figure 1 in the Henley et al., 2015) and thus provides robust evidence supporting the strong relationship between the WNP TC migration and IPO.

To further compare the contributions of various natural variabilities to the WNP TC migration, we calculate the residuals of annual-mean φ LMI regressed onto the JASON-mean Niño-3.4, PDO, and TPI indices (i.e., ε _Niño-3.4, ε _PDO, and ε _TPI), respectively. The residuals of the trivariate regression of φ LMI onto the three indices (i.e., ε _ALL3) are also calculated. For instance, ENSO is removed from the φ L M I time series by regressing the φ LMI series onto the Niño-3.4 index averaged over the peak TC season (July-November), and the residuals are analyzed. When the PDO and IPO variabilities are removed respectively, the residual TC migration rates (i.e., ε _PDO, and ε _TPI) are found to decrease substantially for various best-track data (Table 1; Figure 3). This suggests that the PDO and IPO considerably contribute to the TC migration rate. Moreover, as expected, the regression coefficients (i.e., R ²) are statistically significant, and the JASON-mean ENSO, PDO, IPO, and all these three indices together respectively explain 34.5%, 24.0%, 44.7%, and 44.7% of the φ LMI variance in the IBTrACS best-track data (Table 1). These results further suggest that natural variability, especially IPO (Wang et al., 2013; Henley et al., 2015), whose contribution to the change of φ LMI is similar to the total contribution of the three indices, is playing a critical role in the observed TC migration. Furthermore, when the IPO variability (i.e., time series of the TPI) is removed, the residual TC migration rate decreases significantly from 0.29 ± 0.38 latitude decade⁻¹ (i.e., φ LMI trend which is significant at 76% confidence level) to 0.03 ± 0.29 latitude decade⁻¹ (see φ LMI and ε _TPI in Table 1). This indicates that the poleward migration trend of the WNP TC track does not occur outside the IPO variability based on the observational records (1980–2020), and thus provides evidence against the proposed conclusions of previous studies on the contribution of anthropogenic factors to the TC poleward migration ((Kossin et al., 2014; Kossin et al., 2016). In this study, rather than studying the φ LMI trend, we put emphasis on the close relationship between TC migration and IPO, especially the modulation of IPO on the WNP TC migration.

Despite the strong agreement between the two indices of PDO and IPO (r = 0.75), the IPO index (i.e., TPI) is based on the difference between SST averaged over the central equatorial Pacific and that over the northwestern and southwestern Pacific, and thus is substantially different from the PDO index which is mainly determined by SST in the North Pacific (north of 20ºN) (Mantua et al., 1997). Thus, compared with the PDO index, the IPO index is more related to meridional and zonal gradients of SST over the WNP. The interdecadal changes of the SST gradients and related spatial SST distribution pattern make large contributions to the WNP TC migration, which will be discussed in detail in the next section.

As mentioned above, our observational results imply that the anthropogenic forcing does not contribute significantly to the observed WNP TC migration as expected. This also does not necessarily mean that the anthropogenic factors (e.g., global warming induced by greenhouse gases) will not significantly contribute to the poleward shift of the WNP TC in the future, since it is difficult to distinguish the contributions of anthropogenic forcing and natural variability based on the short time-period records (Allen et al., 2014; Lucas et al., 2014; Boo et al., 2015).

Thanks to the relatively long period (2006–2100) of RCP4.5 and RCP8.5 scenarios which is much longer than the period during which the known modes of natural variability are derived, it is relatively robust and reliable to draw the conclusions on the influence of anthropogenic forcing by using the results of these 14 CMIP5 models under the RCP4.5 and RCP8.5 scenarios. In the RCP4.5 projection simulations, the projected rate of migration of φ LMI (0.11°±0.12° decade⁻¹) is insignificant just under the 95% confidence level (Figure 4). The difference in φ LMI between RCP8.5 and RCP4.5 scenarios is hardly distinguishable, although the φ LMI in RCP8.5 scenario (0.12°±0.11° decade⁻¹) is just at the 95% confidence level (Figure 4). Thereby, the projected rates of WNP TC migration are closely related to the carbon emission intensity, as the migration rates, in terms of φ LMI trends, increase with the carbon emission intensity (Figure 4). In RCP4.5 scenario, the global warming remains unable to produce significant trend, indicating that the observed insignificant trend is confirmed in the model.

The φ LMI of a TC is determined by TC genesis positon, TC intensity evolution and TC track pattern. The first determines the starting position of a TC during its lifetime, the second evolves before LMI is reached and controls when the LMI would occur during the TC lifetime, and the third determines the path of the TC and thus controls where the LMI would possibly occur during the TC lifetime. In this study, these three factors can be analyzed from the aspects of GPI, PI, and large-scale environmental flow respectively (Sun et al., 2020). Firstly, the GPI is widely used to analyze spatial distribution of TC genesis potential and the influence of environmental factors on tropical cyclogenesis (Camargo et al., 2007b). Secondly, PI is a known major factor modulating the TC intensity evolution (Emanuel, 1999). Thirdly, WNP TCs are usually steered by large-scale environmental flow, especially the flow in the main TC activity region (Chan & Gray, 1982; Sun et al., 2015), which modulates the geographical properties of TC track over the WNP. Therefore, the change of φ LMI can always be attributed to changes in the GPI, PI and large-scale environmental flow.

Consistent with the previous studies (Burgman et al., 2008; Chen et al., 2008; Lee & McPhaden, 2008; Meehl & Arblaster, 2012; Henley et al., 2015), the IPO transferred from positive to negative phase in the late 1990s after the 1997–98 El Niño event. Thus, the period of 1980–1998 (1999–2015) can be considered as the positive (negative) IPO phase on interdecadal timescale (Figure 2B). Moreover, the time-averaged φ LMI in the negative IPO phase (21.2°N) is notably higher than that in the positive IPO phase (20.4°N), which is significant above the 95% confidence level by the Two-sample t-test. The positive and negative phases of IPO are characterized by different regional patterns (Power et al., 1999; Folland et al., 2002; Henley et al., 2015), which are expected to play an important role in determining the GPI and PI spatial distributions and the large-scale environmental flow in the main TC activity region, and eventually affect φ LMI of TCs over the WNP. However, based on what we understand, little is known about the impact of the regional pattern of IPO on the TC track migration. To further investigate the linkage between IPO and the WNP TC migration, we compare SST, GPI, PI, large-scale circulation, and TC exposure during the negative and positive IPO phases (Figure 5 and Figure 6). Note that, TC exposure is defined as the number of days per year per 2×2 latitude–longitude grid box when a TC center is located in that grid box.

During the transitional period from the positive IPO phase to negative IPO phase over the WNP, there are two main regions (Region 1 25°N-45°N, 130°E-180°E; Region 2 0°N-20°N, 135°E-155°E) where positive SST anomalies are significant at/above the 99% confidence level (Figure 5A). At first, the higher SST in Region 1 and Region 2 contribute to the significant increase of GPI in the three ocean areas (i.e., an ocean area to the south of Japan, an ocean area near the equator, and an ocean area near the central North Pacific). Compared with the first area that is to the north of the main WNP TC genesis area (see Figures 6C,D), the latter two areas are relatively farther from the main genesis region of WNP TC and thus contribute relatively smaller to the migration of TC genesis position (Figure 5B). Under the joint effects of GPI changes in the three areas, the time-averaged TC genesis position shifts poleward by 0.47 latitude (Figure 5E). Moreover, the higher SST in Region 1 leads to a relatively higher PI that favors further development of TCs at the high latitudes over the region (Figure 5C and Figures 6E,F), and eventually leads to the poleward shift of φ LMI . More importantly, the higher SST in Region 2 induces positive geopotential height anomalies in this region as the ocean warming is favorable for the development of the deep warm high (i.e., WPSH) in the non-equatorial region (>5°N), which favors westward expansion of the WPSH over Region 2 (Figure 5D and Figures 6G,H). The anticyclonic anomaly in the western boundary of the westward-expanded WPSH contributes to increases in the northward steering flow over the area of 10°N-22°N, 120°E-135°E, which tends to promote northward turning of westward moving TCs over the area (Figure 5D). In addition to the poleward shift of φ LMI , the intensified northward steering flow also contributes to the increase in TC exposure over the regions of Taiwan, Ryukyu Islands, the southern coast of Japan, and the southeastern coast of China, but decreases TC exposure in the regions of Hainan, Vietnam, and Philippines (Figure 5E). Note that the response of equatorial circulations to higher SST is substantially different to that of mid- and low-latitude circulations. Over the equatorial area in Region 2, the prevailing easterly winds basically agree with the theory of Gill model (Gill, 2010), which indicates that an equatorial heating (e.g., SST warming) produces easterly flow over the eastern region of the heating (Figure 5D). Overall, due to the response of GPI, PI and large-scale flow to IPO, the WNP TCs tend to originate and intensify at high latitudes and turn poleward over certain regions, which contributes to the poleward shift of the location of LMI during the IPO transition period from positive to negative phase in the late 1990s.