Geographical Trapping of Synchronous Extremes Amidst Increasing Variability of Indian Summer Monsoon Rainfall

Concurrent extreme rainfall events, or synchronous extremes, during Indian Summer Monsoon Rainfall (ISMR), cause significant damage, but their spatiotemporal evolution remains unclear. Using the event synchronization approach to examine the synchronicity of extreme rainfall events from 1901 to 2019, we find that Central India consistently hosts strongly connected synchronous extreme hubs with localized connections, indicating the geographical trapping of these concurrent events in the region. We observe a moderate positive correlation between network cohesiveness and El Niño Southern Oscillations (ENSO), and a negative correlation between ENSO and link lengths, suggesting localized synchronicity during El Niño dominant decades and opposite patterns in La Niña periods. Despite increasing ISMR variability and spatial nonuniformity, the persistence of hubs and network attributes could offer insights for predicting synchronous extremes, informing effective adaptation and risk management strategies during the monsoon season.

weather events and their potential impacts is crucial to managing risks and reducing the societal impacts of these events.
With the Indian Summer Monsoon (ISM) exhibiting increased variability in recent decades (Ghosh et al., 2016(Ghosh et al., , 2012;;D. Singh et al., 2014), the understanding of the intensity, duration, frequency, and spatial patterns of extreme precipitation during ISM remains a subject of ongoing debate (Ali et al., 2014;Ghosh et al., 2009Ghosh et al., , 2012;;Goswami et al., 2006;Mishra et al., 2012;A. Mondal & Mujumdar, 2015;Nikumbh et al., 2019;Rajeevan et al., 2008;Roxy et al., 2017;D. Singh et al., 2014).Traditional methods, such as spatial regression analysis and empirical orthogonal functions, do not investigate the spatial characteristics of extreme rainfall co-occurrences (Gupta et al., 2022) and other complex statistical interrelationships of higher order (S.Mondal et al., 2023).Thus, a robust nonlinear approach is necessary to evaluate extreme event time series and comprehend synchronization, particularly when examining the dependence between event time series across different locations.
Complex network (CN) theory, combined with non-linear synchronization measures, offers a new approach to tackle the problems described earlier.Researchers have used CN-based algorithms to analyze climate and weather extremes such as extreme rainfall (Boers et al., 2019;Gupta et al., 2022;Malik et al., 2012;Stolbova et al., 2014), heat waves (S.Mondal & Mishra, 2021), and droughts (Konapala & Mishra, 2017;S. Mondal et al., 2023;J. Singh et al., 2021).Furthermore, CNs have also helped map the influence of various teleconnections such as El Nio/South Oscillation (ENSO) (Tsonis & Swanson, 2008), North Atlantic Oscillation (Guez et al., 2013), and Rossby Waves (Boers et al., 2019) on regional climate systems.Over the Indian subcontinent, researchers (Malik et al., 2012) applied the CN method to investigate the spatial structure of rainfall extremes during the ISM and identified certain regions that receive rainfall only during the most active phase of ISM.In the study of large-scale extremes, researchers have used correlation networks and synchrony-based methods to understand the statistical properties of the underlying networks for extremes in rainfall, temperature, and droughts (Boers et al., 2019;Gupta et al., 2022;Konapala & Mishra, 2017;Malik et al., 2012;S. Mondal et al., 2023;J. Singh et al., 2021;Stolbova et al., 2014).However, little effort has been made to understand the evolution of spatial extent and temporal trends in the properties of the underlying networks.This understanding is crucial, particularly in light of claims that "stationarity is dead" (Milly et al., 2008) and that the assumption of stationarity in assessing hydrological extremes should be abandoned due to substantial anthropogenic changes in the climate of the Earth (Haasnoot et al., 2020;Lins & Cohn, 2011;Milly et al., 2008).
In the realm of understanding monsoon rainfall extremes, our study presents a comprehensive and long-term analysis of the evolution of network properties over the past century.While many existing studies provide a snapshot of these networks using composite data from all years, our approach delves into the dynamic nature of these networks, studying their evolution across 110 network instances over the last century.This approach provides new insights into the changing dynamics of Indian Summer Monsoon Rainfall (ISMR), offering a more nuanced understanding of this complex climatic phenomenon.Our findings are then corroborated with the known drivers of ISMR extremes, providing a coherent narrative for the observed trends and variability patterns in network attributes.
Specifically, our study aims to address several important questions related to synchronous extreme events, including (a) How has the synchronicity of rainfall extremes (both in terms of strength and scale) during ISMR evolved over the past century?(b) What is the spatial location of the most connected nodes or hubs?(c) How has the spatial pattern of these hubs evolved over the last century?(d) How do synchronous rainfall extremes statistically correlate with large-scale climatic patterns such as ENSO and sea surface temperature (SST) of the Tropical Indian Ocean (TIO)?Our results show that despite the increasing nonuniformity in the attributes of extreme rainfall events during the ISM over the last 119 years, we find a persistent and localized synchronization of extreme events in CI, which underscores the urgency for region-specific adaptation and risk management strategies.We show that CI regularly becomes a hotspot for simultaneously interconnected extreme weather events, indicating a geographical concentration of such synchronized extremes over the past century.Our research also uncovers significant connections between these synchronous extremes and global weather phenomena such as the El Nio Southern Oscillations and the TIO Sea-Surface Temperature.By providing insights into the correlation between synchronous rainfall extremes and these climate patterns, our study aims to shed light on the informed choice of covariates and, therefore, predictability of the synchronized extremes, whose management remains elusive (Jongman et al., 2014).

Rainfall and Climate Data Analysis
In this study, we analyzed daily 0.25° × 0.25° gridded rainfall data from the Indian Meteorological Department (IMD) (Pai et al., 2014) for the summer monsoon period (June to September or JJAS) from 1901 to 2019.We define extreme rainfall events as those exceeding the 90th or 95th percentile of all wet days during the summer monsoon months.We then calculate the corresponding annual frequency and intensity of daily rainfall exceeding these percentile thresholds.Figure S1 in Supporting Information S1 shows the 90th and 95th percentile threshold of daily rainfall during ISM for the period 1901-2019.We consider successive decadal periods over this period to understand the spatio-temporal evolution of synchronous events.
To estimate the association between SST and synchronous rainfall extremes, we use the Nino 3.4 index and SST data in the TIO region from 1901 to 2019.We spatially averaged the SST data over the region (40°E-120°E, 30°S-30°N) and took a decadal moving average for the summer monsoon months (JJA).To further confirm the generality of our findings, we also analyzed daily gridded rainfall data using the APHRODITE data set (Yatagai et al., 2012) available at 0.25° × 0.25° spatial resolution for the time period 1951-2015.
Finally, we study zonal and meridional wind components over the Indian subcontinent during the premonsoon (March-May) and monsoon (June-September) seasons using hourly resolution data at 500 and 850 hPa from the ERA5 reanalysis data set (Hersbach et al., 2023).We calculate the magnitude of the wind (V mag ) as the square root of the sum of the squares of the components of the zonal and meridional wind.

Statistical Tests
The Kolmogorov-Smirnov (KS) test is a widely used nonparametric method for comparing two samples to determine if they originate from a population with a particular distribution, while the Mann-Kendall (MK) trend test is frequently used in climatology and hydrology to determine the statistical significance of a trend in the time series of rainfall extremes.Details of these tests can be found in Huth and Pokorná (2004).The p-values, sample size, significance level, and slope for various parameters for which statistical tests are conducted, as summarized in Table S1 in Supporting Information S1.To measure the divergence of probability distributions, we employed Jeffreys divergence, which is a symmetric version of the Kullback-Leibler (KL) divergence (Kullback & Leibler, 1951).It is a type of statistical distance that measures the divergence (difference) between the two probability distributions (say P and Q).Specifically, the KL divergence from P to Q can be given as similarly, the KL divergence from Q to P is given as The Jeffreys divergence is defined as the sum of D KL (P‖Q) and D KL (Q‖P) given as Jeffreys divergence, like all divergences or distance measures, is a nonnegative quantity.A minimum value of 0 signifies identical distributions, with increasing divergence indicating an increasing difference between the distributions.Although lower divergences imply more similarity, there is no universal threshold for "low" divergence.Comparisons to a reference often provide a meaningful interpretation.To exemplify the Jeffreys divergence "low" and "high", we established a baseline Normal Distribution,   (0, 1) .We also tested four other distributions:   (10, 1) ,   (0.1, 1) ,   (0, 5) , and   (0, 1.1) (Figure S15 in Supporting Information S1).These yielded Jeffreys divergences of 1,663.33,0.17, 185.94, and 0.30, respectively, demonstrating the spectrum of divergence when distributions vary considerably.
To determine the synchronization between two grid sites i and j, we consider events l and m that occur at these grid points at times     and     , respectively, and are above a certain threshold α.We consider them to be synchronized if their time difference is within a certain range of dynamical lag  ±    , which is determined as follows: Here, l ranges from 1 to s i , and m ranges from 1 to s j , where s i and s j are the total numbers of events that occurred at the grid site i and j, respectively.This dynamical lag  ±    helps to distinguish independent events, accounting for different atmospheric processes that lead to precipitation at different time scales (Malik et al., 2012).We count the number of times an event appears at the grid site i after it appears at j and vice versa by defining c(i/j) and c(j/i), respectively, where and Similarly, we can define c(j/i).Then we calculate the synchronization strength Q ij between the points of the grid i and j as follows: where Q ij measures the strength of the synchronization of the extremes of rainfall between the points of the grid i and j and is normalized to 0 ≤ Q ij ≤ 1.Here Q ij = 1 means complete synchronization.After repeating this procedure for all pairs of grid points, we obtain a correlation matrix of size 4,331 × 4,331, where 4,331 is the number of grid locations in India.This correlation matrix is a square symmetric matrix that represents the strength of synchronization of extreme rainfall events between each pair of grid sites.

Construction of Synchronization Network
The connection between two nodes or points on the grid is formed when their synchronization falls within the top 5% of all synchronization values, signifying a high level of synchronization and an adequate number of rainfall events for comparison (Stolbova et al., 2014).Thus, the θ value is identified as the 95th percentile of all non-zero synchronization values.The binary adjacency matrix is constructed as follows: where A Q is a symmetric matrix with     = 1 representing a link between the grid sites ith and jth, and 0 represents otherwise.The adjacency matrix is used to quantify the network coefficients to characterize rainfall extremes.

Network Coefficients
Degree Centrality (DC) is the simplest CN measure, which is given as 10.1029/2023GL104788 5 of 14 where N is the total number of grid locations.C(Dj) measures the number of connections a node j has.Nodes with higher DC values have greater significance in network performance.During periods of heavy rainfall, areas with high DC values may synchronize with several other locations in the network, leading to a substantial impact on the spatiotemporal progression of synchronized rainfall extremes.
Previous climate network studies (Boers et al., 2014a(Boers et al., , 2014b(Boers et al., , 2015(Boers et al., , 2019;;Gupta et al., 2022;Malik et al., 2012;Stolbova et al., 2014), heat waves (S.Mondal & Mishra, 2021) have extensively used CN distance measures such as "average length of links" or "geographic distance."The Link Length (LL) L ij between two connected grid locations i and j is calculated by using the formula for spherical earth projected onto the plane as given where δϕ ij and δλ ij are the differences in latitude and longitude in radians between the grid locations i and j, ϕ m is the mean of the latitudes of i and j, and R is the radius of the earth.The average geographical LL     is the average geographical distance of ith node's links.
The average clustering coefficient (ACC) is calculated for the whole network and measures the overall network cohesiveness (Konapala & Mishra, 2017).Average clustering coefficient is the average of the local clustering coefficient (LCC) of a network's nodes, and can be expressed as where N is the total number of grid locations in a network.LCC i is the LCC of the ith node and is given as where L i represents the number of links between the k i network neighbors of node i.

Spatiotemporal Variability and Trends in Extreme Rainfall During ISMR
First, we analyze the statistical attributes of extreme precipitation (90th and 95th percentile) to gain insights into the spatiotemporal trends and variability in extreme rainfall over India during the ISMR, based on long-term data analysis.We use 0.25-degree gridded daily rainfall data obtained from the IMD from 1901 to 2019.Our analysis reveals a high variation in the magnitude and frequency of 90th percentile extremes, with higher values observed in the Western Ghats and Northeastern India (Figure 1a-1e).Furthermore, the frequency of such extremes is higher in east-CI adjoining the Bay of Bengal.However, only a few regions have a statistically significant trend in intensity (frequency) at the significance level of 0.01 (Ghosh et al., 2012), with 12.5% (13%) and 3.5% (11%) of the grid points showing increasing and decreasing trends, respectively (Figures 1c and 1f).
The spatial average of the intensity of extreme rainfall events in India shows an increasing trend.On the contrary, the standard deviation does not show a significant trend (Figure 1g).This suggests that, at aggregated levels, extreme rainfall events are becoming more intense on average with no significant increase in temporal variability.On the contrary, the variability of the frequency of extreme events exhibits an increasing trend when calculated based on all grids in India (Figure 1h).For more extreme rainfall considering 95th percentile, our findings are consistent with those of the 90th percentile, with the added knowledge that intensity variability also shows an increasing trend (Figures S2a-S2c for intensity and Figures S2d-S2f for frequency in Supporting Information S1).
While our results are in agreement with previous studies that highlight the lack of uniform trends in Indian rainfall extremes (Ali et al., 2014;Ghosh et al., 2009Ghosh et al., , 2012;;A. Mondal & Mujumdar, 2015), other studies have also reported that there is an increasing tendency for rainfall extremes to aggregate (Nikumbh et al., 2019(Nikumbh et al., , 2020;;Roxy et al., 2017), which may not be discernible through traditional statistical analyzes.Hence, it is crucial to investigate whether synchronous extreme events occur despite the lack of uniform trends across the region and if certain regions are particularly vulnerable to such events over long-term periods.Synchronous extreme event is defined by the statistically significant (typically defined by a threshold on the strength of synchrony) co-occurrence of extreme precipitation between two grid points.In contrast, widespread extreme events highlight the extensive spatial coverage of extreme precipitation without regard to the simultaneity of such occurrences.

Spatiotemporal Evolution of Network Attributes
We construct synchronization networks of extreme rainfall events during the ISMR at decadal scales to understand the space-time evolution of synchronous events.We use the ES method, a nonlinear correlation technique (Malik et al., 2012;Quiroga et al., 2002;Stolbova et al., 2014), to measure the synchronization strength between two locations.We generate 110 instances of these networks over the study period using a moving window with Intensity shows an increasing trend at the 1% significance level; however, its variability does not exhibit a statistically significant trend.(h) All India's average frequency of rainfall extremes and its standard deviation.The frequency of rainfall extremes shows no significant increasing trend at the 0.01 significance level; however, its variability exhibits an increasing trend.
a size of 10 years.We then examine topological properties, such as DC (number of network neighbors of a particular node), LL (geographical distance between a connected pair of nodes), and average LCC (a measure of network cohesiveness, which represents the extent of local coherence in extreme occurrences among the nodes with respect to a given node) (S.Mondal et al., 2023).
To ensure that each moving decade has a sufficient number of extreme events, we calculate the mean and standard deviation of extreme events at each grid point (Figure S3 in Supporting Information S1).Our analysis demonstrated that, while the number of extreme events varies spatially, most of the grid points exhibit a low standard deviation in the number of extreme events.This consistency in our data is coupled with the inherent robustness of the ES matrix methodology, which can handle varying numbers of extreme events across different windows.Furthermore, each 10-year period inherently encompasses a diverse range of climatic phenomena, including years with varying degrees of El Nino dominance.These variances, in addition to other factors that affect mean and extreme rainfall during the ISMR, ensure that each decade presents a sufficient number of events above the (percentile) thresholds under consideration.
The time average of the network attributes suggests that a considerable fraction of high DC nodes are concentrated over CI for both 90th (Figure 2a) and 95th (Figure S4a in Supporting Information S1).To understand whether high synchronization over CI is due to the interaction of nearby or far-placed nodes, we select 50 nodes/ grids of high-degree centrality and track their links with the rest of the network.We observe that synchronous rainfall extremes are highly localized (Figure 2b and Figure S4b in Supporting Information S1).The physically plausible explanation for this localization lies in the fact that CI receives a fair share of extreme rainfall during the active phase of the Indian monsoon (positive phase of intraseasonal monsoon oscillations (MISO)), with 60%-80% of these extremes being attributed to low-pressure systems (Krishnamurthy & Ajayamohan, 2010;Thomas et al., 2021) that cluster over this region during the active phase (Goswami et al., 2003).These extremes also arise from dynamic contributions from the Bay of Bengal and the Arabian Sea that generate large-scale extreme rainfall events (Nikumbh et al., 2020).Furthermore, the absence of a high DC node and therefore the effects of synchronization in South India and Northeast India (Figure 2a and Figure S4a in Supporting Information S1) can be attributed to the fact that these regions receive their share of rainfall during the negative phase of MISO, that is, the break period in the monsoon.While the rest of the Indian subcontinent receives maximum rainfall during the active phase, the northeast and southern regions receive maximum rainfall during the break phase (Rajeevan et al., 2010).The Western Ghats, although showing the highest frequency of rainfall extremes, show a low DC, possibly due to the barrier of high mountains along the Western Ghats, which inhibits the interaction of rainfall processes with the rest of India.Interestingly, the persistence of these nodes over the last century in a specific region with sustained localization hints toward stability in underlying generational processes of synchronized extremes.
Given recent debates around trends and variability in means and extremes of ISMR, and anthropogenic changes to the Earth's climate, understanding time variation in the context of network attributes of underlying synchronous networks of extremes is important.Time variance, in this case, denotes the changes in both the mean and variance of the network attributes that characterize these synchronous networks of extremes over the past century.
Our analysis of decadal variations in the DC for both (threshold) percentiles shows consistent synchronicity over CI with a slight shift in the high-degree nodes from east to west (Movie S1).Further, we generate the empirical probability degree distribution functions for degree distributions of the 110 variants of networks (Figure 2c).The spatial average of degree, calculated for all grids in India, does not show a statistically significant trend in mean and standard deviation despite oscillations around the mean (Figure 2d).In addition to examining the time invariance in mean and standard deviation, it is crucial to explore whether the underlying probability distribution functions of the DCs generated for different moving decades are also similar to each other.The physical significance of this investigation lies in its potential to reveal any changes in the overall behavior of the spatiotemporal structure of ISMR extremes.If the underlying probability distribution functions remain consistent over time, it suggests that the fundamental processes and patterns governing the synchrony of extremes are relatively stable.Conversely, if the distributions diverge or differ significantly, it could indicate that certain aspects of the system are undergoing substantial changes, potentially driven by external factors or intrinsic variability (Upadhyay et al., 2021).
To investigate the divergence between samples that do not come from the same parent distribution as identified statistically using the Kolmogorov-Smirnov test at 0.01 significance level (Figure 2e (lower triangle)), we employ Jeffreys divergence as a measure to evaluate the extent of the distance between these distributions (Figure 2e; upper triangle).Our findings indicate that even when the underlying distributions are not the same, the distances among the distributions remain relatively "low" (Figure 2e; see Section 2 for explanation of "low" divergence).Similar conclusions are drawn from the analysis of more extreme rainfall events (95th percentile) (Figures S4c-S4e in Supporting Information S1).To test the robustness of these results to the choice of data sets, we performed a similar analysis with different thresholds on the same data set (Figure S5 in Supporting Information S1), as well as on another data set (APHRODITE) in this case.Our analysis reveals a consistent pattern of extreme event localization in CI across both data sets (Figures S6 and S7 in Supporting Information S1) and for different thresholds.However, it is pertinent to highlight a unique pattern observed in the APHRODITE data set.We find a distinct decreasing trend in the average DC, alongside an increasing trend in its standard deviation for extreme rainfall events.This finding deviates from the patterns observed in the IMD data set.Furthermore, the DC value in a specific location is likely influenced by the DC values in nearby locations due to the spatial autocorrelation structure.Therefore, to address this, we perform an analysis using the residuals of DC.These residuals are obtained after regressing the DC values at each location against the DC values of surrounding grids up to 250 km in either direction.We observe similar patterns in PDFs, and Jeffreys divergence, thus supporting the residuals' time-invariance (Figure S8 in Supporting Information S1 for both percentiles).

Localized Connections of High-Degree Nodes
We analyze the scale and evolution of the geographical distances between nodes over a moving decadal window for all connections in the network.The time average of the average link lengths, calculated for all the networks in India during the entire period (Figure 3a, Figure S9 in Supporting Information S1), shows that the nodes with relatively high DC tend to have shorter link lengths (e.g., locations in CI) than the nodes with lower DC (e.g., areas in western and eastern ghats).Interestingly, correlations of degree-link length have been observed in disparate networks, such as social networks (Newman, 2001) and transportation networks (Sevtsuk & Ratti, 2010), which also generalize to synchronization networks.One reason why nodes with a low DC may have longer link lengths is that in networks with spatial or geographical constraints, nodes with a low DC may have fewer opportunities to connect to nearby nodes, resulting in longer link lengths.This effect is more pronounced when the geographical barriers are substantial (e.g., stretch along the Western Ghats), as it limits the opportunities for the nodes to form local connections.Furthermore, the synchrony of extremes in these regions may only occur during mesoscale events (Vijaykumar et al., 2021), influencing the larger area beyond geographic barriers.These events can create a synchronizing effect among the nodes, resulting in a higher degree of connectivity among the nodes that were previously disconnected.
Furthermore, we analyzed the stationarity of the topological attributes of synchronous rainfall extremes from 1901 to 2019 by generating empirical pdfs of the average LL over moving decadal windows (Figure 3 for the 90th percentile and Figure S9 for the 95th percentile in Supporting Information S1).For average link lengths, we again employ the KS test and Jeffreys divergence to determine whether two samples originate from different underlying distributions and quantify the extent of the divergence between them.Although statistical tests suggest that link lengths across moving decades are not always drawn from the same distribution (Figure 3; lower triangle), the low value of Jeffreys divergence indicates that the distributions across the study period exhibit minimal dissimilarities (Figure 3; upper triangle), yielding similar insights to those obtained from degree distributions.
Our results highlight the consistency of topological features in the network of synchronous rainfall extremes over the past 119 years.
Building on our analysis of changes and divergences in the distributions of network properties, we further explored the consistency of relationships between these network properties of synchronous rainfall extremes across India.Pairwise correlations between various network properties used throughout this study confirm this consistency, as evidenced by all 110 instances of networks generated (Figures S10a-S10c for 90th percentile extremes and Figures S10d-S10f for 95th percentile extremes in Supporting Information S1).In each case, the sign of the trend remains consistent across the network instances, underscoring the stable patterns and further substantiating the minimal temporal divergences observed in the network properties.
To gain further insight into the plausible reasons for the time invariance in topological properties of synchronous extreme networks, we analyzed changes in wind speeds at upper level (500 hPa) and near-surface (850 hPa) in our study area (Nikumbh et al., 2020).Our analysis reveals that during the period from 1940 to 2020, there are no widespread trends in wind speed magnitudes at either level during either the premonsoon (March to May or MAM) or the monsoon (June-September or JJAS) periods in the regions where high synchronization is observed when wind vectors were analyzed across India (y Figures S11 and S12 in Supporting Information S1).However, we observe increasing trends in wind speed magnitudes at specific grid points during the premonsoon season in the northern and northeastern regions, and during the monsoon season in the Tibetan region (Figure S11 in Supporting Information S1).In particular, in CI, where a persistent cluster of high DC nodes is located for the entire period, no significant trends were observed in the magnitude of wind speeds.This lack of trends could indicate dynamic stability in the circulation patterns in the region dominated by hubs (Christidis & Stott, 2015), contributing to the time invariance of the topological attributes of synchronization networks.

Network Attributes and Ocean-Atmospheric Processes
To better understand the impact of ocean-atmospheric processes like El Niño-Southern Oscillation (ENSO) and SST of the TIO on the intensity and frequency of ISMR (Ummenhofer et al., 2011;Yun & Timmermann, 2018), we contrasted various network attributes and the indices of these large-scale ocean-atmospheric phenomena (time-series in Figures 4a and 4b and corresponding scatterplots with 95% CI in Figures S13a-S13d in Supporting Information S1).Our findings show that the Nino 3.4 index correlates significantly with network cohesiveness, measured by the ACC (Figure 4a).Despite the SST anomaly of the TIO showing a strong correlation (at a significance level of 0.01) with the ACC (Figure 4b), visibly discernible differences in oscillations are observed.The sudden drop in the network ACC from 1940 to 1970 coincides with the cooling phase of Niño 3.4 and the SST anomaly of TIO (England et al., 2014).Furthermore, we observe that for the majority of the gridpoints, the average link lengths are comparatively shorter (Figure 4c) along with increased cohesiveness (Figure S13e in Supporting Information S1) in El Nino dominant decades than in La Nina dominant decades (Figure 4d; Figure S13f in Supporting Information S1), indicating the localization of extremes in the predominantly positive ENSO phase (Figures 4c and 4d).The most probable explanation is that during a positive phase of ENSO, the Indian monsoon circulation weakens, leading to reduced overall rainfall over the Indian subcontinent.However, this weakening could result in more localized extreme rainfall events, probably with reduced magnitudes compared to the La Niña years, particularly in regions more susceptible to convective rainfall (Yoshida et al., 2007), thus increasing the cohesiveness in networks of synchronous extremes.Conversely, during La Niña, the walker circulation over the tropical Pacific strengthens, causing the Indian monsoon trough to shift westward (Roxy et al., 2015).Furthermore, the western Indian Ocean tends to be warmer than normal (Roxy et al., 2015), drawing moisture away from the central and eastern parts of India and reducing rainfall in regions with high concentrations of high DC nodes.As a result, the rainfall patterns become more fragmented and less cohesive.The results for 95th%ile are shown in Figure S14 in Supporting Information S1.

Implications and Conclusions
Our study provides crucial information on the synchronization and spatial evolution of extreme rainfall events during the ISM over the past century.The persistent and localized synchronization of extreme events in CI, despite the increasing variability of the IMSR, underscores the need for targeted adaptation and risk management strategies in this region.Our findings suggest that CI consistently hosts strongly connected synchronous extreme hubs, indicating the geographical trapping of synchronous extremes.Despite changes in the magnitude and frequency of rainfall extremes, as demonstrated by the minimal dissimilarities in the distributions according to Jeffreys divergence, these networks can provide a robust framework for understanding the spatial patterns of extreme events.Furthermore, our study reveals significant associations of synchronous extremes with El Nino Southern Oscillations and TIO Sea-Surface Temperature, emphasizing the importance of considering coupled ocean-atmosphere processes when predicting concurrent extreme events.Despite fluctuations in the intensity and frequency of extreme rainfall, the relatively stable nature of network properties demonstrates that these networks can offer a reliable framework for understanding the spatiotemporal patterns of extreme events.
Our study has some limitations that need to be addressed in future research.This study does not consider the relative influence of atmospheric processes, tropical oceans, and other oscillations, such as the Indian Ocean Dipole and Rossby waves, on the synchronization structure and strength of rainfall extremes.The implications of these dynamic features on the structural stationarity of synchronous extremes need further exploration.Further studies are needed to obtain details of the dynamical processes leading to the geographic capture of extremes and their aggregation over the Indian subcontinent.
Understanding the complex interplay of extreme weather events and their driving processes is essential for skillful predictions and reliable projections to provide guidance to policymakers on mitigating risks and allocating resources.Our findings on structural stability present a promising framework for improving the predictions of extreme precipitation events during the ISMR.The complexities involved in refining coupled forecast models are well-acknowledged.In this context, our research offers a targeted methodology: we identify geographical regions prone to synchronous rainfall extremes.Focusing on these areas is an effective approach for sensitivity studies aiming to improve forecast models.By doing so, modeling approaches can channel their attention toward refining a subset of critical parameters within these regions, notably topography and land-atmosphere interactions.This focus could streamline these models' calibration and validation process, leading to more accurate and timely predictions of synchronized extreme rainfall events.The approach provides a structured pathway for enhancing the process and predictive understanding for of forecast models, specifically in regions that are most susceptible to extreme weather phenomena.
Moreover, the persistent 'hubs' of extreme rainfall events, as revealed in our study, also have practical applications in enhancing the predictability of such extreme events.Particularly, the time-invariant nature of these hubs can be utilized effectively by emerging machine learning models such as Graph Neural Networks, which are adept at learning patterns over heterogeneous data structures like graphs and networks.One of the prerequisites for applying these models in regression problems is a time-invariant graph structure while accounting for temporal dynamics at each node, which is provided by the persistent 'hubs' in our study.These models can learn the relationships between the attributes of the target variable and predictors (network neighbors that influence the results at locations of interest and exogenous variables, including indices of ocean-atmospheric coupled processes).Predictive models typically rely on geographical neighbors in space and time in the context of extremes.However, our research suggests the potential benefits of incorporating network neighbors, regions interconnected within our extreme event network, in addition to exogenous variables, into these predictive frameworks.Given the influential role of these network neighbors on the outcomes at locations of interest, integrating the network structure revealed in our study into these models can encapsulate the inherent spatial synchronizations among extreme events In summary, our study could help develop a framework to improve the predictability of synchronized extreme weather events.This could facilitate effective adaptation and management of risks related to synchronous rainfall extremes.This becomes increasingly crucial in the face of growing temporal variability and spatial non-uniformity associated with the Indian Summer Monsoon Rainfall.

Figure 1 .
Figure 1.Spatiotemporal distribution of frequency and intensity of rainfall extremes and their trends during the Indian Summer Monsoon.(a) Spatial distribution of the temporal mean of rainfall extremes (90th percentile) and (b) its standard deviation, with high values over the northwestern region and parts of NEI.(c) Spatial distribution trends of rainfall extremes.(d) Spatial distribution of frequency of rainfall extremes (90th percentile), with high values over NEI and Western Ghats.(e) Standard deviation of frequency, with high values over the Western Ghats and parts of NEI.(f) Spatial distribution of trends at each grid location, indicating decreasing trends in the frequency of rainfall extremes along the Indo-Gangetic plains and increasing trends along the east coast of India.(g) All India average of intensity of rainfall extremes and its standard deviation.Intensity shows an increasing trend at the 1% significance level; however, its variability does not exhibit a statistically significant trend.(h) All India's average frequency of rainfall extremes and its standard deviation.The frequency of rainfall extremes shows no significant increasing trend at the 0.01 significance level; however, its variability exhibits an increasing trend.

Figure 2 .
Figure 2. Degree Centrality (DC) and synchronization of rainfall extremes in India during 1901-2019.(a) All India average DC for 1901-2019 with nodes having higher degree centralities clustered over Central India (in Red).(b) Top 50-nodes with the highest DC (in Black) and the spatial extent of their connections.For these 50 nodes, the connections are primarily confined with the nodes in close proximity.(c) Empirical Probability Distributions of degree distributions (non-normalized total number of connections)) for 110 network instances with each decade shown in a different color.Three non-overlapping decades are shown as dotted PDFs (d) The spatially averaged degree and its standard deviation reveal no significant mean and standard deviation trend despite fluctuations around the mean.(e) Quantification of differences in generational processes and divergences among the distributions depicted in (c): The lower triangle displays the KS divergence at the 0.01 significance level, with blue regions indicating that samples are not coming from the same generational process.The upper triangle indicates the magnitude of Jeffreys divergence among the degree distributions.

Figure 3 .
Figure 3. Average Link Lengths and Probability Distribution Functions for 90th percentile extremes: (a) Spatial Distribution of Average LL: Notably high average link lengths are observed in regions with geographical barriers, such as the Western Ghats and sections of the Eastern Ghats (depicted in red).In contrast, Central India, or the core monsoon zone, exhibits relatively low values of average LL (represented in blue).(b) The lower triangle illustrates the Kolmogorov-Smirnov divergence at a 0.01 significance level, where gray regions signify that samples do not originate from the same generational process.The upper triangle displays the magnitude of Jeffreys divergence across all possible combinations of distributions.(c) Empirical Probability Distribution Functions for Average Link Lengths across 110 Networks: Each decade is represented by a unique color, with three non-overlapping decades portrayed as dotted PDFs.

Figure 4 .
Figure 4. Relationship between Network Properties and Ocean-Atmospheric Indices: (a) During decades with a positive phase of the Nino Index 3.4 (°C), a moderate yet statistically significant negative correlation (p < 0.01) is observed between average Link Length (LL) and the index, along with a moderate positive correlation with the average local clustering coefficient (LCC) (Average LCC).(b) Analyzing the correlations of average LL (red) and average LCC (blue) with Tropical Indian Ocean Sea-Surface temperature reveals patterns similar to the Nino Index 3.4(°C), with distinct differences emerging after 1980.(c) Comparing the average link lengths in decades dominated by positive Nino Index 3.4 (°C)  phases to the entire study period, an increase in average link lengths is apparent, with 57.5% of gridpoints exhibiting length increments.(d) A parallel analysis is conducted for decades dominated by a negative Nino Index 3.4 (°C) phase with 64% of nodes exhibiting reduced link lengths during these decades.