Borneo Vortices in a warmer climate

Borneo Vortices (BVs) are weather systems that are responsible for devastating hydro-climatic extremes and significant losses of life and property in Southeast Asia. The typical resolution of most current climate models is insufficient to resolve these high-impact, synoptic-scale weather systems. Here, an ensemble of high-resolution models projects that future BVs may become less frequent and more stationary, driven by the weakening of the Northeast monsoon flow and associated cold surges across North Borneo. However, substantial increases in both the intensity and the total amount of precipitation from BVs are projected. Such changes are driven by the more humid and convectively unstable lower troposphere. As a result, the contribution of BVs to the accumulation of both total precipitation and extreme precipitation is projected to increase considerably in the vicinity of the southern South China Sea, making individual BVs more threatening to the adjacent coastal regions.


INTRODUCTION
Borneo Vortices (BVs) are cyclonic weather systems that are often associated with extreme precipitation and flood events over Southeast Asia. They contribute significantly to extreme winter precipitation over northern Borneo (up to 60-70% of the total extreme precipitation above the 95th percentile daily precipitation rate) and the east coast of Peninsular Malaysia (15-25%) during the Northeast monsoon season 1 . Such weather systems feature a shallow cyclonic circulation in the lower troposphere 2 . Their motions are usually quasi stationary 3 while their equatorward propagation can lead to torrential rainfall and flooding over populated cities in Indonesia, which include Jakarta 4 and Padang 5 . Most BVs occur in boreal winter, with a peak in December-January, 1,6,7 as the formation of BVs is facilitated by the northeasterly cold surge within the Northeast monsoon flow due to its meso-scale convergence 8 and interactions with local topography 8,9 . Their maintenance is considered to be associated with the release of mixed barotropic and baroclinic instability 10 . In extreme cases, BVs can develop into tropical cyclones (e.g. Tropical Storm Greg in 1996 and Typhoon Vamei in 2001) and have further devastating effects on life and property 7,9,11 . This has been highlighted by the recent Tropical Depression 29 W which developed from a BV that formed to the northwest of Borneo on 13 th December 2021. In a week, this storm made landfall on the east coast of Peninsular Malaysia and triggered extreme flooding, leaving approximately 70,000 people displaced and severely disrupting infrastructure 12 .
Although the socio-economic impacts of BVs can be substantial, there remains relatively little research seeking to understand their response to global warming. There are two reasons for this knowledge gap. First, the relationship between BVs and precipitation extremes stems from the vigorous rainfall-favoring environments near the vortex center, including the intense lowtropospheric convergence of cyclonic winds and the development of cumulus convection at a fine spatial scale of 20-200 km 8 . Simulating these fine-scale features using climate models requires a sufficiently fine horizontal resolution of at least~65 km, as the frequency and wind structures of BVs are poorly resolved when coarser model resolutions are used 1 . Thus, current general circulation models (GCMs), such as many of the GCM members (at resolutions corresponding to~100-300 km) 13 of the Coupled Model Intercomparison Project Phase 6 (CMIP6) and used by the latest IPCC Assessment Reports, cannot adequately resolve these synoptic processes. Second, the dynamical characteristics of BVs are usually distinctly observed when an equatorward cold surge exists within the Northeast monsoon flow over the South China Sea (SCS) 14 . The development of such a winter circulation relates to winter high-pressure disturbances over the upstream midlatitudes via a teleconnection associated with meridional dispersion of the zonal Rossby wave group 15 . Precise detection and analysis of these multi-scale characteristics are essential for understanding the changes in BVs and their dynamical drivers, but are technically challenging.
Recently, state-of-the art climate model outputs from the High Resolution Model Intercomparison Project (HighResMIP 16 ), an experimental protocol aiming to deliver a multi-GCM ensemble at horizontal resolutions finer than the current CMIP6 experiments, have demonstrated an ability to simulate the historical climatology 17,18 and natural variability 19 of regional precipitation over Southeast Asia. HighResMIP also plays an important role in understanding near-future projected changes in regional precipitation and their uncertainties due to changes in climate model resolution 16 . This implies a necessity, but also an opportunity, to fill the knowledge gap on possible future changes in BVs using the state-of-the-art high-resolution GCMs.
Here, we undertake the first study to project the near-future response of BV characteristics to a high-level emission scenario of greenhouse gases as used for the current IPCC AR6 assessment report. This is achieved by using six high-resolution GCMs (referred to as HR-GCMs and summarized in Supplementary Discussion) of HighResMIP and an advanced objective identification algorithm 1,20,21 (see "Methods") to detect the full-lifetime trajectories of BVs in these models. To demonstrate the importance of the use of high-resolution climate models for BV projection, a validation of the models and their low-resolution counterparts, compared to one of the latest high-resolution climate reanalysis datasets (ECMWF Reanalysis v5, ERA5 22 ), is discussed in detail in Supplementary Discussion. This study investigates changes in BV characteristics (including frequency, intensity, precipitation and associated large-scale environmental controls) for the near future climate (2030-2050) under the Shared Socioeconomic Pathway 5-85 (SSP5-85) scenario 23 relative to the historical baseline period . This analysis provides insight into the hydrological changes associated with high-impact weather systems over the next few decades, which require relatively urgent strategies of mitigation and adaptation to deal with changes in these events. Figure 1 shows anomalies of BV frequency relative to the mean frequency for the historical baseline period . The historical time series for the ERA5 reanalysis dataset exhibits a multi-decadal variability of BV frequency, as the period 1986-1999 shows relatively high frequency and a lower frequency is shown for 2000-2014 for all-regions (Fig. 1a). The cause of such variability remains unclear, with a possible association with decadal climate variability in the tropical Pacific 24 Table 2), indicating a limited but noticeable model uncertainty. However, all the projected trends in BV frequency are within the magnitude of historical variations for both ERA5 and the GCM simulations, thus they should be interpreted with caution due to uncertainties from internal climate variability.

Changes in vortex features
The spatial distribution of BVs in the historical baseline period ( Fig. 1d-f) and projected future changes ( Fig. 1g-i) are analyzed using the multi-GCM ensemble mean. Figure 1g indicates a projected decrease (p value < 0.05 according to the Wilcoxon ranksum test) in the genesis density of BV (up to 1.8 vortices per 3.5°s pherical area per season) to the north of Borneo which is the primary BV formation region (Fig. 1d), a change consistent with the decreasing trend in BV frequency (Fig. 1c). For the mean track density of BVs, a significant reduction (p value < 0.05) by up to 2.5 vortices per 3.5°spherical area per season is projected over the south of Peninsular Malaysia and the east of Sumatra. Some decrease (p value < 0.05) is seen over the most active region near the western coast of North Borneo. For the lysis density, Fig. 1i shows a significant decrease by up to 2.1 vortices per 3.5°s pherical area per season near the Karimata Strait (105°E, 0°N). This reveals a decrease in the southwestward propagating BVs; hence, BVs are projected to be confined within the offshore region of North Borneo and become more stationary. Figure 2 shows the spatial distribution of total precipitation during ONDJFM and the contribution from BVs. There are significant increases (greater than 90 mm per season, p value < 0.05) in BVassociated precipitation (defined in Methods) over most of the southern part of the SCS (Fig. 2e). Although most of the statistically significant changes are over ocean, increases in BVassociated precipitation are found in most of Borneo, with the most significant increases to the north and east of the region. Also, the fractional contribution of BVs to the total precipitation during ONDJFM increases significantly over these regions (Fig. 2f). The largest increase (up to 13.5%) is seen off the coast of North Borneo (110°E, 5°N), associated with both increased BV-associated precipitation and some slight decreases in the seasonal total precipitation (Fig. 2d). To the southwest of Peninsular Malaysia and east of Sumatra, there is some decrease in BV precipitation and the fractional contribution to the seasonal precipitation . d-f Ensemble mean of simulated historical densities during ONDJFM (unit: number per 3.5°spherical area per season) for BV geneses, tracks and lyses respectively. g-i As (d-f) but for ensemble mean of the projected absolute changes for SSP5-85 relative to the historical period. Stippling indicates statistical significance at confidence levels of above 90% (black) and 95% (red).

Changes in BV-associated precipitation
associated with the decrease in southwestward propagating BVs over these regions, though these changes are small compared to the general increases for most of the study region. As the occurrence of BVs decreases in most of the study region, significant increases in the total amount of BV-associated precipitation are primarily driven by the significant enhancement of the mean precipitation intensity of BVs ( Supplementary Fig. 6), i.e. the mean daily precipitation rate within a geodesic distance of 3.5°of the BV centers. BV-associated precipitation has been observed in other studies to be enhanced by the occurrence of northeasterly cold surges 11,14 . However, the multi-GCM ensemble mean shows no apparent changes in the precipitation enhancement by cold surges relative to the seasonal mean BV precipitation intensity. This suggests a limited effect of the changes in cold surges on the changes in BV-associated precipitation.
During the Northeast monsoon season, a fraction of BVs can trigger extreme precipitation over Borneo 1 and the Peninsular Malaysia-Sumatra-Java region 25 . Hence, we analyze projected changes in composite precipitation averaged over 50 of the most intense BVs (ranked by their central relative vorticity at 850-hPa) in each multi-GCM ensemble member. The ensemble mean projects a significant increase in precipitation, by up to 12-16 mm per day to the north (about 1°of geodesic distance) of the BV center (Fig. 3a). Such an intensified near-center precipitation is projected by all the selected GCMs except EC-Earth3P. This can be explained by the projected intensification of extreme BV dynamics, in terms of significant increases in the low-level wind speed (g of Supplementary Fig. 7) driven by the sharper horizontal gradient of mean-sea-level pressure (h of Supplementary Fig. 7) and lowlevel convergence (i of Supplementary Fig. 7). For the changes in the near-center monsoon flow in the northern vortex sector, intensification of the meridional wind component is much more pronounced than that of the zonal component, which indicates a stronger inward radial flow that favors the near-center convergence. These structural changes also partly explain the increased mean precipitation intensity of BVs ( Supplementary Fig. 6). Analysis of changes in the probability density function for BV intensities indicates an increase of the mean BV intensity (i.e. increases in the proportion in the relatively high intensity range) when BV intensities are measured by the central relative vorticity at 850-hPa ( Supplementary Fig. 8), which is consistent with the above-mentioned intensification of extreme BVs.
To quantify extreme precipitation from BVs, we compute the number of Heavy Rain Days (HRDs) occurring within a 3.5°g eodesic radius from the identified BV center. A HRD (also referred to as R50mm 26 ) is defined as a day with daily precipitation >50 mm, a threshold suggested by the World Meteorological Organization 27 . We found a significant increase (by up to 0.8 days per season) in BV-associated extreme precipitation over most of the southern SCS during ONDJFM. An increase (up to above 12%) in the fractional contribution of BVs to the total HRDs (all days with precipitation >50 mm per day) during ONDJFM is also seen near the offshore region of North Borneo. This suggests that BVs will be more closely linked to the occurrence of extreme precipitation due to near-future climate warming.

Environmental drivers of BV responses
The projected reduction in frequency and increased precipitation intensity of BVs can be explained by two environmental mechanisms: (1) changes in the Northeast monsoon flow; (2) changes in convective instability. As BVs are observed in the lower troposphere and are mainly dependent on the low-level environmental fields 8,11 , we analyze the historical distribution and the projected future changes of environmental variables at 925-hPa, including the wind field, relative vorticity, moist static instability and vertical velocity averaged during ONDJFM (Fig. 4). In Fig. 4e, a weakening of the Northeast monsoon flow is seen across 2.5-10°N, with the most significant weakening occurring upstream of the monsoon flow near the Philippines (p value < 0.05). Significant intensification of the monsoon flow is seen to the north of 10°N. This creates a less favorable environment for BVs over their most active region near the western coast of northern Borneo. A similar result is seen for the changes in the mean synoptic fields extracted 24-h prior to the genesis of the 50 most intense BVs (Supplementary Fig. 9). Note that the analysis of the synoptic fields for the most extreme BVs indicates a significant enhancement of high-pressure disturbances expanding from the Indochina Peninsula to southern China. This indicates a more active cold air mass southward extending from the continental Fig. 2 Spatial distributions of the total precipitation amount, accumulated BV-associated precipitation and the fractional contribution of BVs to total precipitation during ONDJFM. a-c Ensemble mean of historical simulations . d-f Ensemble mean of the near-future simulations under SSP5-85 (2030-2050). g-i Ensemble mean of the projected absolute changes for SSP5-85 relative to the historical period. Stippling indicates statistical significance at confidence levels of above 90% (black) and 95% (red). high-pressure system in a warmer climate (d of Supplementary  Fig. 9), consistent with the increase in the magnitude of winter temperature variability in southern China, as projected by the CMIP5 GCMs 28 . These changes lead to a sharper pressure gradient that drives a significant intensification of low-level winds across 10-20°N. This makes the Northeast monsoon flow and the associated cold surges more intense across their northern flanks and weaker across the southern flank, making it more difficult for these systems to propagate southward and trigger BVs through orographic effects 8,14 , an observation which is consistent with the decreased occurrence of BVs. This is further confirmed by the projected northward shift in the identified BV tracks during the days when cold surges occur, according to the mean track density and occurrence fraction ( Supplementary Fig. 10).
Analyses of the inertial (in terms of relative vorticity, Fig. 4f) and moist static instability (moist-adiabatic lapse rate, Fig. 4g) show significant increases in the lower-tropospheric instability (p value < 0.05). These changes favor an increase in BV-associated precipitation and extreme precipitation events. An analysis of vertical profiles shows that increased relative vorticity dominates 10°N throughout the much of the troposphere and is accompanied by increased vertical ascent (e of Supplementary Fig. 11), a result of the northward shift in the Intertropical Convergence Zone during winter 29 . The increase in moist static instability extends throughout the lower troposphere from 700 to 925-hPa (f of Supplementary Fig. 11), which is driven by the larger magnitude of the increase in specific humidity (g of Supplementary Fig. 11) while the decrease in the tropospheric lapse rate, associated with a greater magnitude of the projected warming in the upper troposphere, has an opposite effect (h of Supplementary Fig. 11). As the maintenance of BVs relates to the low-level barotropic and baroclinity instability near the southern SCS 10,11 , further analysis of the changes in the decomposed 925-hPa barotropic and baroclinic stabilities ( Supplementary Fig. 12), in terms of the barotropic (MPV 1 ) and baroclinic moist potential vorticity (MPV 2 , both defined in Methods) suggests a significant increase in lower-tropospheric barotropic instability. In contrast, baroclinic instability decreases significantly to the north of 15°N. Although this region is controlled by an intensification of the monsoon flow and associated cold surges (which act to increase the baroclinic instability through intensified cold advection), a pronounced decrease in the meridional temperature gradient (h of Supplementary Fig. 11) is projected, which tends to baroclinically stabilize the low troposphere. Furthermore, changes in largescale ascent in terms of changes in the 925-hPa vertical velocity (Fig. 4h) show a significant increase to the north of 10°N and decrease in the vicinity of the Peninsular Malaysia-Sumatra-Java region, which agrees well with the changes in the BV-associated precipitation (Figs. 2e and 3b).

DISCUSSION
We have investigated projected near-future changes in BVs and their associated precipitation. The main conclusions are summarized as follows.
(1) An ensemble of high-resolution GCMs projects an overall decrease in the occurrence of BVs under near-future humaninduced climate warming. A reduction in storm lysis density is projected over the downstream area of the main BV propagation pathway, suggesting more stationary BVs.  with orographic effects to the south of the SCS. Increased BV-associated precipitation is primarily driven by increased barotropic instability associated with a more humid lower troposphere. In contrast, the changes in cold surges and the temperature distribution in the lower troposphere play a minor role in the projected BV precipitation changes.
This study indicates that BVs may become more impactful in a warmer climate due to projected increases in the intensity of BVassociated precipitation. The intensified BV-associated precipitation is accompanied by a general increase in the low-level relative vorticity of BVs. This finding is consistent with a BV case study by Koseki et al. 8 suggesting that the occurrence of the rainfall peak follows the growth of low-level vorticity near the vortex center. Alternatively, Braesicke et al. 30 suggested that BV-associated precipitation can be intensified by the occurrence of warm inner-cores that enhance specific humidity between 800 and 200 hPa. However, such a dynamical characteristic is not investigated here due to limited availability of the mid/high-tropospheric data from the selected high-resolution GCMs. The limited data availability also poses difficulties in presenting the relationship between low-level convective instability and the maintenance of BVs during different stages of the vortex lifecycle. The analysis presented here should be supplemented in future work using appropriate model experiments in order to understand the internal mechanisms related to the changes in the vortex precipitation and the vortex interaction with the more convectively unstable environment in a warmer climate. In addition, as the GCM experiments are atmosphere-only simulations, two-way atmosphere-ocean coupling processes are not considered in this study. Although the GCM simulations consider the warming of the sea surface temperatures (SSTs) using an SST projection from the CMIP5 ensemble mean 16 , the absence of two-way coupling processes may lead to uncertainties in the simulated BVs, especially for those interacting with the winter SST cold tongue over the southern SCS [31][32][33] and those developing into tropical storms, which may be particularly dependent on the local SST pattern 34,35 .
BVs are closely associated with hydrological extremes in Southeast Asia. This study has shown evidence that anthropogenic climate warming will alter these high-impact synoptic processes in the near future, making them more impactful as they alter the frequency and intensity of rainfall. This study is expected to motivate more research on the future impact of BVs; in particular, the hydrological impacts of BVs in present and projected future climates that have not been adequately quantified in hydrologically vulnerable regions within Southeast Asia. Hence, we call for more research in these aspects of BV behavior using advanced hydrological modeling techniques, a scientific requisite for achieving better management of future storm-induced hydrological risks.

Models
This study is based on six high-resolution versions of the HighResMIP GCMs 16 . The GCM outputs are from atmosphereonly (i.e. AMIP-style experiments of Tiers 1 and 3) simulations for the present-day and future climate. Further information on the model configuration, horizontal and vertical resolutions as well as for their low-resolution counterparts in the HighResMIP ensemble are summarized in Supplementary Table 1. The horizontal resolutions of the GCM simulations correspond to a range of 10-50 km grid spacing, which is finer than most of the non-HighResMIP CMIP6 simulations. Note that the selection of the GCM experiments for this study is limited by the current availability of six-hourly outputs for pressure-level winds and mean-sea-level pressures as well as the daily output of precipitation (required for the diagnoses of BV characteristics and the associated precipitation) in the HighResMIP data repository. The  (u, v). b, f Relative vorticity (ζ p ). c, g Moist static instability (i.e. the vertical gradient of pseudo-equivalent potential temperature, δθ se /δp). d, h Vertical velocity (ω). a-d Ensemble mean of historical simulations . e-h Ensemble mean of the projected absolute changes for SSP5-85 (2030-2050) relative to the historical period. Stippling indicates statistical significance at confidence levels of above 90% (black) and 95% (green).
selected HighResMIP GCMs were the only models satisfying this requirement at the time this study was undertaken. Outputs from the coupled atmosphere-ocean (Tier 2, named hist-1950 and highres-future) experiments are also provided by HighResMIP; however, they are not used here due to the limited availability of all required six-hourly variables.
The present-day climate simulations cover the period 1979-2014, and are forced by time-dependent SSTs and sea-ice (named highresSST-present) from the HadISST2.2.0 observational dataset 36 . These experiments start on 1st January 1950 using the initial atmospheric and land-surface conditions from the ERA-20C reanalysis data 37  It is noted that some HighResMIP Tier-1 experiments provide future climate simulations of a longer period till the end of the 21st century, while the data required for the detection and tracking of BVs is currently not available. Information on the applied external climate forcings, such as the anthropogenic aerosol forcing, volcanic stratospheric aerosol, natural aerosol forcing and solar forcing, can be found in the HighResMIP description paper 16 .
This paper is based on high-resolution GCM simulations instead of regional climate model (RCM) simulations given the authors' limited data storage capacity and the known issues of RCMs, including the combined model uncertainties from the RCM and the driving models providing lateral boundary conditions 39 , uncertainty related to the choice of model domain 40 and the lack of two-way interactions between the simulated regional climate and the external global climate 41 . However, as RCMs can run at relatively high spatial resolution with relatively low computational costs compared to GCMs at equivalent spatial resolutions, it is worth using RCM simulations (e.g. Coordinated Regional Climate Downscaling Experiment-Southeast Asia 42 ) to investigate the projected future changes in BVs and the uncertainties related to RCM configurations in future work.

Detection and tracking of BVs
In the ERA5 reanalysis data and the HighResMIP simulations, BVs are detected and tracked using the objective tracking algorithm (TRACK) developed by Hodges 20,21,43 , which is a synthesized program combining the tracking of targeted synoptic processes with statistical analyses of characteristics of the targets. The effectiveness of TRACK has been demonstrated in BV analysis by Liang et al. 1 as well as in studies of extratropical cyclones [44][45][46][47][48] , tropical cyclones 44,[49][50][51][52][53] and tropical easterly waves 54,55 .
The procedures for the detection and tracking of BVs are as Liang et al. 1 . First, a spectral filter provided by TRACK is applied to the 6-hourly 850-hPa relative vorticity fields to remove signals with total wavenumbers less than 5 and truncate the field wavenumbers to T63, using spherical harmonics and spectral tapering of the spectral coefficients based on the method of Hoskins and Sardeshmukh 56 . This is to obtain more coherent BV tracks by removing the non-BV high-frequency noise and the lowfrequency background in the vorticity fields. In the filtered fields, features are initially identified as grid point maxima that exceed 5 × 10 −6 s −1 in the tropics. These are then refined by searching for the off-grid vorticity maxima using B-spline interpolation and steepest ascent maximization as described by Hodges 20 . The tracking of BVs is then performed by initially finding the nearest neighbors and then using minimization of a cost function for track smoothness 20,57 . As BVs are quasi-stationary systems near Borneo 3 , we filter out the non-BV disturbances by removing the tracks that have a pressure minimum for less than one day within a targeted area (2.5°S-7.5°N, 102.5°-117.5°E). The genesis and lysis locations of BVs are then defined as the first and last feature points of the identified tracks.
The intensities (in terms of the local maximum 850-hPa relative vorticity) of BVs are determined following the tracks. The minimum mean-sea-level pressure near the detected BV centers is also identified using B-spline interpolation and steepest descent for searching the minimum within 4°of the detected center. The precipitation associated with BVs is accumulated on the model grids for all precipitation occurring within a 3.5°geodesic radius of the identified BV centers, according to Liang et al. 1 . Then, the fractional contribution of BVs to the total precipitation during ONDJFM is calculated. The search radius of 3.5°is approximately 400 km, covering the intense precipitation near the BV center as well as the broader precipitation footprint typically observed during the developing and decaying stages 8 . To examine simulations of BV-related extreme precipitation, on each model grid the occurrences of extreme precipitation for daily precipitation greater than 50 mm (e.g. heavy rain days, HRDs) is determined. Then we define the BV-associated extreme precipitation as any HRDs occurring anywhere within a 3.5°geodesic radius of the identified BV centers. Based on this, we further compute the fractional contributions of BVs to the seasonal total of HRDs during ONDJFM.
Determination of cold surge occurrence in the SCS Analyses of projected changes in the mean precipitation intensity and track distribution of BVs during the occurrence of the northeasterly cold surges in the SCS are shown in Supplementary  Figs 6 and 10. The occurrence of cold surges are defined as the days when the daily northeasterly wind speed at 850-hPa regionally averaged within the region of 5°-10°N, 107°-115°E exceeds 0.75 standard deviations above its climatological mean, according to Lim et al. 15 and Xavier et al. 58 .

Analysis of BV-associated large-scale environment
The analysis is based on monthly fields of horizontal wind field (u, v), relative vorticity (ζ p ), moist static instability in terms of the vertical gradient of pseudo-equivalent potential temperature (δθ se /δp) and vertical velocity (ω). As in Liang et al. 1 , the averaged fields (925-hPa winds and mean-sea-level pressure) for the timesteps 24 hours prior to the genesis phases of the 50 most intense BVs are analyzed to show the synoptic environments which lead to the formation of BVs, e.g. the southward extending high-pressure disturbances and northeasterly cold surges 11 . In addition, to further understand the changes in the atmospheric instability, in Supplementary Fig. 12 we analyze the large-scale low-level barotropic stability (to consider the changes in inertial and moist static instabilities) and baroclinic stability (to consider the changes in the horizontal thermal advection by the cold surges) in terms of the barotropic (MPV 1 ) and baroclinic moist potential vorticity (MPV 2 ) at 925-hPa. These are defined by Qian et al. 59 as: where g is the acceleration of gravity, f is the planetary vorticity, ∂θse ∂x and ∂θse ∂y are respectively the zonal and meridional gradients of pseudo-equivalent potential temperature. ∂v ∂p and ∂u ∂p are the vertical shears of the meridional winds and zonal winds respectively.

Statistical significance tests
Statistical significance tests were conducted to evaluate the difference between the samples between the present-day and future climate periods. These are based on the p value of the Wilcoxon rank-sum test 60 (also referred to as the Mann-Whitney U test 61 ), which indicates the probability of accepting the null hypothesis that the present-day and future climate samples are from the same population.

DATA AVAILABILITY
The CMIP6 HighResMIP data are downloaded from the data node website of the Lawrence Livermore National Laboratory (https://esgf-node.llnl.gov/projects/cmip6/). The ERA5 climate reanalysis datasets can be downloaded from the website of the Copernicus Programme (https://cds.climate.copernicus.eu/cdsapp#!/dataset/ reanalysis-era5-pressure-levels?tab=overview). The TRACK outputs for the identified BV features based on the datasets above are available upon request from the corresponding author Ju Liang: J.Liang@exeter.ac.uk.