The annual mean global mean surface temperature (GMST) broke the previous records 3 years in a row during 2014–2016 and by a large margin. The magnitude and duration of this record‐breaking event were highly unusual (Mann et al., 2017 ) and exceeded previous events associated with strong El Niño. As a consequence, global warming has passed the 1°C milestone relative to the preindustrial period (1851–1880), with some monthly GMST values in 2015 and 2016 even close to the 1.5°C threshold. This jump of record warm GMST also represents a sizable fraction of the 2°C stabilization goal set by the 2015 Paris Agreement (Paris Agreements, 2015 ) to avoid dangerous climate change. The record high monthly and annual mean GMSTs in 2015 and 2016 were coincident with outbreaks of extreme weather worldwide (Herring et al., 2016 , 2018 ), extensive melting of polar ice (Nicolas et al., 2017 ; Simpkins, 2017 ), and the third ever global coral bleaching event (NOAA, 2015 ). Given these significant impacts, it is extremely important and timely to understand the mechanism for this rapid warming and provide perspectives for future large record‐breaking events of GMST. Here we examine various observational data and the simulations and projections by 40 climate models to address these issues.

Climate sensitivity values in these models depend on representations of different feedback processes in the atmosphere. Different approaches have been used to estimate climate sensitivity in climate models (Yoshimori et al., 2016 ). CMIP5 used the approach of Gregory et al. ( 2004 ) to calculate effective climate sensitivity. The values range from 2.1° to 4.7°C across the CMIP5 models (Table S2 ), therefore covering a large uncertainty space. As the transient climate response values show a narrower spread across the CMIP5 models, we find that the use of effective climate sensitivity can enhance correlation analyses in this study.

In terms of modeling data, we use 40 Coupled Model Intercomparison Project Phase 5 (CMIP5) models and their piControl and 20th century simulations, as well as 21st century projections under four Representative Concentration Pathway (RCP) emission scenarios (Table S1 ) (Flato et al., 2013 ; Taylor et al., 2012 ). Based on the 200 year piControl simulations, the CMIP5 models show a range of internal variability of GMST (as the standard deviation of the annual time series) from 0.06° to 0.16°C, with an ensemble mean of 0.11°C (Table S2 ). For observational estimates, we apply the Ensemble Empirical Mode Decomposition method (Wu et al., 2011 ) to remove the nonlinear secular trend from the observed GMST time series (Figure S5 a). We also subtract the forced GMST response in the CMIP5 historical simulations from the GMST observations (Figure S5 b). All of these methods consistently indicate a 0.11°–0.12°C standard deviation of internal GMST variability. Some observational data such as HadCRUT4 are spatially incomplete and blend surface air temperature and SST. This may slightly influence the estimated internal variability of GMST and therefore the data‐model comparison (Cowtan et al., 2015 ). ENSO variability in the tropical Pacific contributes to GMST variability.

The Multivariate El Niño–Southern Oscillation (ENSO) Index (MEI) is a comprehensive index for quantifying ENSO strength (Table 1 ) (Wolter & Timlin, 2011 ). According to the MEI time series (Figure S3 ), the 2015/2016 El Niño was strong but not as strong as the 1982/1983 and 1997/1998 El Niño. Although there were signs of El Niño development in 2014 (Hu & Fedorov, 2016 ), 2015 was the real onset of the recent strong El Niño. As another widely used index, Niño3.4 confirms the strong El Niño in 1982/1983, 1997/1998, and 2015/2016 (Figure S3 ) (Rayner et al., 2003 ). To identify decadal and longer time scale variability in the Pacific, we apply a 9 year low‐pass filter to the annual mean index of the Pacific Decadal Oscillation (PDO) (Figure S4 ). We quantify the contribution of PDO to the ocean heat accumulation and sea level rise in the western tropical Pacific (WP) during 1993–2012. The use of the Interdecadal Pacific Oscillation (IPO) index gives similar results.

We choose two tide gauge data in the northwestern tropical Pacific (NWP) with relatively long and complete records (Table 1 ) (Holgate et al., 2013 ). The relative sea level data at Kwajalein (167.74°E, 8.73°N) and Guam (144.65°E, 13.44°N) span from 1950–2015 and 1948–2016, respectively. For satellite altimetry data, we use the Archiving, Validation, and Interpretation of Satellite Oceanographic data (AVISO) (Table 1 ) (Ducet et al., 2000 ). The data show dynamic sea level (i.e., departure of sea surface height from the geoid) at an eddy‐resolving scale and over a global domain. We find that the tide gauge and altimetry data are highly correlated at Kwajalein and Guam during 1993–2016 ( r = 0.98) and show similar interannual changes and decadal/interdecadal trends (Figure S2 ). The long‐term global mean sea level reconstruction is from Church and White ( 2011 ).

As closely related fields, we use three data sets for ocean heat content (OHC) and temperatures including Levitus (Levitus et al., 2012 ), Ishii (Ishii & Kimoto, 2009 ), and Argo (Roemmich & Gilson, 2009 ) (Table 1 ). The OHC anomalies are relative to the WOA09, 1981–2010 and 2004–2015 climatology in the three data sets, respectively. The values are integrated in the upper 700 m and across each grid box. It should be noted that the three data sets are not independent. The Argo coverage is usually not considered adequate until 2005.

Table 1 summarizes the observational data sets used in the present study. Six data sets are available for GMST including Hadley Center and Climate Research Unit analysis version 4 (HadCRUT4) (Morice et al., 2012 ), Goddard Institute for Space Studies (GISS) Surface Temperature Analysis (GISTEMP) (Hansen et al., 2010 ), NOAA Global Surface Temperature (NOAAGlobalTemp) (Vose et al., 2012 ), Japanese Meteorological Agency (JMA) (Ishihara, 2006 ), Berkeley (Rohde et al., 2013 ), and Cowtan&Way (Cowtan & Way, 2014 ). The detailed data descriptions and uncertainty quantification can be found in the references (Figure S1 in the supporting information ). All six GMST data sets were downloaded as of 7 March 2017.

3 Results

3.1 Observed Big Jump of Record Warm GMST Over 2014–2016 The mean of the six GMST data sets indicates that the past 3 years (2014–2016) taken together broke the previous record warm year set in 2010 by a large amount of 0.24°C (0.22°–0.27°C, mean and range of the six data sets) (Figure 1). Year 2014 alone was a record‐breaking year in five of the six data sets, with a difference from the previous record that is within the observational uncertainty range (Figure S1). However, the 3 year jump of 0.24°C represents a large change compared with GMST variability. Such events with similar or greater record‐breaking magnitudes over consecutive record‐breaking years (0.24°C+ events hereinafter) are highly unusual in the instrumental records (Figure 1) (Mann et al., 2017). Most of the record increase in GMST over the past 3 years occurred in 2015 and 2016 during a strong El Niño. Years 2015 and 2016 surpassed 2014 by 0.15°C (0.11°–0.19°C) and 0.21°C (0.18°–0.24°C), respectively. In 2015–2016, forecasts were published (Peyser et al., 2016; UK Met Office, 2015) predicting this jump in GMST using either operational forecast systems (Folland et al., 2013) or Pacific sea level anomalies along with a dynamical‐statistical method (Peyser et al., 2016). The subsequent observations confirm these GMST predictions. Figure 1 Open in figure viewer PowerPoint Observed GMST anomalies and record‐breaking events. (a) Time series of six GMST data sets and the mean of the six data sets. (b) Cumulative record‐breaking magnitude of GMST during consecutive record‐breaking years after 1980 (mean and range of the six data sets). For the 2014–2016 event, the cumulative record‐breaking magnitude of GMST is the difference between 2016 and 2010 (the previous record prior to the event). All anomaly values are relative to the 1981–2010 period. The mean in Figure 1 a is calculated for the common period of the six data sets (1891–2016). The red dots denote record‐breaking years of GMST after 1980. See Figure S1 for plots of each data set and associated uncertainty. In three longer data sets (HadCRUT4, Berkeley, and Cowtan&Way), the 1877–1878 GMST spike was the only other possible 0.24°C+ event in the instrumental period (Figure S1). This 19th century event was prior to large anthropogenic climate forcing and warming, but it is also associated with relatively large uncertainty. The large drop of GMST in the following year (1879) suggests that this GMST spike was likely caused by natural and transient processes, particularly a possible strong El Niño in the Pacific (Aceituno et al., 2009). By contrast, recent predictions of GMST based on operational systems suggest that 2017 will be among the warmest years and the GMST, although unlikely to exceed the 2016 record, will remain at a high level after the 2014–2016 jump (UK Met Office, 2016).

3.2 Unusually Large Oceanic Heat Releases From WP/NWP during 2014–2015 The 2014–2016 GMST jump is a combined effect of the long‐term increase in the atmospheric GHG concentration and the short‐term strong El Niño in the Pacific Ocean. The warming associated with the El Niño should not be considered as purely natural variability because, in addition to warming the SSTs along the eastern equatorial Pacific, the El Niño also released excess heat sequestered into the subsurface layer (50–300 m) of the WP, especially the NWP east of the Philippines (Figures 2, 3, and S6). This heat had built up over the previous two decades (Figures 2a and S7), leading to rapid sea level rise in the WP (up to four times faster than the global mean) since the 1990s due to the dominant thermosteric effect (Figures 3c and S2) (Griffies et al., 2014; Peyser et al., 2016). Figure 2 Open in figure viewer PowerPoint 18 J, Levitus data). (a) Total heat accumulation during 1993–2012 based on the linear trend. (b) Ocean heat release (negative values) during 2013–2015. Notice the different scale. (c) The 1993–2012 ocean heat accumulation due to the transition of PDO/IPO to their negative phase. (d) Difference between Figures Interdecadal accumulation and rapid release of ocean heat in the upper 700 m (10J, Levitus data). (a) Total heat accumulation during 1993–2012 based on the linear trend. (b) Ocean heat release (negative values) during 2013–2015. Notice the different scale. (c) The 1993–2012 ocean heat accumulation due to the transition of PDO/IPO to their negative phase. (d) Difference between Figures 2 a and 2 c indicating the ocean heat accumulation due to factors external to the Pacific. The green boxes in Figures 2 a, 2 c, and 2 d indicate WP (120°E–180°E, 20°S–20°N), SWP (145°E–195°E, 0°–15°S), and NWP (120°E–180°E, 3°N–18°N), respectively. The circles mark Kwajalein and Guam (Figure 3 c). See Figure S6 for similar analysis with the Ishii data. Figure 3 Open in figure viewer PowerPoint Consistent changes of observed ocean temperature, heat content, and sea level in the WP and NWP. (a) Time series of the upper 700 m OHC anomalies in the WP from three data sets. The dashed lines indicate 20 year linear trends prior to the 1981–1982, 1996–1997, and 2013–2015 ocean heat release events associated with strong El Niño. (b) Anomalies of the area mean ocean temperatures in the WP as a function of depth and time (Levitus data). The vertical dashed lines indicate the 1981–1982, 1996–1997, and 2013–2015 ocean heat releases from the subsurface layer (horizontal dashed lines). (c) Two long‐term high‐quality tide gauge data of relative sea level at Kwajalein and Guam in the NWP. The dashed lines show the linear trend during 1950–1992 and 1993–2012. The reconstructed global mean sea level rise is also plotted with uncertainty. (d) The 1993–2012 ocean heat accumulation in the NWP (based on the linear trend) and the contributions from PDO (Figure S8 ) and other factors external to the Pacific (total minus PDO), and the magnitudes of the large drops during 2013–2015 and 1996–1997. The tide gauge data at Kwajalein and Guam are analyzed similarly. The linear trends of the observed upper 700 m OHC indicate a total heat accumulation of 3.5–4.5 × 1022 J in the WP during 1993–2012 (Figure 3a). As shown in the low‐pass filtered time series, the gradual buildup of ocean heat in the WP/NWP started from around 1993–1994 (Figure S7). The heat accumulation in the WP is an important component of global ocean heat uptake during the same period (Gleckler et al., 2016). Both internal variability and external forcing contributed to this unprecedented heat accumulation in the WP. On decadal/interdecadal time scales, the PDO or IPO is the dominant variability mode in the Pacific region that can influence OHC and GMST (England et al., 2014; Mantua et al., 1997; Medhaug et al., 2017; Meehl et al., 2016; Nieves et al., 2015; Steinman et al., 2015). It should be noted that the mechanisms of the PDO/IPO are still under investigation about whether it is an intrinsic variability mode, a mixture of tropical and higher latitude phenomena (Newman et al., 2016), or involves a forced component (Smith et al., 2016). The regression of the upper 700 m OHC onto the observed PDO or IPO indices since the 1950s suggests that the transition of PDO/IPO to their negative phase during 1993–2012 had caused heat accumulation primarily centered in the southwestern tropical Pacific (Figures 2c, S4, and S6). By contrast, the buildup of excess ocean heat in the NWP east of the Philippines is mainly attributable to factors external to the Pacific, including GHG forcing and possible remote effects of other ocean basins (Figures 2d, 3d, and S6) (Luo et al., 2012; McGregor et al., 2014). As a good indicator of OHC, two long‐term high‐quality tide gauge data from the NWP confirm the unusual heat accumulation during 1993–2012 (Figure 3c). The relative sea level data at Kwajalein and Guam since 1950 are well correlated, suggesting common mechanisms. Sea level trends show a clear change during the 1990s at both locations. According to a clear linear trend, the total sea level rise over 1993–2012 was 194 and 186 mm at Kwajalein and Guam, respectively. The regression of the detrended tide gauge data onto the PDO index suggests that natural variability can explain 53 and 75 mm sea level rise at the two sites, respectively (Figures 3d and S8). Thus, other factors including GHG forcing and remote oceanic effects likely played a more important role in causing the rapid sea level rise in the NWP during 1993–2012 (Figures 3c and S2) (Hamlington et al., 2014; Han et al., 2014). The excess subsurface ocean heat in the WP/NWP rapidly resurfaced and was then released to the atmosphere during 2014–2015 (Figures 2, 3, and S9). Here 2013 is used as the reference level for the heat release during 2014 and 2015. Compared to 2013, the WP OHC shows a significant drop by about 4.0–5.2 × 1022 J in 2015 (Figure 3a), with the maximum heat release at the NWP (Figures 2b and S6). Our results indicate that in addition to the natural and normal amount of ocean heat release, the strong 2015/2016 El Niño also completely released the ocean heat accumulated in the NWP during 1993–2012 (Figure 3d). The intense heat release resulted in a drastic sea level fall in the NWP by up to 300 mm during 2013–2015 (Figure S2). As a consequence, the once fastest sea level rise east of the Philippines as evidenced in altimetry data tends to subside after the 2014–2015 event (Figure S2). EOF1 of the OHC in the tropical Pacific reveals an east‐west seesaw pattern in the tropical Pacific (Figure S10). This seesaw also reflects vertical redistribution of ocean heat in the tropical Pacific associated with ENSO and PDO/IPO as well as induced by external forcing (Figure S9). Following a gradual decline during 1993–2012, PC1 shows a large jump during 2013–2015 exceeding the previous events associated with strong El Niño (Figure S10). It indicates an eastward movement and rapid upward emergence of massive ocean heat previously sequestered in deeper oceans (Figure S9), thereby directly boosting the Pacific SST and GMST (Figure 1). During the large 1997–1998 El Niño, GMST also broke the previous record warm year set in 1995 by a large margin of 0.20°C (0.18°–0.22°C) (Figure 1). Conventional indices such as MEI and Niño3.4 suggest that the 1997–1998 El Niño was at least as strong as the 2015–2016 one (Figure S3). In terms of ocean heat release from the WP, and especially the NWP, however, the 2013–2015 event was clearly stronger and more pronounced than the previous ones during 1981–1982 and 1996–1997 (Figures 3 and S10). In fact, after two decades of unprecedented heat accumulation, the OHC in the WP dropped all the way back to the 1997 level in 2015. The ocean heat release started from 2014, coincident with the increase in GMST.