Extreme rainfall in Northern China in September 2021 tied to air–sea multi-factors

The September rainfall over Northern China (NC) in 2021 was the heaviest since 1961 and had unprecedented socioeconomic impacts. Holding the hypothesis that the drivers of extreme climate events usually contain extreme factors, we firstly propose the Ranking Attribution Method (RAM) to find the possible air–sea multi-factors responsible for this rainfall event. Via the atmospheric bridges of zonal-vertical circulation and Rossby wave energy propagation, the remote factors of warm sea surface temperature anomalies (SSTA) over the tropical Atlantic, cold SSTA over the tropical Pacific, Southern Annular Mode-like pattern in the Southern Hemisphere and North Pacific Oscillation-like pattern in the Northern Hemisphere jointly strengthened the Maritime Continent (MC) convection and Indian monsoon (IM). Through meridional-vertical circulation, the intensified MC convection enhanced the subtropical high over southern China and induced ascending motion over NC. The local factor of extreme air acceleration in the east Asian upper-level jet entrance region further anchored the location of the southwest-northeast rain belt. The strengthened IM and subtropical high over southern China induced considerable moisture transport to the rain belt via two moisture channels. The combined effect of these extreme dynamic and moisture conditions formed this unprecedented rainfall event. This study suggests that the RAM can effectively reveal the factors that contributed to this extreme rainfall event, which could provide a new pathway for a better understanding of extreme climate events.


Introduction
Autumn rainfall events in China contribute most strongly to the annual mean rainfall in the northern area, and the autumn rainfall in central China has a decadal transition from a decreasing trend to an increasing trend in the late 1990s (Liu and Fan 2013;Wang and Zhou 2019). Extreme Northern China (NC) rainfall events in autumn have been shown to result in great socio-economic impacts (Xu et al. 2016). In the early autumn of 2021, 1.76 million residents suffered from the extreme rainfall event in Shanxi Province, with the collapse of more than 17,000 houses and about 190,000 hectares of crops damaged (XINHUA). Moreover, the continuous downpours damaged the ancient town of Pingyao, which is a world cultural heritage site with a history of 2700 years (XINHUA). The severe consequences of this extreme rainfall event urged us to investigate its inner mechanism.
Understanding the processes that lead to extreme events is a key goal of climate science (Stephenson et al. 2008). Previous studies regarded the El Niño-Southern Oscillation (ENSO) as an important factor of autumn rainfall anomaly in NC (Li and Zeng 2013;Liu and Fan 2013;Xu et al. 2016). Extreme autumn rainfall days in NC are more frequent in La Niña phases than in El Niño phases during 1951-2004(Li et al. 2011. Recent studies analyzed specific extreme autumn rainfall events in NC at synoptic timescale (Zhao et al. 2016;Li et al. 2018). Less attention has been paid to thoroughly investigating a specific monthly-to-seasonal extreme autumn rainfall event in NC. Therefore, this study focuses on the attribution of NC extreme early-autumn (September) rainfall in 2021.
Attribution refers to tracing the causation (NAS 2016). Event attribution is an integral component of climate services to inform adaptation and mitigation programs around the world and to support climate risk management (Stott et al. 2016). Attribution of a specific extreme event is a pressing challenge for the climate science community (WRCP 2021). Generally, a specific extreme climate event can be attributed to multiple causes which usually contain extreme factors. For instance, there have been mixed conclusions regarding the severe air-sea causes of record-breaking June rainfall over the Yangtze River Valley in 2020 recently. The causes include the significant quasi-biweekly oscillation of East Asian monsoon circulation system (including the western Pacific subtropical high, the upper-level East Asian westerly jet, and the low-level southwesterly jet) (Ding et al. 2021), North Atlantic Oscillation sub-seasonal phase transition (Liu et al. 2020), exceptionally persistent Madden-Julian Oscillation activity ), excessive atmospheric blockings over East Siberia , persistent positive Pacific-Japan pattern, rapid development of La Niña (Qiao et al. 2021), extremely strong Indian Ocean Dipole event (Takaya et al. 2020;Zhou et al. 2021), and positive sea surface temperature (SST) anomalies (SSTA) in Atlantic (Zheng and Wang 2021).
Finding possible factors of an extreme climate event is not only the precondition of understanding its formation, but also the first-step towards physical verification through model simulations and experiments. Correlation analysis is a frequently-used method to find relavant factors, with the correlation coefficient representing general relationships between the rainfall amount and relevant factors. Since the causality of each extreme event is unique, correlation analysis is not sufficient to figure out the causes of an individual record-breaking climate event. According to previous studies, we hold a relatively justified assumption that serious consequences are usually associated with severe control factors. For instance, super El Niño events gave rise to severe summer floods in China (Li 1999;Yuan et al. 2017), record-breaking warming and extreme drought in the Amazon rainforest (Jiménez-Muñoz et al. 2016); the strongest Australian high staged the record-breaking rainfall over South China in June 2017 (Sun et al. 2018); and recordbreaking northward shift of the western Pacific subtropical high led to record-breaking heat event over Northeast Asia in summer 2018 (Ren et al. 2020). Since severe factors occupy special historical statues, and ranks can not only act as the direct measurements of anomalies but also represent their relative statuses in history. We firstly proposed the Ranking Attribution Method (RAM) to make an attempt for objectively searching for possible responsible air-sea factors associated with this extreme September rainfall event in NC. To exclude the high (low) rankings that just simultaneously occur without a specific causal relation, we established the physical linkages between the NC rainfall in September 2021 and the potential contributors, which may provide responsible causality.
In Sect. 2, we describe the data used in this study, RAM, calculations of significance test between the means of two samples where one sample contains only one value, the air acceleration in the jet at monthly time scale, the unified Dynamic Normalized Seasonality (DNS), Indian monsoon (IM) index, and the Rossby wave ray trajectory. Section 3 depictes the record-breaking rainfall in NC. Section 4 investigates the combined effect of record-breaking Maritime Continent (MC) convection, Indian monsoon (IM) and upper-level air acceleration on this rainfall event. Section 5 explains the possible causes of extreme MC convection and IM anomalies. Section 6 is the conclusion and discussion.

Data
We obtain the 699 in-situ rainfall observations covering the period 1961-2021 from the National Meteorological Information Center of China. Monthly reanalyzed rainfall, surface air temperature (SAT), SST, sea level pressure, horizontal wind, vertical velocity, geopotential hight (GHT), and vertical integrated moisture flux in 1979-2021 are taken from the ERA5 reanalysis dataset on a 0.25° × 0.25° grid (Copernicus Climate Change Service 2017; Hersbach et al. 2018Hersbach et al. , 2019. According to the data length, September anomalies of observed rainfall (reanalyzed data) are derived relative to the climatological mean of 1961-2021 (1979-2021). The NC rainfall in this study refers to the rainfall within the southwest-northeast rain belt (marked by the blue contour in Fig. 1).

Ranking attribution method (RAM)
According to the hypothesis that extreme consequences are usually associated with extreme control factors and the fact that extreme factors occupy special historical statuses, the RAM employs top/bottom ranks to search for possible extreme factors. For an individual extreme event Y and the year T when Y occurs, the RAM obtains the historical statue of a relevant variable X at the year T by calculating the rank R T(Y) is a time series for ranking, n is the number of samples (number of years). Positive ranks represent the ranking from the maximum and negative ranks from the minimum, with 1 referring to the historical high and − 1 the historical low. Detrending is not needed for X before ranking since the extreme event Y is defined by the raw data, and detrending may change its historical statuses, which leads to mismatched results with the extreme event Y. The result based on the RAM does not depend on the selection of climatological baseline period. Besides extreme climate events, the RAM can also be used for analyzing extreme weather events. In this study, we obtain the historical status of September rainfall in 2021 on every spatial grid [named as R 2021 (Rainfall) for short] through the RAM to identify the specific region of extreme rainfall in NC. Since it is critical to understand the multiple air-sea causes of an extreme rainfall event (Otto et al. 2015), the RAM is used to provide the rankings (R 2021 (X)) of relevant air-sea variables X on every spatial grid. Possible extreme factors, such as striking anomalies of SST and atmospheric circulations, are nominated by the top (bottom) rankings. Finally, responsible causality is provided by establishing the physical connections between possible extreme factors and drastic vertical motion as well as sufficient water vapor transport which are directly associated with the NC extreme rainfall. The ranking of observed rainfall (reanalysis data) in September 2021 is calculated from 1961 (1979) onward.

Significance test
To rigorously test whether the September statistics in 2021 are indeed extreme and significantly different from those in other years, we set a null hypothesis that the variable being tested (X 0 ) has the same true mean as those in other years (X 1 , …, X n ). Let then t has a distribution of Student's t with n − 1 degrees of freedom. This is a special case of significant test between the means of two samples where one sample contains only one value (Trenberth 1984). The statistical significance at the 95% confidence level are marked by white dots in this study.

Air acceleration in the jet entrance region
In the jet entrance region of upper-level atmosphere, the equation of motion for the zonal wind (Rochette and Market 2006) on synoptic time scale (neglecting the viscosity term) is where u, v, and z are the zonal wind, meridional wind, and geopotential height at 200 hPa, respectively, with the subscript p indicating on the isobaric surface; f is the Coriolis force; v g and v ag are meridional geostrophic wind and geostrophic departure; g is the gravitational acceleration. Equation (4) exhibits the relationship between the meridional ageostrophic wind and zonal air acceleration in the jet entrance region.
To quantify the air acceleration in the East Asian jet at monthly time scale, we promote Eq. (4) into monthly (or seasonally) time scale by decomposing the variables into monthly-mean term and disturbance term and obtaining: where the variables with the bar are the mean, and the ones with the prime ′ are the disturbance. For the large-scale motion, the scales of the variables in Eq. (5) are as follows: After omitting the terms of smaller magnitudes, Eq. (5) turns into: Thus, the relationship between the meridional ageostrophic wind and air acceleration also works at monthly (or seasonally) time scale.

The unified dynamic normalized seasonality
The DNS has been widely practiced to depict the interannual variability of monsoons (Li and Zeng 2000, 2003, 2005, and the calculation is based on the seasonal alternation in wind direction over monsoon regions. For a given pressure level, the September DNS index in a year on the grid point is as follows: where ��� ⃗ V 9 is the wind vector of September, ��� ⃗ V 1 and ��� ⃗ V 7 are the climatological wind vectors of January and July, respectively. The norm A is ( ∬ |A| 2 dS) 1∕2 where S denotes the i n t e g r a t i o n d o m a i n .
is the calculation at a point, with a and j indicating the mean radius of the earth and latitude of the point, respectively. Δ and Δ in radians are resolutions in the meridional and zonal directions. We employ the area-mean (pink box in Fig. 5c) DNS index over the India Peninsula to indicate the IM index (IMI). The dynamic IMI is also calculated by the zonal wind difference between a northern region (20°N-30°N, 70°E-90°E) and a southern region (5°N-15°N, 40°E-80°E) around India Peninsula at 850 hPa (Wang et al. 2001, blue boxes in Fig. 5c).

Rossby wave ray tracing
The wave ray trajectory is based on the dispersion relation of the barotropic nondivergent vorticity equation on a time-mean slowly varying basic state. The stationary Rossby wave ray method Zhao et al. 2019) is employed to depict the wave energy propagation in the mean flow of September 2021. On a horizontally non-uniform flow, the dispersion relation of a barotropic Rossby wave is  where K = √ k 2 + l 2 is the total wave number, with the k and l representing the zonal and meridional wave numbers, respectively. A Mercator projection of the sphere is used for a more convenient analysis (Hoskins and Karoly 1981). u M , v M = u, v ∕ cos is the Mercator projection of horizontal winds, where is the latitude. is the frequency and q = 2Ω sin + ∇ 2 is the absolute vorticity. The Ω and are the earth rotation rate and the stream function of the basic-state. The local group velocity c g = u g , v g are: represents the material derivative moving with the group velocity. Determined by the kinematic wave theory (Whitham 1960;Shaman et al. 2012), the wavenumbers k and l vary with the position along the way ray: Equations (9a, b) and (10a, b) provide the complete set of equations for Rossby wave ray tracing for stationary waves ( = 0 ). For a given starting point and an initial zonal wave number k, the initial local meridional wave number l is determined by the dispersion relation (Eq. 8). Thus, numerically integrate the above equations can derive the Rossby wave ray trajectory. The integrations were terminated to omit the small-scale Rossby waves once the local meridional wave length was less than 1000 km.

Record-breaking rainfall in Northern China
The gauge observation in September 2021 depicts a distinct southwest-northeast rain belt over NC (Fig. 1a, b). Accumulated September rainfall at many stations significantly exceeded their 1961-2021 mean values by more than 300% or 300 mm. Meanwhile, southern China was drier than the climatological mean. The ERA5 rainfall data reliably reflects the spatial distribution of the anomalous rainfall compared with the gauge observations in September 2021 (Fig. 1c, d). Consequently, we mark the rain belt over the land by the 50 mm contour of excess rainfall from the ERA5 data, and this also indicates the NC domain where the rainfall anumalies are statistically significant at 95% confidence level. The interannual variation of the September rainfall anomalies in NC from ERA5 is also closely correlated with that of gauge observations, with the correlation coefficient reaching 0.95 (Fig. 2a). The accumulated rainfall anomaly in September 2021 in NC is 140 mm (equivalent to 5.5 SDs), which breaks the record that has stood since 1961. We use the RAM to present the extreme September rainfall in 2021 through the gauge observation, and find that rainfall data from 74 (165) stations in NC rank in the first place (top five) since 1961 (Fig. 2b). It is notable that the RAM can not only depict the extreme heavy rainfall in NC but also reveal the minimal September rainfall in southern China in 2021. The RAM output based on the ERA5 rainfall resembles that of observed data (Fig. 3a). Ranking of ERA5 rainfall also shows that aside from NC, there are recordbreaking September rainfall over the MC, west of the Indian Peninsula, south of the Indo-China Peninsula, and to the east of Baikal in 2021 since 1979.
Via the spatial distribution of the ranks that locate near the top or bottom, the RAM can detailly describe the spatial pattern of extreme situation. Then, we employ the RAM to search for possible extreme signals that may contribute to this record-breaking rainfall event from the simultaneous ranking conditions of relevant fields and consider possible mechanisms as well.

Dynamic effects of extreme maritime continent convection and upper-level jet stream
The ranking of the vertical velocity ( Fig. 3b) shows that all of the above-mentioned record-breaking rainfall regions ( Fig. 3a) experienced extremely strong upward motion at 500 hPa. The significantly vigorous convection activity over the MC played an important role in stimulating the anomalous ascending motion over the rain belt in NC via the atmospheric bridge (Lau and Nath 1996;Li et al. 2019aLi et al. , b, 2022. The zonal-averaged meridional circulation between 110°E and 125°E (Fig. 3c) shows that the convection over the MC favored strengthened regional Hadley cells (Lau and Chan 1983a, b) with anomalous downward motion between 15°N and 25°N, which enhanced the subtropical high over southern China and reduced the rainfall there (Fig. 1a). The strengthened subtropical high further prompted the ascending motion spanning 35°N-40°N through significant southerly anomalies in the lower troposphere, which dynamically nurtured the formation of the rain belt. However, distinct southerly anomalies also appear in the upper troposphere from 25°N to 60°N, which may be linked with the anomalous upper-level jet stream and distribution of trough and ridge.
Apart from the meridional circulation originating from the MC, the local factor of East Asian upper-level jet stream also contributed to upward motion over NC rainfall during September 2021. There are strong associations between China rainfall and the upper-level jet at seasonal-to-interannual scales (Liang and Wang 1998). First, the air acceleration in the upper-level jet entrance region is related to the mass confluence and upward motion, and dominated by the ageostrophic motion (Namias and Clapp 1949;Reiter 1963;Blackmon et al. 1977;Keyser and Shapiro 1986). The crossjet ageostrophic circulation determines the rainfall location relative to the jet (Liang and Wang 1998). Thus, we calculate the air acceleration in September 2021 using the meridional ageostrophic wind data at 200 hPa based on the monthlyscale equation of motion (Eq. 6). The air acceleration in the entrance region of the jet reaches its highest value along the west flank of the rain belt since 1979 (Fig. 4a). This favored the southwest-northeast record-breaking divergence at the upper atmosphere (Fig. 4b), which further facilitated the convergence at the lower atmosphere and formed the southwest-northeast rain belt. Second, a poleward shift of the jet is also associated with the heavy rainfall (Yeh 1959;Trenberth and Guillemot 1996). Anomalous westerlies (easterlies) appeared around the north (south) of the rain belt at 200 hPa, and the zonal wind speed in September 2021 around the south of the rain belt was the lowest since 1979 (Fig. 4c). This reversal of the winds indicates a northerly shift of the jet. The anticyclonic circulation induced by the northward jet shift contributed to the anomalous upper-level southerly anomalies around the rain belt (Figs. 3c, 4d) and Ranking of a September rainfall R 2021 ( ) (shaded) and b September 500 hPa vertical velocity R 2021 − 500 (shaded) since 1979. c Regionally anomalous vertical velocity (shaded, − 10 -2 Pa s −1 ) by averaging − between 110°E and 125°E, with right and up arrows indicating southerly and upward motion. The solid (dashed) white contour indicates the 95% (90%) confidence level. The gray vectors represent the anomalous meridional-vertical circulation (m s -1 ; − 30 Pa s −1 ), and the black ones are above the 90% confidence level was conducive to the development of cyclones at lower atmosphere. Thus, the upper-level anticyclonic circulation favored the formation of rainfall by inducing the lower-level cyclonic circulation, and the upper-level air acceleration further anchored the southwest-northeast rain belt (Fig. 4d).
In short, from the perspective of dynamic condition, the record-breaking convection over the MC, significantly strengthened air acceleration in the East Asian upper-level jet entrance region together with the jet shift contributed to this extreme rainfall event. We then focus on the moisture transport which is another essential condition (Westra et al. 2014) for extreme rainfall events.

Moisture transport of extreme Indian monsoon
In the climatology for September (vectors in Fig. 5a), the vertical-integrated water vapor transported by the easterly wind on the southern flank of the subtropical high and the westerly wind associated with the IM converges over the Indo-China Peninsula at first. Then it extends along the northern flank of the subtropical high to feed the rainfall over NC. This process was distinctly strengthened in 2021 (vectors in Fig. 5b). In addition, ranking of the eastward vapor flux (color in Fig. 5b, with positive and negative values indicating extreme eastward and westward moisture transport) in September 2021 reflects that the vapor transported by both the easterly and westerly winds to the Indo-China Peninsula was the most. Ranking of the northward vapor flux (color in Fig. 5a, with positive and negative values indicating extreme northward and southward moisture transport) shows that abundant vapor was further transported northward to supply the rain belt over NC (Fig. 5b). According to the condition of extreme moisture transport, the above two processes may result from the enhanced subtropical high over southern China (Fig. 3c) and record-breaking IM.
However, is the IM in September 2021 record-breaking in reality? We obtain the spatial distribution of IM intensity ranks at 850 hPa (color in Fig. 5c) through the unified Dynamic Normalized Seasonality (DNS) index (Eq. 7; Li and Zeng 2000, 2003, 2005. The DNS indexes over northern and southern India reached their highest values since 1979, resulting from the anomalous northwestward and northeastward flows over the north India and the Arabian Sea, respectively (vectors in Fig. 5c). The IM index (Fig. 5d) calculated by the zonal wind difference at 850 hPa (blue boxes in Fig. 5c; Wang et al. 2001) is highly correlated with the area-mean DNS index over the India (deeppink boxe in Fig. 5c). Both the two IM indices show a strongest September IM in 2021.
Besides the aforementioned pathway of moisture delivery, the anomaly and ranking of vapor flux reveals another water vapor channel. Induced by the extreme IM, the vapor from the South China Sea and the northern Indian Ocean was firstly delivered northwestward along the southern edge of Tibetan Plateau to the Iranian plateau (color in Fig. 5b). Then the vapor over the Iranian plateau keeps going northward until the anomalous westerlies turned it to the east, causing extremely eastward vapor delivery over the Mongolia (Fig. 5a). This moisture channel also fed the rainfall in NC. The anomalous westerlies may have been related to the northerly shift of the jet (Fig. 4c). Thus, the joint impact of the extreme IM and the enhanced subtropical high over southern China delivered sufficient water vapor to the rain belt.
As mentioned earlier, record-breaking September rainfall also occurred over the MC, the western Indian Peninsula, south of the Indo-China Peninsula, and east of Baikal in 2021 (Fig. 3a). The rainfall over the MC can be interpreted by the extreme local convection. Rainfall over the western Indian Peninsula and south of the Indo-China Peninsula was related to the strong IM, which transported sufficient water vapor to feed the rainfall there. This suggests that the record-breaking MC convection and IM contribute to the dramatic regional rainfall while inducing the distant rainfall in NC.

Causes of record-breaking maritime continent convection and Indian monsoon
The above analysis of the dynamic and moisture conditions for the rainfall in September 2021 implies that the extreme MC convection and IM were essential. At interannual timescale, the variation of MC convection is primarily controlled by the ENSO (Philander 1983;Hamada et al. 2002;Hendon 2003). A negative SSTA in the central to eastern tropical Pacific could strengthen the zonal SST gradient over the Pacific and enhance the Walker cell, which leaves the MC region controlled by anomalous ascending motion (Lau and Chan 1983a, b). The interannual variation of the IM intensity is also significantly linked with the simultaneous ENSO (Rasmusson and Carpenter 1982;Webster and Yang 1992). In general, a stronger (weaker) Indian summer monsoon is associated with a cooler (warmer) phase of ENSO, resulting from the upward (downward) branch of anomalous Walker circulation (Webster and Palmer 1997).

Sea surface temperature anomalies
To identify the causes of the strong MC convection and IM in September 2021, we firstly examine the simultaneous SST and SAT anomalies. The anomalous SST and SAT over the central to eastern tropical Pacific are negative but not significant (Fig. 6a, b). Although the cold SSTA of the tropical Pacific favors the enhancement of MC convection and IM intensity, this seems to be insufficient on its own. Thus, the record-breaking MC convection and IM may have been generated by other factors rather than the weak La Niña. Meanwhile, there are significantly warm regions over the tropical Atlantic, tropical Indian ocean, mid-latitude Pacific, and mid-latitude Atlantic (Fig. 6a, b). The high ranking of SST and SAT also successfully identifies these regions (Fig. 6c, d).
The extremely warm SSTA of the tropical Atlantic may have played an important role in influencing the MC convection and IM in September 2021. The tropical Atlantic SST connects the tropical Pacific SST and the atmospheric circulation above via the inter-basin SST gradient (Wang 2006). Warm tropical Atlantic and cold tropical Pacific SSTAs increase the zonal SST gradient, which results in a zonal-vertical circulation with an ascending branch over the Atlantic and a descending branch over the central to eastern Pacific (Rodríguez-Fonseca et al. 2009). Although the negative SSTA over the tropical Pacific in September 2021 was not extreme, the strengthened descending motion induced by a warm tropical Atlantic could have enhanced MC convection and IM intensity via the atmospheric bridge of anomalous Walker cell over the Pacific. The record-breaking warm SST over the tropical Indian ocean helped to maintain a stronger subtropical high over southern China through equatorial Kelvin waves and the resultant divergence in the subtropics (Xie et al. 2009(Xie et al. , 2016. Moreover, the SST around the MC ranks the top five, which also helped to facilitate the development of the convection (Fig. 6c).
Besides the tropical SSTA, the SSTA over the mid-latitude Pacific in both hemispheres also broke their record since 1979. The SSTA over the North Pacific may be associated with the negative phase of the Pacific Decadal Oscillation, which is connected with more extreme rainfall events in NC at inter-decadal time scale (Wei et al. 2021). The SSTA over the South Pacific may be related to the water transport induced by abnormal atmospheric circulation.

Atmospheric circulation anomalies
Atmospheric circulation anomalies also contributed to the MC convection and IM. The anomalous atmospheric circulation pattern (Fig. 7) exhibits an equivalent barotropic structure over the middle-high latitudes, with a significantly positive North Pacific Oscillation (NPO)-like pattern in the Northern Hemisphere and a significantly positive Southern Annular Mode (SAM)-like pattern in the Southern Hemisphere. Ranking of the GHT at 500 hPa and 850 hPa highlights the particularly distinct anomalous subtropical high over the northern and southern Pacific (Fig. 8a, b). Besides the Pacific, the GHT over the Mascarene High, the Australian high and the subtropical high over the southwestern South Atlantic ranks in the top five at both 500 hPa and 850 hPa (Fig. 8a, b). It is worth noting that the extremely high GHT around the eastern tropical Pacific suggests the downward motion and low-level wind divergence there. This probably originated from the upward motion over the warm tropical Atlantic. Furthermore, the cold SSTA in the tropical Pacific also facilitated this process. In addition, the recordbreaking negative signal over the Arabian Sea is evident only at 850 hPa, and this is associated with the extreme IM.
Possible mechanisms of the SAM effects on the MC convection and IM have been explored in previous studies. A positive SAM is concurrent with strong MC convection via the meridional circulation along the central South Pacific and two meridional teleconnection wave train patterns, with one over the southern Indian Ocean at the lower level and the other along the central South Pacific at the upper level (Sun et al. 2009). The positive SAM favors the development of IM by inducing the enhanced westerly flow and anomalous cyclone at 850 hPa near the Indian Peninsula (Viswambharan and Mohanakumar 2012).

Linkages between maritime continent convevtion, Indian monsoon and extreme factors
To dynamically link the MC convection and IM with these possible air-sea factors through atmospheric bridges, we use the Rossby wave ray tracing theory (Eqs. 9 and 10) to describe the wave energy propagation at higher (200 hPa) and lower (850 hPa) atmosphere Zhao et al. 2015). The raw monthly-mean horizontal wind of September 2021 is employed as the background state to calculate the wave ray trajectories, which are tangential to the group velocities (Lighthill 1978). Based on the regions with extreme tropical SST (Fig. 6c) and GHT (Fig. 8a, b), we set the initial points at the tropical Atlantic, the equatorward lobes of the SAM-like and NPO-like patterns, respectively (Fig. 8c, d). Since that the effect of ENSO on MC convection and IM is generally acknowledged, we also set initial points at the negative SSTA of the eastern tropical Pacific. The wave ray trajectories at 200 hPa (Fig. 8c) show that the wave energy from the tropical Atlantic, tropical eastern Pacific, southern lobe of the NPO-like pattern and mid-latitude southern Pacific can propagate to the MC region (the right black box) and influence the convection there. At 850 hPa, except for the southern lobe of the NPO-like pattern and south of Australia, the wave energy from all the initial regions can propagate to the IM region (the left black box, Fig. 8d) and affect the IM intensity. Moreover, the wave ray from the tropical Pacific, the southern Pacific, and the southeastern South Atlantic pass through the MC region before arriving the IM region. This suggests that the wave energy could fasten the upward motion over the MC before enhancing the IM intensity. The wave rays from south of the Australia also reach the MC region. We then set the wave sources over the MC region at 850 hPa, and the wave energy from the MC propagates directly to the IM region. In short, from the perspective of wave energy propagation, the MC convection was influenced by multiple air-sea factors from both the upper and lower atmosphere, while the IM can not only be directly affected by these factors but also indirectly influenced by the MC convection at lower atmosphere.

Conclusion and discussion
Based on a relatively reasonable hypothesis that drastic consequences are usually associated with drastic factors, we propose the RAM to search for possible factors that contributed to the heaviest September rainfall over NC in 2021. We find several concurrent record-breaking factors which may dominate this event, and establish the physical linkages between multi-factors and the rainfall (Fig. 9). First, these air-sea factors intensified the MC convection and IM through the atmospheric bridges of anomalous zonal-vertical circulations and Rossby wave propagation. Then, the MC convection enhanced the subtropical high over southern China and gation from the mid-latitude northern Pacific, mid-latitude southern Pacific, tropical Atlantic and tropical Pacific to the MC region, respectively. At 850 hPa, the pink, green, purple, blue, and orange lines are the wave rays from the tropical Pacific, tropical Atlantic, southern India Ocean, mid-latitude southern Pacific, and southwestern South Atlantic to the IM region. Cyan lines are the wave rays from the southern Australia to the MC, and the red lines are those from the MC to the IM region induced ascending motion over NC. The East Asian upperlevel jet stream anchored the southwest-northeast rain belt in NC by significantly strengthened air acceleration and northward jet shift. The joint effect of the southern China subtropical high and extreme IM induced considerable moisture transport to the rain belt via two moisture channels. Finally, the combined effect of extreme dynamic and moisture conditions formed this unprecedented rainfall over NC. The reason of extreme air acceleration in the jet entrance region needs further investigation. We also did not calculate the relative contributions of these factors, and the method of statistically quantifying the effects of multi-factors on a single event still needs development. Further idealized modeling studies are needed to reveal the relative roles of these factors. While it is challenging to design the experiment for quantifying the effects of atmospheric forcing and multiair-sea factors. Suitable experiment should be designed in the future.
Additionally, the heavy rainfall east of the Baikal (Fig. 3a) locates between the low trough and the high ridge, where the upward motion ( Fig. 3b and color in Fig. 3c) and northward moisture transport (color in Fig. 5a and vectors between 40°N and 60°N in Fig. 3c) are distinct. This could be interpreted by the anomalous SAT over Eurasia, which is related to the anomalous GHT in the upper troposphere. The SAT anomaly over Eurasia exhibits a meridional warm-coldwarm sandwich pattern (Fig. 6b). Similarly, the anomalous GHT at 500 hPa (200 hPa) shows a high-low-high pattern (Fig. 7a, b). This agrees with the theory that the high ridge (low trough) axis inclines to the warm (cold) regions in an asymmetric thermal system (Li and Zhu 2010;. The existence of the extremely high ridge and low trough over the Eurasia (Fig. 8a) were very likely or in part induced by the extreme regional meridional asymmetric warming (Fig. 6d). The enhanced South Asian High over the Iranian Plateau (Fig. 7a) may also associated with the extreme IM. A stronger IM is associated with more condensation heat release in the upper troposphere over the northern Indian peninsula. The anomalous heating stimulates positive GHT anomalies to its northwest in the upper troposphere over the Iranian Plateau, causing a stronger South Asian High (Wei et al. 2014).
The results based on the RAM emphasize the synergistic effect of multi-factors and highlight the remote effects of mid-latitude atmospheric variability on the simultaneous climate events in addition to tropical SSTA. We verified the robustness of the results of ERA5 (0.25° × 0.25°) using Hadley Centre Global Sea Ice and SST (1° × 1°) and NCEP/ NCAR reanalysis (2.5° × 2.5°), and found the above results based on the RAM are not sensitive to the choice of dataset and spatial resolution (not shown). This study suggests that the RAM is an effective and promising approach to search for possible factors that are responsible for this tremendous rainfall event.
Like any data-driven approach, the RAM is limited by the assumption that extreme climate events are linked with extreme factors. Some extreme events may also result from a combination of several non-extreme factors, while it is a great challenge to pay attention to these no-significant Fig. 9 Schematic for physical processes causing the record-breaking NC rainfall in September 2021. The abnormal atmosphere patterns in the Northern and Southern hemispheres are similar to the positive North Pacific Oscillation (NPO) and Southern Annular Mode (SAM), respectively. First, the strengthened mid-latitude subtropical high (light blue shading) in both hemispheres, warm SSTA over the tropical Atlantic (red shading) and cold SSTA over the tropical Pacific (blue shading) enhance the MC convection (red circle) and IM (orange circle) through zonal-vertical circulations (red arrows) and Rossby wave energy propagation (lines formed by triangles at 200 hPa and 850 hPa). Then the MC enhances the subtropical high (black circle at 850 and 500 hPa) over the southeastern China and induces strong ascending motion over NC. The joint effect of ascending motion, upper-level air acceleration (yellow shading) in the jet (blue vector) entrance region and anticyclonic circulation (black circle with vectors at 200 hPa) anchors the southwest-northeast rain belt (deep blue shading) over NC. The southern China subtropical high and IM induce considerable moisture transport to the rain belt via the two moisture channels (green arrows). Finally, the combined effect of extreme dynamic and moisture conditions genertates this unprecedented rainfall over NC in September 2021. The blue dashed lines indicate the distribution of ridge and trough, which contribute to the record-breaking rainfall east of the Baikal (deep blue shading). The black circle at 200 hPa is the upper-level anticyclonic circulation over the MC. The orange shading in this schematic represents the poleward lobes of the NPO-like and SAM-like patterns potential factors. For instance, the cold SSTA in eastern tropical Pacific contributed to the extreme event in NC in this study, while is not significant and can not be revealed by the RAM. In a warming climate, the frequency of worldwide extreme climate event occurrence is increasing (Rahmstorf and Coumou 2011;Cai et al. 2014a, b;Cai et al. 2014a, b;Chen and Sun 2015). The RAM makes an attempt on shedding light on a better understanding of these events. Besides extreme rainfall events, the RAM is also suitable for diverse types of extreme climate and weather events, such as record-breaking heat waves, droughts, and hurricanes. The potential of RAM for different types of events will be tested in future studies.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.