Emergent constraint on the projected central equatorial Pacific warming and northwestern Pacific monsoon trough change

The northwestern Pacific monsoon trough (NWPMT) deeply impacts socio-economic development and human security over East Asia by supplying moisture to the summer monsoon rainfall and modulating tropical cyclone activities. However, considerable inter-model spreads in the coupled model inter-comparison project phase 6 models make the future projection of the NWPMT less reliable. Here, we find that the inter-model spread of the NWPMT change is significantly correlated with the central equatorial Pacific sea surface temperature change, and mainly determined by the equatorial thermocline sharpness in the historical simulations. According to the emergent constraint method, the central equatorial Pacific SST would warm up about 6% slower than the multi-model mean with 56% uncertainty reduced. Correspondingly, the NWPMT would slacken westward with 36% uncertainty reduced. Results here emphasize the importance of examining and reducing systematic model biases in simulating thermocline sharpness that have been overlooked in past literatures, before achieving more reliable future projections.


Introduction
The northwestern Pacific monsoon trough (NWPMT) is associated with convergence of the southwesterly monsoon winds and the prevailing northeasterly trade winds over the northwestern tropical Pacific [1][2][3][4] in boreal summer.The induced weak vertical wind shear and cyclonic vorticity provide a favorable environment for the genesis and development of tropical cyclones (TCs) in boreal summer [5,6].Therefore, anomalous shifts in the NWPMT location have great socio-economic impacts on East and Southeast Asia where billions of people live in [7][8][9][10].
Although state-of-the-art climate models cannot simulate key features of TCs to investigate the future change of TCs activities directly [11][12][13], it has been suggested that projected extension (withdrawal) of the NWPMT is a good indicator of more (less) strong TCs over the northwestern Pacific in boreal summer under global warming [10].However, the projected zonal shifts of the NWPMT under global warming remain inconclusive due to the considerable intermodel spreads in CMIP6.Actually, the majority of the uncertainty in the projected shifts of the NWPMT has been attributed to the inter-model spread in the projected changes of the Pacific sea surface temperature (SST) in coupled model inter-comparison project phase 5 [10].
Here, we suggest for the first time that the intermodel spreads in the projected central equatorial Pacific SST and NWPMT changes are significantly correlated with the simulated equatorial Pacific thermocline sharpness in the CMIP6 models.Based on the emergent constraint method [30][31][32], we show that the SST change averaged over the central equatorial Pacific is projected to warm up about 6% slower than the raw projection with about 56% reduction in uncertainty.Correspondingly, the NWPMT is projected to slacken westward under global warming with about 36% reduction in uncertainty.

Definitions of changes of SST, thermocline depth, precipitation, relative vorticity, horizonal wind speed at 850hPa and the NWPMT in the future projections
The change of each variable mentioned above is computed by taking the difference between the future projections (2059-2098) and historical simulations .Additionally, the change of each variable is normalized by the globally averaged SST change in each model to account for the different sensitivity of each model to global warming.

Inter-model EOF analysis
The leading modes of the inter-model SST change over the equatorial Pacific (10 • S-10 • N, 120 • E-60 • W) are obtained by the conventional EOF method, applied to model-spatial dimensions: Here, ∆SST denotes SST change, m denotes the model number, s represents the spatial grid, and n denotes the mode number.Prime means the deviation from the multi-model mean.

Emergent constraint method
The emergent constraint method has been widely used to reduce uncertainty in future projections [30][31][32].By establishing a physical linkage between the target future projection Y and the historical simulation X, the observational constraint method can be described as follows: Here, b is the intercept, a = σY σX p is the linear regression coefficient, p is the correlation coefficient, σ Y , σ X is the standard deviation of Y, X across CMIP6 models, respectively.In this study, Y is the PC1 of the equatorial Pacific SST change among CMIP6 models, or the NWPMT index change, and X is the thermocline sharpness in the historical simulations.
Since the thermocline sharpness in the observed or assimilated datasets (i.e.X obs ) is used to constrain the Y, the uncertainty in the observations should be considered.With an additive-noise model under Gaussian assumptions that relates the observations to current climate [32], the signal-noise ratio (SNR) in the observed thermocline sharpness is derived to correct the linear regression coefficient (i.e.a) by multiplying x represents the inter-model spread of the simulated thermocline sharpness among CMIP6 models, and σ 2 obs denotes the spread of the observed or assimilated datasets.
After taking the SNR into account, the observational constraint can be expressed as: Correspondingly, the residual inter-model spread of Y (i.e.Y c ) is: Therefore, the reduction of total inter-model spread (1 − Once the SNR is large enough (i.e.SNR ≫ 1), the reduction of total inter-model spread would depend only on the linear correlation coefficient (i.e.p), or the impacts of SNR cannot be neglected.

Correction on the projected central equatorial Pacific SST change
After obtaining the constrained PC1, the central equatorial Pacific SST change in future projection can be corrected based on an EOF reconstruction as follows: Here, the subscript 'c' denotes the results constrained by the optimal PC1, mode number i is set to 1, and an overbar denotes the multi-model mean of the CMIP6 models.
Next, the difference of equatorial Pacific SST change between original multi-model mean and the projection after correction can be achieved by subtracting one from the other: The percent of difference in the projected SST change is estimated as follows:

Data
The monthly SST (tos), precipitation (pr), multiplelevel zonal (ua) and meridional wind (va), and ocean water potential temperature (thetao) in the 22 CMIP6 models (see figure 1(d)) are used here [33].In this research, seasonal-mean during boreal summer (i.e.July-October) are analyzed.A 40 year period from 1981 to 2020 is used for the historical simulation.Since the historical runs mostly end in 2014, outputs from 2015 to 2020 are from the shared socioeconomic pathway (SSP) 245 scenario [34].On the other hand, outputs from the last 40 year (2059-2098) of the SSP 245 scenario are used for the future projections.The MME mean is defined as the equal weight average of the 22 models.Only one ensemble (i.e.r1i1p1f1) of each CMIP6 model is used in this research.
The observational SST and precipitation in the ERA5 datasets [35] are used.The subsurface temperature for ocean includes the EN4.2.2g10 [36], EN4.2.2 c14 [37] and EN4.2.2l09 [38] from the Hadley Centre, the ocean data assimilation product (i.e.GODAS) [39] from the National Centers for Environmental Prediction (NCEP) and the Simple Ocean Data Assimilation Version 3 (SODA3) [40] from National Oceanic and Atmospheric Administration (NOAA).Three datasets from the Hadley Centre are averaged and labeled as the 'EN-group' to ensure equal weighting with GODAS and SODA3.The mean of the EN-group, GODAS and SODA3 is used in the emergent constraint method to yield the optimal constraint.The same 40 year period from 1981 to 2020 are used as the model baselines.All model outputs and reanalysis datasets are interpolated to the common 1 • × 1 • grid.

Results
The cross-equatorial flows turn into the southwesterly wind over the South China Sea (SCS), and converge with the prevailing easterly winds over the northwestern tropical Pacific in boreal summer (July-October).The induced positive relative vorticity extends southeastward from the SCS to the equatorial Pacific with negative relative vorticity lying on both sides (figure 1(a)).Following previous works [7,8], the NWPMT index is defined as the easternmost longitude of zero-zonal wind contour within positive vorticity region at 850 hPa (figure 1(a)).The spatial pattern of the NWPMT in the CMIP6 multimodel mean (MME) in both the historical simulation (figure 1(b)) and the SSP 245 scenario (figure 1(c)) are similar to that in the observation.
The change of NWPMT index between SSP245 (2059-2098) and historical simulation (1981-2020) is used to represent the zonal shift of NWPMT.Although the CMIP6 MME suggests that the NWPMT will extend slightly eastward in future, the inter-model spread is quite large.More specifically, half of 22 CMIP6 models project that the NWPMT would retreat westward, while other models project that it would extend further eastward (figure 1(d)).
Based on the high-quality observational datasets over the past 40 years, the zonal shift of the NWPMT seems to be not determined by the local SST anomaly, because the negative correlations between precipitation and SST anomaly signify that the SST anomaly is the response, rather than the forcing of the convective anomaly over the northwestern Pacific in boreal summer (figure S1(a)) [41].In fact, the interannual and interdecadal shifts of Moreover, the projected shift of the NWPMT is significantly correlated with the principal component of the leading mode (PC1, figure 2(c)).According to the Matsuno-Gill response [42], the CMIP6 models that project a faster SST warming over the central equatorial Pacific, by strengthening the westerly winds, positive relative vorticity, and rainfall over the northwestern tropical Pacific (figure S2), would thus project a larger eastward shift of the NWPMT (figure 2(c)).
Accordingly, a more reliable projection of the NWPMT can be achieved by making clear the SST change over the central equatorial Pacific under global warming.To unravel the root cause of the inter-model uncertainty in central equatorial Pacific SST change, correlation coefficients between the PC1 and the vertical potential ocean temperature gradient (i.e.dT dz ) in the historical simulations are examined.It is found that both meridionally and zonally averaged correlation coefficients are significantly negative along the thermocline (figures 3(a) and (b)).More specifically, the CMIP6 models with larger vertical ocean temperature gradient across the thermocline in the historical simulation correspond to a smaller PC1, and would project a slower central equatorial Pacific SST warming in the future.Hereafter, the vertical temperature gradient across the thermocline in each model is measured by the depth between 20 • C and 16 • C isotherms, and called the thermocline sharpness.
According to the emergent constraint method, the thermocline sharpness in the observed datasets can be used to derive a more robust and reliable central equatorial Pacific SST change than that in the original projection in the CMIP6 MME (figure 4(a)).It is found that only 6 out of 22 CMIP6 models locate within the range of the observed or assimilated datasets and could thus realistically simulate the thermocline sharpness, while four models simulate excessively sharp thermocline and 12 models simulate too diffuse thermocline.The PC1 constrained by the observed thermocline sharpness equals about −1.40, while the raw one projected by the CMIP6 MME equals to 0. Accordingly, the projected central equatorial Pacific SST change in boreal summer after emergent constraint would warm slower than that projected by the MME (figure S3(a)).
The physical mechanism underpinning the emergent constraint may be explained as follows: The CMIP6 models with sharper thermocline have larger vertical temperature gradient across the thermocline (figure 3), correspond to more pronounced shoaling of thermocline over the central equatorial Pacific under global warming (figure S4(a)), and thus could cool the SST more efficiently by the upwelling damping effect (figure S4(b)) [17,43,44].Compared with the CMIP6 models with more diffuse thermocline, the SST warming due to the anthropogenic forcing over the central-eastern equatorial Pacific would be greatly damped in the CMIP6 models with realistic thermocline sharpness.The induced faster SST warming over the western equatorial Pacific than that over the central-eastern equatorial Pacific would further benefit the slower SST warming over the centraleastern equatorial Pacific according to the Bjerknes feedback (figure 4(a)) [45].
In addition, the uncertainty in PC1 of equatorial Pacific SST change can also be greatly reduced after emergent constraint [30,32].The probability density function of standardized PC1 can be assumed to have a Gaussian distribution with a mean of zero and a standard deviation of 1 (figure 5(a)).Since the intermodel correlation coefficient between the projected PC1 and the simulated thermocline sharpness reaches 0.75 (figure 4(a)), about 56% of total inter-model spread of PC1 can be explained.The residual uncertainty in PC1 is 44% of the total inter-model spread, and the standardized PC1 after emergent constraint still has a Gaussian distribution but with a mean of −0.28 and a standard deviation of 0.66 (figure 5(a)).Thus, the SST change averaged over the central equatorial Pacific is projected to warm up about 6% slower than the MME (figure S3(a)), and its uncertainty is reduced by 56% after constraining the thermocline sharpness along the equatorial Pacific.
Due to the significant correlation between NWPMT change and the equatorial Pacific SST change (figure 2(a)), the more diffuse thermocline simulated in the CMIP6 MME would also exaggerate the projected expansion of the NWPMT under global warming (figure 4(b)).More specifically, the CMIP6  Note that the spread of thermocline sharpness among five observed or assimilated datasets is far less than that among 22 CMIP6 models (i.e.SNR ≫ 1), results above would keep almost the same after taking the differences of thermocline sharpness among five observed or assimilated datasets into account [30,32].Moreover, the thermocline sharpness index used in this study is free from the potential impacts of internal variabilities, such as inter-decadal Pacific oscillation (IPO) or Pacific decadal oscillation (PDO), because it keeps almost unchanged in each 40 year sliding period from 1951 to 2020 in five observed or assimilated datasets and 22 CMIP6 models (figure S5).

Conclusions
Changes in the zonal location of the NWPMT have great socio-economic impacts on East Asia.However, without a reliable future projection of the NWPMT change, no proper measures for disaster prevention and mitigation can be undertaken.Our results here suggest that the SST over the central equatorial Pacific would warm up slower than the CMIP6 MME projection with about 56% uncertainty reduced, and the NWPMT would slacken westward with about 36% uncertainty reduced by constraining thermocline sharpness along the equatorial Pacific.Considering the anti-phase shift of NWPMT and northwestern Pacific Subtropical high (NWPSH) [7], our results here are in agreement with the previous work that projected the future strengthening of the NWPSH based on the CMIP5 models [30].Moreover, a slackened NWPMT might imply a stronger summer monsoon rainfall over east Asian, but fewer typhoon landfalls and a reduction of strong TCs over the northwestern tropical Pacific under global warming [10,30], which aligns with statistical results derived from reanalysis datasets [46,47].
Although the mechanism underpinning the emergent constraint method is straightforward and easy to understand, the cause of the systematic bias in simulating the thermocline sharpness requires further investigations [48,49].

Figure 1 .
Figure 1.Geographical location of the northwestern Pacific monsoon trough (NWPMT) and its projected shift in the coupled model intercomparison project phase 6 (CMIP6) future projections.(a), (b), (c) Climatological locations of the NWPMT (green dot) in the observed datasets (1981-2020), and the multi-model mean of the historical simulations (1981-2020) and future projections (2059-2098) in boreal summer, respectively.The color shading denotes the relative vorticity (10 −5 s −1 ).The black and red solid curves denote the zero-zonal wind contour (m s −1 ) and the zero relative vorticity contour, respectively.(d) Rank of the NWPMT index change (longitude • C −1 ) in the CMIP6 future projections.The black bar represents the multi-model mean.Note that the zonal shift of NWPMT has been normalized by the corresponding globally averaged SST change in each model to account for different sensitivity of each model to global warming.

Figure 2 .
Figure 2. Correlation between the projected NWPMT and the central equatorial Pacific sea surface temperature (SST) change.(a) Inter-model correlation coefficients between the NWPMT change and the tropical Pacific SST change in the 22 CMIP6 models.Hatching indicates correlation coefficients significant at the 95% confidence level with a Student t-test.(b) The first empirical orthogonal function (EOF) mode of equatorial Pacific SST change among CMIP6 models.The spatial pattern correlation between panel (a) and panel (b) is indicated in the upper right of panel (a), and the variance contribution of the leading mode is indicated in the upper right of panel (b).(c) Scattor plot constructed by the principal component of the leading mode (PC1, X axis), and NWPMT change (Y axis) in each CMIP6 model.The red solid line denotes the regression line, while the gray dashed lines represent the 95% confidence range of the linear regression.Their correlation coefficient is indicated in the upper left.
NWPMT are mainly determined by variabilities of the central equatorial Pacific SST (figures S1(b) and (c)).Consistently, the projected shift of the NWPMT is also significantly correlated with the projected change in the central equatorial Pacific SST (figure 2(a)).Interestingly, the spatial pattern of the intermodel correlation coefficients between the projected shift of the NWPMT and the change in the equatorial Pacific SST can be well captured by the leading empirical orthogonal function (EOF) mode of the inter-model spread of the equatorial Pacific SST change; the pattern correlation coefficient reaches 0.92 (figures 2(a) and (b)).The leading EOF mode accounts for about 45% of inter-model variance among CMIP6 models, and shows the SST warming center over the central equatorial Pacific (figure 2(b)).

Figure 3 .
Figure 3. Correlation coefficients between PC1 and vertical ocean temperature gradient.(a) Zonal cross section of correlation coefficients averaged along the equatorial (5 • S-5 • N) Pacific.(b) As in (a), but for the meridional cross section averaged over 150 • E-150 • W. Stippling indicates correlation coefficients significant at the 95% confidence level with a Student t-test.The red and black solid curves denote the 20 • C and 16 • C isotherms in the mean of the observed datasets and the multi-model mean of historical simulations in the 22 CMIP6 models, respectively.

Figure 4 .
Figure 4. Inter-model relationship between simulated thermocline sharpness and future projections.(a) Inter-model correlation between thermocline sharpness (m) that is measured by the difference between the depth of 20 • C and 16 • C isotherms over the central equatorial Pacific (5 • S-5 • N, 150 • E-150 • W) and PC1 based on 22 CMIP6 models.The red solid line denotes the linear regression line, while the gray dashed lines represent the 95% confidence range of the linear regression.The blue, purple and cyan vertical dashed lines denote the thermocline sharpness in the EN-group, GODAS and SODA3, respectively.The red solid vertical line denotes the mean of the observed datasets, and the induced optimal projection according to emergent constraint is represented by the red horizontal solid line.(b) As in (a), but for thermocline sharpness and the projected NWPMT shift.The inter-model correlation coefficient is indicated in the upper left of each panel.

Figure 5 .
Figure 5. Probability density function of original and constrained future projections.(a) Probability density function (PDF) of the standardized principal components (PC1) of the leading mode of the equatorial Pacific SST change differences among 22 CMIP6 models are generated under Gaussian assumption.The blue (red) curve denotes the PDF distribution of the original (constrained) PC1.The values in parentheses at the top are mean and standard deviation of the Gaussian distribution.The red and blue dots at the bottom represent six CMIP6 models with realistic thermocline sharpness and other CMIP6 models, respectively.(b) As in (a), but for the PDF distributions of the standardized NWPMT index change in CMIP6 models.