Hydrological impacts of altered monsoon rain spells in the Indian Ganga basin: a century-long perspective

The Indian Ganga basin (IGB) is one of the most valuable socioeconomic regions in the Indian subcontinent. The IGB supports more than half a billion people due to an abundant supply of freshwater for agro-industrial purposes, primarily through Indian Summer Monsoon (ISM) rainfall contributions (∼85%). Any alterations in ISM characteristics would significantly affect freshwater availability, and as a result, socioeconomic activities would be affected. Therefore, in this study, we have attempted to assess how the monsoon rain spell characteristics, i.e. peak, volume, and duration, altered historically between 1901 to 2019. We further analyzed the specific IGB regions where monsoon rain spell changes are more prominent and their hydrological implications. Our estimates reveal that short-duration high-magnitude rain spells have significantly increased across the major regions of the IGB after 1960, which implies the increased probabilities of flash flood hazards. At the same time, the rain spell volumes have been depleted across the IGB after 1960, especially in the eastern Indo-Gangetic plains and southern IGB regions, indicating increased drought frequencies. Further, Himalayan regions, i.e. upper Ganga, upper Yamuna, and upper Ghaghra, have demonstrated increasing magnitudes of rain spell peaks, volume, and duration post-1960. In addition, the continuous warming and anthropogenic alterations might further exaggerate the current situation. Thus, these inferences are helpful for river basin management strategies to deal with the extreme hydrological disasters in the IGB.


Introduction
Hydrological extremes, i.e. floods and droughts, refer to those events with a relatively lower probability of occurrence; however, they significantly impact socioeconomic activities and ecosystems across the globe (He et al 2020).These extreme events are predominantly influenced by climatic variability, changes in land surface responses (i.e.land use land cover alterations), and the development of hydrological infrastructures, i.e. dams & reservoirs (Dankers and Feyen 2009, Joshi et al 2023a, Dey et al 2024).Warming climates have increased the water-holding capacity of the atmosphere (Trenberth 2011), and thereby, increased occurrences of these extreme precipitation events have become frequent worldwide (Field andBarros 2014, Bloschl et al 2019).As a result, high-magnitude flooding events have become frequent, impacting infrastructures, causing biodiversity loss, and over 6.8 million losses of human life in the 20th century (WMO 2018).In the past few decades, nearly 50% of floods occurred in Asia's densely populated and highly flood-prone regions, resulting in the highest mortality compared to the other areas worldwide (Jonkman 2005).In addition, warming climates have amplified low-magnitude hydrological extremes, i.e. more frequent and severe drought conditions worldwide, which have significantly impacted agricultural productivity, freshwater scarcity, and loss of biodiversity (Mishra and Singh 2010, Bennett et al 2015, Van Loon 2015, Kim et al 2019).Compared to floods, droughts are slow-onset disasters and, therefore, can create complex socioeconomic problems for extended durations.The summer monsoon fed Asia, one of the highest drought-prone regions, and has experienced longer and more intense droughts since 1950 (Field et al 2012, Miyan 2015).
India is prone to both the hydrological extremes, i.e. floods and droughts, and these extreme events pose significant social, economic, and environmental impacts (Mujumdar et al 2020).Typically, droughts in Indian regions can be linked with extended duration abnormally low rainfall over large spatial extent (Sikka 1999).However, floods are characterized as extreme hydrological disasters that occur in a relatively smaller region for a shorter duration due to excessive rainfall and bank full conditions of rivers (Dhar and Nandargi 2003, Sharma et al 2018, Swarnkar and Mujumdar 2023).Estimates from previous studies suggest that nearly 40 million hectares of Indian regions are flood-prone, and around 8 million hectares of Indian land are affected by floods every year (Ray et al 2019).Estimates also indicate that approximately 18 million Indians are directly impacted by increased flood magnitude and frequency yearly (Rupa and Mujumdar 2019).In particular, Uttarakhand (2013), Jammu & Kashmir (2014), Tamil Nadu (2015), Gujarat (2017), andKerala (2018) were some of the catastrophic flood events that caused severe socioeconomic disbalance, loss of biodiversity, and human life (Mishra and Shah 2018, Ray et al 2019, Upadhyay et al 2020, Swarnkar et al 2021a, 2021b, Swarnkar and Mujumdar 2023).Similarly, drought events have significantly impacted socioeconomic and agricultural activities in the Indian regions.For example, around 1% of India's gross domestic product decreased in 2002 due to severe drought conditions (Gadgil et al 2003).Between 1901 and 2010, around 17% of the years were recognized as drought years, affecting agriculture, freshwater availability, biodiversity, and other socioeconomic activities in India (Kumar et al 2013).
The northern Indian river basins, i.e. the Ganga, Indus, and Brahmaputra Basins, are strongly influenced by the summer monsoon, which brings significant rainfall (Immerzeel et al 2009).These river basins support over 700 million people, and 60% of the population lives in the Ganga basin (Nepal and Shrestha 2015).The Indian summer monsoon (ISM) period (June to September) contributes approximately 80% of rainfall and is very crucial for direct and indirect freshwater demands for the Ganga basin (Gadgil andGadgil 2006, Auffhammer et al 2012).In particular, the Indo-Gangetic Plains (IGP) is critical for the agro-industrial activity that receives over 80% rainfall during ISM months (Bhatla et al 2015).Any significant deviations from the average monsoon rainfall can accentuate hydrological extremes across the country, including the densely packed Ganga basin.Therefore, understanding monsoon rainfall using long-duration observed rainfall records is vital.Most of the previous studies have assessed the trends and seasonal behavior of rainfall magnitude across the different parts of the country (Singh and Sontakke 2002, Ghosh et al 2009, Konwar et al 2012, Sinha et al 2015, Jin and Wang 2017, Preethi et al 2017, Swarnkar et al 2021b, Falga and Wang 2022).However, for a complete understanding of monsoon rainfall, the information related to monsoon rain spell characteristics, such as volume, peak, & duration, and their spatio-temporal evaluation is critical for assessing the impacts of hydrological extremes affected regions and planning flood mitigation strategies in the Ganga basin.Therefore, we tried assessing the spatio-temporal monsoon rain spell characteristics across the Ganga basin in this work to answer the following crucial research questions using century-long rainfall records.
(1) How have the monsoon rain spell characteristics, i.e. rain spells peak, volume, and duration, changed historically?(2) How have historical changes in monsoon rain spell characteristics amplified the hydrological extremes across the Ganga basin?

Study area
Indian Ganga basin (IGB) is selected as the study area.The IGB covers around 8 61 452 km 2 , which is geographically almost one-fourth, approximately 26.3% of the area of India, making it India's most extensive river basin (figure 1(a)).The IGB is home to a large human population, with over 500 million people in the region.The IGB is also a critical basin because it covers large agricultural, industrial, rural, and urban areas, with inhabitants relying on the freshwater from summer rainfall during monsoon and snow & glacial meltwater during pre-monsoon periods for irrigation, drinking water, and industrial purposes.Average annual rainfall has been recorded as about 1100 mm, with more than 85% of the total rainfall received during the monsoon period (June to September; figures 1(b) and (c)).Ganga has some major rivers as its tributaries, like Yamuna, Chambal, Betwa, Ghaghara, and Koshi, and minor rivers like Gomti, Tons, Sone, Damodar, and Gandak, which collectively form the entire Ganga basin.For this study, we have divided the whole Ganga basin into 13 sub-basins, i.e. numbered 1-7 in the northern region and 8-13 in the southern region of the basin, as shown in figure 1(a).The Level-4 HydroBASINS shapefiles from the HydroSHEDS dataset (Lehner and Grill 2013) were used to divide the IGB into sub-basins from northern and southern Indian regions.The IGB was divided into 13 sub-basins by merging several smaller sub-basins using the 'Union' tool in ArcGIS 10.8 (figure 1(a)).Chambal is the largest sub-basin with an area of 142 886 km 2 , whereas Sunderban delta sub-basin is the smallest with an area of 12 818 km 2 .In addition, the shuttle radar topography mission (SRTM; Jarvis et al 2008) digital elevation model (DEM) was used to depict the topographic variation across the IGB (figure 1(a)).The online sources of all the listed input datasets in this study are provided in the manuscript's 'data and code availability' section.

Selection of gridded rainfall dataset
We used an open-source daily gridded rainfall dataset of resolution 0.25  Rajeevan et al (2006Rajeevan et al ( ), (2008Rajeevan et al ( ), (2009)), global gridded rainfall dataset (APHRODITE) proposed by Yatagai et al (2012) and their gridded data with IMD operational dataset constructed using area weightage of district rainfall dataset.The monsoon rainfall time series derived from the 0.25 • (Pai et al 2014) dataset closely resembles the IMD operational dataset from 1901 to 2010.Around 0.26 mm daily rainfall bias for monsoon months (June-September) was observed between 0.25 • (Pai et al 2014) and the IMD operational rainfall dataset.Therefore, we argue that the 0.25 • gridded century-long daily rainfall dataset proposed by Pai et al (2014) is one of the most robust datasets available for rainfall analysis in India.Based on the aforementioned arguments, this study selected the Pai et al (2014) dataset for characterizing monsoon rain spells using 1157 0.25 • resolution IMD rainfall grids in the IGB.The IMD rainfall data was available in binary format.We used the C code provided by IMD Pune to convert each year's rainfall data into an ASCII file.Further, the generated ASCIIs were imported in MATLAB 2022 for structuring into days, months, and years.Finally, the structured data in '.mat' format was used for further rainfall analyses.It should be noted that this study does not include the mountainous regions that fall in Nepalese and Chinese territories because the observed daily IMD 0.25-degree rainfall dataset was unavailable for the selected 119 years in these regions.

Identification of monsoon rain spell characteristics
The ISM contributes a significant amount (more than 85%) of rainfall from June to September between 1901-2019 in the IGB (figure 1(c)).Therefore, the rainfall analysis was primarily focused on the ISM period (a total of 122 d) in the IGB.The monsoon rain spell rainfall characteristics, i.e. rain spell volume, rain spell duration, and peak rainfall, were estimated by identifying the rain spell events that occurred during the monsoon periods between 1901-2019 for 1157 grids present in the IGB.In a particular year, zero values (non-rainy days) and above zero values (rainy days) are possible at each grid in the daily monsoon rainfall dataset (figures 1(b) and (c)).Hence, the rainy days between two zero rainfall events were identified as a single monsoon rain spell event at a particular grid for a given monsoon year.The accumulated rainfall received during each rain spell was summed and considered as rain spell volume (in mm).The total duration of non-zero rainy days in an identified rain spell event was considered rain spell duration (in days).The daily maximum rainfall value from the determined rain spell events was considered peak (or extreme) rainfall (in mm).Multiple rain spells, and therefore, multiple values of rain spell volume, rain spell duration, and peak rainfall, were possible for a selected grid in a given monsoon year.Thus, the maximum values of rain spell volume, rain spell duration, and peak rainfall were estimated for a selected grid in the IGB in a particular monsoon year.These steps were applied to all the grids in the entire IGB to identify the maximum rain spell characteristics.The time series of maximum monsoon rain spell volume (hereafter V max ), maximum monsoon rain spell duration (hereafter D max ), and maximum monsoon peak rainfall (hereafter R max ) were identified for the entire period (119 years, 1901-2019) in the IGB (1157 grids).Finally, the 119 values of V max , R max , and D max for all the grids were used for further statistical analysis.Jenkinson (1955) introduced the concept of generalized extreme value (GEV) distribution.Since then, it has been widely used in flood frequency estimation, rainfall frequency analysis, wind analysis, and several other applications in natural and social sciences (Martins and Stedinger 2000, Katz et al 2002, Engeland et al 2004, Ailliot et al 2011, Rulfová et al 2016, Ballarin et al 2022).Thus, the GEV distribution was selected to analyze the probability densities of extreme value time series of V max , R max, and D max .The cumulative distribution function (CDF) of the GEV distribution can be defined as follows:

Extreme value analysis
where, 'x'is the extreme value time series (maxima), 'µ'is the location parameter, 'σ'is the scale parameter, and 'ξ 'is the shape parameter.The shape parameter (ξ ) defines the overall shape and tail behavior of the GEV distributions.For example, suppose ξ can be equal to zero, greater than zero, or less than zero.In that case, the GEV distribution can be classified as Type I (Gumbel distribution; ξ = 0), Type II (Fréchet distribution; ξ > 0), and Type III (Weibull distribution; ξ < 0) GEV distributions.The Fréchet (Type II) distribution is considered heavy-tailed compared to the other two distributions, i.e.Type I & II (Ailliot et al 2011).Independence and stationarity are essential criteria that must be satisfied before applying extreme value analysis.Therefore, the Wald-Wolfowitz (WW; Wald and Wolfowitz 1943) test was applied to assess the independence and stationarity of the extreme value time series.The null (H 0 ; p > 0.05) and alternate hypotheses (H a ; p < 0.05) of the WW test suggest fulfillment and nonfulfillment of independence and stationarity criteria, respectively.Further, the Anderson-Darling (AD) test was used to assess if any statistically significant dissimilarities were present when comparing two probability distributions (box plots).The null (H 0 ; p > 0.05) and alternate hypotheses (H a ; p < 0.05) in the AD test suggest significant and insignificant differences present in the two probability distributions.
The GEV distributions were fitted for the maximum rain spell characteristics, i.e.V max , R max, and D max, at each grid in the IGB.The parameters of CDF (equation ( 1)) were estimated using the 'maximum likelihood estimation (MLE)' method (Prescott and Walden 1980).Further, the performances of the fitted GEV distributions were analyzed by comparing the theoretical quantiles (obtained from fitted distribution) and empirical quantiles (obtained from empirical CDF) using non-parametric Spearman correlation (ρ).Further, the obtained shape parameters (ξ) from the fitted GEV distributions were analyzed to understand the heavy-tailed and light-tailed distribution of extreme rain spell characteristics at the individual grids in the IGB.
The extreme value magnitudes of V max , R max, and D max at 2 year (median), 10 year, and 100 year return periods were estimated using the fitted GEV distributions.In addition, 95% confidence bounds around the return period magnitudes were also estimated to analyze the accuracy of the fitted GEV models.It should be noted here that the extrapolation of frequency analysis is generally allowed up to twice the extreme value time series record.Therefore, to estimate the 100 year return period magnitudes for extreme variable time series, data points around 60 years were considered at each grid in the IGB.Further, the percentage differences in the post-and pre-1960 at 2-,10-and 100 year monsoon rain spell characteristics return levels were quantified with respect to the respective pre-1960 return levels.All the mentioned statistical analyses were performed using the 'extRemes' & 'trend' R packages and the inbuilt function 'cor' in R Studio (Version 4.2.0).The 'data and code availability' section lists the online sources of all the packages used in this study.
The percentage difference is calculated by subtracting the pre-1960 values from the post-1960 values and then dividing the result by the pre-1960 values.The percentage difference variations are shown by boxplots.Each boxplot displays the median (Q2) as the middle mark, with the 25th (Q1) and 75th (Q3) percentiles represented by the bottom and top margins of the box, respectively.A whisker is drawn from the upper quartile (Q3) to the biggest data point within 1.5 times the interquartile range (IQR; Q3-Q1) above Q3.Similarly, a distance of 1.5 times the interquartile range (IQR) is marked below the lower quartile (Q1), and a line is drawn down to the lowest recorded data point within this range.Due to the requirement that the whiskers terminate at a certain data point, the lengths of the whiskers may appear uneven despite the fact that both sides have the identical length of 1.5 times the IQR.Any data points that fall outside the whiskers are considered outliers and are shown separately.

Physical connections between pre-and post-1960 monsoon rainfall
Past studies (Kripalani and Kumar 2004, Kumar et al 2009, Rao et al 2012, Abish et al 2018, Todmal et al 2022) reported that the sea surface temperature of western Indian Ocean is non linearly increasing, particularly after 1960.It is argued that the warming in the western Indian Ocean is perhaps influenced by surplus heat transported from the Pacific Ocean to the Indian Ocean (Dong et al 2016), increased solar flux absorption due to less cloud cover (Mayer et al 2013) and increased greenhouse gases concentration (Du and Xie 2008).Kumar et al (2013) explained that increased warming of the equatorial Indian Ocean may affect the (a) enhanced convection and (b) land-sea temperature contrast.Meehl (1994), Li and Yanai (1996) and Chou (2003) further explained that the intensity of the Asian monsoon is directly related to the strength of the land-sea temperature gradient.Increased evaporation (latent heat loss) over the Indian Ocean possibly increases the land-sea temperature contrast (Sutton et al 2007).Roxy et al (2015) indicated that increased warming of the Indian Ocean drives the weakening of monsoon circulation and rainfall over north and central Indian regions.Mishra et al (2016) estimated that the entire Indian region has witnessed a significant reduction in monsoon rainfall in the post-1960 period compared to pre-1960 rainfall.Further, Mishra et al (2016) also pointed out that the Indo-Gangetic plains have experienced a very significant decline in the monsoon rainfall volume post-1960.Our 30 years normalized rain spell patterns also reveal that there is significant reduction in rain spell volumes and duration post-1960 in the IGB (figures 2(a) and (c)).
Furthermore, several studies (Krishnan et al 2016, Guhathakurta et al 2017, Joshi et al 2023b) have pointed out that the local surface warming of the Indian subcontinent and the rise in humidity levels have a direct role in increased extreme rainfall events.Roxy et al (2017) revealed a threefold increase in widespread extreme rainfall across the Central Indian regions primarily because of the significant contribution of moisture transport from the Arabian Sea since the 1950s.Our monsoon rain spell characteristics results reveal that there is significant alteration in the maximum rainfall variability post-1960 in the IGB (figure 2(b)).Hence, based on the above-mentioned studies, it is evident that significant changes occurred in the monsoon rainfall pattern post-1960 due to increased sea surface warming of the Indian Ocean, high amount of moisture carrying capacity, weakening of monsoon and land-sea temperature gradient contrast.Thus, the pre-and post-1960 division of rainfall dataset can reveal significant changes in multi-decadal monsoon rain spell characteristics in the IGB, which includes vast Indo-Gangetic regions in the North and peninsular regions of Central India in the South.We divided the monsoon rain spell characteristics time series of the IGB into two periods, i.e. pre-1960 (1901-1960) and post-1960 (1961-2019).In addition, the 60 year extreme value time series is a good record length for predicting 100 year return levels (lower frequency extreme events) in the IGB monsoon rain spell characteristics.

Historical data analysis and application of GEV distributions
The annual V max , R max , and D max time series show decreasing magnitude for V max & D max and higher variability for R max post-1960 in the IGB (figures 2(a)-(c)).To visualize these changes further, we have normalized the 30 year moving mean of monsoon rain spell characteristics (V max , R max , and D max ) time series with respect to their mean and standard deviations to determine whether any statistically significant differences exist between pre-and post-1960 rain spell characteristics (figures 2(d)-(f)).The 30 year moving window is considered an ideally large sample size for regions that are greatly influenced by monsoon rainfall (Li et al 2018).It is evident from these figures that the monsoon rain spell characteristics have changed after 1960 (figures 2(d)-(f)).For example, the spatial mean of 30 year moving average V max and D max time series anomalies show a consistently decreasing trend after 1960 (figures 2(d) and (f)).On the other hand, the spatial mean of 30 year moving average R max shows an increased variability after 1960 (figure 2(e)).As a result of a decrease in the amplitude of R max peak values (figure 2(b)) and a decrease in the number of samples (with an average value of less than 30 years moving window) for recent years, the 30 year moving average R max shows a falling trend after 1980.The AD test results, when applied to 30 years normalized V max , D max , and R max time series, revealed that there are statistically significant (p < 0.05) changes present while comparing boxplots of pre-and post-1960 30 year normalized rain spell characteristics (figures 2(g)-(i)).The normalized R max shows slightly lower statistical significance (p = 0.12; figure 2(h)).However, the post-1960 normalized 30 year moving mean R max boxplot shows an increase in variability compared to the pre-1960 (figure 2(h)).
The pre-and post-1960 annual monsoon rain spell characteristics were further analyzed for independence and stationarity using the WW test.The WW test results revealed that more than 95% of grids present in the IGB fulfill the independence and stationarity criteria (p > 0.05; null hypothesis) for monsoon rain spell characteristics during pre-and post-1960 periods (figures 2(j) and (k)).Thus, the extreme value analysis can be performed on the pre-and post-1960 extreme value time series for the IGB.The GEV was fitted to pre-and post-1960 V max , D max , and R max time series.Their performances were assessed by comparing theoretical and empirical quantiles using the Spearman correlation.The correlation results suggest that GEV distributions performed well for pre-and post-1960 monsoon rain spell characteristics (figure 3(a)), with both pre-and post-1960, V max , and R max correlations being close to one.The correlation values ranged from 0.98 to 1 for the pre-and post-1960 D max .The shape parameter (ξ) values from fitted pre-and post-1960 GEV pdf reveal that more than 70% of grids demonstrate heavy-tailed (ξ > 0; Fréchet distribution) monsoon rain spell characteristics in the IGB, indicating a higher probability of occurrences of extreme values in the time-series (figure 3(c)).The fitted GEV pdf was further analyzed to assess the 2 year (median), 10 year, and 100 year return levels for the pre-and post-1960 monsoon rain spell characteristics.Figure 3(b) shows 2-, 10-and 100 year return levels for pre-and post-1960 monsoon rain spell characteristics with 95% confidence bounds for an example grid (located at 31.25 • North 78 • East) in the IGB.The following sub-sections further describe the percentage differences obtained by comparing pre-and post-1960 monsoon rain spell characteristics return levels and their possible hydrological implications.

How has the monsoon rain spell volume (V max ) changed after 1960?
Maximum monsoon rain spell volume (V max ) analysis was performed for pre-and post-1960 periods (figure 4).Their magnitudes at 2-, 10-and 100 year return periods and percentage differences are shown in figure 4 & summarized in table 1.The upper mountainous regions of the northern IGB sub-basins, i.e.Yamuna, Ganga, Bhagirathi, and lower plateau regions in the southern IGB sub-basins, i.e.Chambal, Betwa, Tons, and Sone, have relatively higher magnitudes (>500 mm) of V max at 2-(figures 4(a) and (b), 10-(figures 4(d) and (e)) and 100 year return periods (figures 4(g) and (h)) during the entire study period.The Indo-Gangetic plains (IGP; low elevation regions; see figure 1(a)) bounded by northern Himalayan (sub-basins 1-6) and southern plateau regions (sub-basins 8-12) have relatively lower magnitudes (<500 mm) of V max at 2-(figures 4(a) and (b)), 10-(figures 4(d) and (e)), and 100 year return periods (figures 4(g) and (h)) during pre-and post-1960.It is further noted that the V max magnitudes moderately increase from western to eastern Indo-Gangetic plains in the IGB at all three frequencies (figure 4).
Figures 5(a)-(d) show the difference between post-and pre-1960 V max and the area covered by positive differences for the northern and southern IGB sub-basins.Our post-and pre-1960 maximum monsoon rain spell volume difference (∆V max ) results reveal that the sub-basins (1-7) originating from the northern IGB show a predominantly decreasing percentage of ∆V max (up to −70%; figures 4(c), (f), (i), 5(a) and table 1).However, there are some positive clusters present in Yamuna (up to 50% sub-basin area), Ganga (up to 50% sub-basin area), and Ghaghara (up to 60% sub-basin area; figure 5(c)) that exhibit increasing percentages of ∆V max (up to 175%) at all three frequencies.Exceptionally, the positive difference cluster located at the upper Ganga basin (figure 4(c), (f) and (i)) shows an increase in ∆V max up to 85%, 90%, and 150% at 2-, 10-and 100 year, respectively.On the contrary, a decreasing ∆V max cluster (up to −70%) shows near the lower mountainous regions of the Ganga and Yamuna sub-basins at 2 year, 10 year, and 100 year frequencies (figures 4(c) and (f)).Further, the northeastern sub-basins in the IGB, i.e.Gomti, Gandak, Koshi, and Bhagirathi, exhibit more than 80% of regions experienced decreasing ∆V max at 2-and 10 year return periods (figures 4(c), (f), 5(a), 5(c) and table 1).
Further, the sub-basins (8-13) originating from the southern IGB regions also show a predominantly decreasing percentage of ∆V max (up to −70%; figure 4(c), (f), (i), 5(b) and table 1).However, some small positive ∆V max clusters are present in the western region of the sub-basins, i.e.Chambal (up to 50% sub-basin area) and Betwa (up to 30% sub-basin area).Further, some positive ∆V max clusters are also present in the eastern region of the sub-basins, i.e.Damodar (up to 58% sub-basin area) and Sunderban (up to 95% sub-basin area; figure 5(d)) in the IGB.The areal extents of the positive ∆V max clusters in the southwestern IGB sub-basins are relatively widespread compared to the southeastern IGB sub-basins.In particular, the Chambal shows a relatively increasing ∆V max in the small patches located at the upstream and downstream regions for 2-(up to 25%), 10-(up to 40%), and 100 year (up to 130%) return periods (see figure 5(b)).In addition, almost all the regions in Tons and Sone have significantly experienced a reduction (up to −50%) at 2-, 10-and 100 year frequency V max (figures 4(c), (f), (i), 5(b) and table 1).Overall, this analysis suggests decreasing ∆V max and negative ∆V max areal extents are relatively more prominent in the southern and the northern IGB sub-basins.

How has the extreme monsoon rainfall (R max ) changed after 1960?
Peak Monsoon rain spell rainfall (R max ) analysis was performed for pre-and post-1960 periods (figure 6).Their magnitudes at 2-, 10-and 100 year return periods and percentage differences (∆R max ) are shown in figure 6 & summarized in table 1.The northwestern to northeastern foothill regions present in the northern sub-basins (1-7) demonstrate a relatively highest magnitude of R max at 2 year (up to 146.8 mm; figures 6(a) and (b)) and 10 year (up to 249.6 mm; figures 6(d) and (e)) frequencies during the entire study period.However, in the southern IGB sub-basins, only Chambal and Betwa show significantly high R max magnitudes at 2-, 10-and 100 year frequencies during pre-and post-1960 (figures 6(a), (b), (d) and (e)).It is also observed that the Ghaghra, Chambal, and Betwa have consistently higher dominance of higher magnitude R max during pre-and post-1960 (figures 6(g) and (h)).On the other hand, the western Indo-Gangetic plains partly covered by Yamuna & Chambal and the upper Ganga basin have a consistently lower magnitude of R max at all three frequencies during the entire study period (figure 6).
Figures 7(a)-(d) show the difference between post-and pre-1960 R max and the area covered by positive differences for the northern and southern IGB sub-basins.Our post-and pre-1960 maximum monsoon rain spell rainfall differences (∆R max ) results reveal that the sub-basins (1-7) originating from the northern IGB regions are majorly covered by increasing magnitudes R max at 2-, 10-and 100 year return periods (figures 6(c), (f), (i), 7(a), (c) and table 1).It is also interesting to note that these positive ∆R max clusters expand and cover almost all the region of the northern IGB for higher magnitude extreme (100-year) monsoon rain spell events (figure 6(c), (f), (i), 7(c) and table 1).In particular, the upper Ganga, whole Ghaghra, upper Gandak, and major regions of the Bhagirathi show more than 70% area covered with positive ∆R max at 2-, 10-and 100 year return periods with increasing ∆R max up to 62.4%, 73.9%, and 182.2%, respectively (figures 6(c), (f), (i), 7(a), (c) and table 1).
Further, the sub-basins (8-13) originating from the southern regions of the IGB are also predominantly covered by increasing ∆R max at 2-, 10-and 100 year return periods (figure 6(c), (f), (i), 7(a), 7(c) and table 1).Interestingly, the total area covered by positive ∆R max is relatively higher in the southern than the northern IGB sub-basins at 2-, 10-and 100 year frequencies (figures 6(c), (f), (i), 7(c) and (d)).Specifically, the Chambal and Betwa regions show higher aerial dominance of positive ∆R max at 2-, 10-and 100 year return periods (figure 6(c), (f), (i), 7(c) and (d)).Further, the most southeastern regions, i.e.Damodar and Sunderban, demonstrate positive ∆R max , and more than 90% of the sub-basin areas are covered with positive ∆R max changes post-1960 (figures 6(c), (f), (i) and 7(b)).Comparatively, the Sone is the only region in the southern IGB that shows a relatively lower increase in ∆R max and lower dominance of positive ∆R max differences areas (figures 6(c), (f), (i), 7(c), (d) and table 1).Overall, this analysis indicates significant IGB regions experienced increasing peak rainfall magnitudes post-1960.

How has the monsoon rain spell duration (D max ) changed after 1960?
Maximum monsoon rain spell duration (D max ) analysis was performed for pre-and post-1960 periods (figure 8).Their magnitudes at 2-, 10-and 100 year return periods and percentage differences (∆D max ) are shown in figure 8 & summarized in table 1.The northwestern mountainous regions located in Yamuna, Ganga, and Ghaghra reveal relatively higher magnitudes of D max at 2-(>25 d; figures 8(a) and (b)), 10-(>40 d; figures 8(d) and (e)), and 100 year (>50 d; figures 8(g) and (h)) frequencies compared to other regions in the northeastern IGB.It is further noted that the western alluvial plains in the northwestern IGB sub-basins receive relatively lower duration rain spells than the alluvial plains of the northeastern IGB sub-basins at all three frequencies (see figure 8).Comparatively, the southern IGB sub-basins show the dominance of more extended duration monsoon rain spell regions compared to the northern IGB sub-basins at all frequencies during pre-and post-1960 (figure 8).In particular, the upper Sone region exhibits a relatively longer duration of monsoonal rain spells (figure 8).Further, the western Chambal is the only region that shows relatively higher dominance of short-duration monsoon rain spells in the southern IGB sub-basins (figure 8).
Figures 9(a)-( d) show the difference between post-and pre-1960 D max and the area covered by positive differences for the northern and southern IGB sub-basins.Our post-and pre-1960 maximum monsoon rain spell durations difference (∆D max ) results reveal that the sub-basins (1-7) originating from the northern regions of the IGB predominately cover (more than 60% sub-basin area; figure 9(c)) by decreasing ∆D max (up to −70%) at all three frequencies (figure 8(c), (f), (i), 9(a) and table 1).Wherein sub-basins 4-7 show significantly depleted ∆D max at 2-, 10-and 100 year return periods (figures 9(c), (a) and table 1).In contrast, some positive ∆D max clusters (figure 8(c), (f) and (i)) are present in the northern sub-basins, i.e. upper Ganga, middle Ganga & Yamuna, and upper Ghaghra.The maximum difference between post-and pre-1960 D max is 254% at a 100 year return period in Yamuna (figures 8(i) and 9(a)).
Further, the sub-basins (8-13) originating from the southern regions of the IGB show major areas are covered with negative ∆D max at 2-, 10-and 100 year return periods (figures 8(c), (f), (i), 9(b), (d) and table 1).Comparatively, the southern IGB sub-basins depict higher dominance of decreasing ∆D max clusters at all three frequencies than the northern IGB sub-basins (figures 8(c), (f), (i), 9(b), (d) and table 1).The southwestern sub-basins, i.e.Chambal, Betwa, Tons, and Sone, have faced significant decrement (up to −67%) in post-and pre-1960 D max at all three frequencies.In contrast, the Damodar and Sunderban sub-basins show relatively increasing ∆D max (up to 100%) covered by slightly higher area (up to 68%) compared to other southern regions of the IGB (figures 8(c), (f), (i), 9(b), (d) and table 1).Overall, these results indicate that significant regions (especially, the southern IGB sub-basins) experienced short-duration rainfall events post-1960.

What are the hydrological implications?
Alterations in the monsoon rain spell characteristics, i.e. peak, duration, and volume, can affect the hydrological processes significantly.For example, high-magnitude short-duration rainfall events can quickly saturate the upper soil surface (Hortonian overland runoff; Subramanya 2017), leading to flash flood hazards (Gaume et al 2009, Archer andFowler 2018).Steep topographic regions, i.e. the Himalayan region in the   northern IGB and steep plateau regions in the southern IGB (figure 1(a)), can further exaggerate the flash flood hazards due to the shortening of overland flow travel time (Kelsch et al 2001, Rao et al 2014, Mahmood and Ullah 2016, Swarnkar et al 2020).In addition, these flash flood hazards not only affect the runoff generation mechanism but also intensify the soil erosion and sediment transport rates that lead to offsite and onsite sedimentation problems (Narayana and Babu 1983, Wasson 2003, Swarnkar et al 2018, 2021c, Kumar et al 2022, Sinha et al 2023).Figures 10(a)-(f) show that R max has increased significantly as compared to V max and D max after 1960 at all three frequencies across the IGB.The relative increments in R max are more prominent for the southern than the northern IGB sub-basins (figures 10(a)-(f)).This comparative analysis indicates that higher magnitude rain spell events increased across the major regions of the IGB post-1960.In particular, the upper part of northcentral and northeastern IGB sub-basins (i.e.Ghaghra, Gomti, Gandak, Koshi, Bhagirathi) and the southern sub-basins (Chambal, Betwa, Tons and Sunderban) have shown the increased occurrence of high-intensity rainfall at the higher (at 2-and 10 year return periods) and lower (100-year return period) frequencies after 1960 (figures 10(a)-(f)).Remarkably, differences between post-and pre-1960 R max are highest for lower frequency events (100-year return period) post-1960 (figure 6(i)), further signifying extreme intensity rainfall events have increased across the Ganga basin (figure 10(c) and (f)).
Our historical analysis further indicates that there have been significant reductions in monsoon rain spell volume across the IGB post-1960 (figures 4, 5, 10 and table 1).As already discussed, most of the regions in the IGB receive more than 85% rainfall during monsoon months, and therefore, monsoon rain spells are critical for the IGB regions (such as the entire Indo-Gangetic Plains and southern IGB plateau regions; figures 1(a)-(c)) that primarily depend on rainfall from southwest monsoon (Sinha andFriend 1994, Bhatla et al 2019).Therefore, any drop in monsoon rain spell volume will further reduce freshwater availability and impose severe stress on agricultural and industrial activities.Figures 4(c), (f), (i), 10(d)-(f), and table 1 show major sub-basins in southern IGB (i.e.Chambal, Betwa, Tons, Sone) have experienced significantly depleted V max post-1960 as compared to the northern IGB.The Sone and Tons sub-basins have experienced the highest reduction in post-1960 V max .Hence, these southern IGB sub-basins will need attention to deal with frequent drought conditions.The alluvial plains in the IGB, also known as Indo-Gangetic plains, support more than 400 million people through the direct and indirect supply of freshwater for agricultural, household, and industrial purposes (Beg et al 2022, Jha et al 2022).Except for the western Indo-Gangetic plains present in the Yamuna and Ganga, the rest of the alluvial plains present in the IGB have shown moderate (up to −10%) to higher (up to −70%) V max reduction after 1960 at all three frequencies (figures 4(c), (f) and (i)).Especially the 2 year ∆V max further indicates that most of the regions in the IGB received relatively depleted amounts of monsoon rain spell volume at a higher frequency than the pre-1960 rain spell volumes.Thus, these inferences reveal that historically reduced monsoon rain spell rainfalls at higher and lower frequencies have significantly impacted the accumulated rain spell volumes received during the monsoon period after 1960.As a result, the freshwater availability was directly affected by the reduced monsoon rain spell volume, especially in the southern IGB sub-basins and Indo-Gangetic plains that have experienced frequent droughts for the past 60 years (Thomas et al 2015, 2016, Nath et al 2017, Lal 2022).
Further, the upper Yamuna, upper Ganga, and upper Ghaghra are exceptional regions where we observed significant increments in V max , R max, and D max at higher (2-year), moderate (10 year), and lower (100 year) frequencies in the post-1960 period.These results indicate that extreme-intensity long-duration monsoon rain spells have increased, accumulating relatively higher amounts of monsoon rain spell volume for extended periods post-1960.Thus, the likelihood of saturation excess runoff generation processes (Bronstert et al 2002) has been amplified in the upper Yamuna, upper Ganga, and upper Ghaghra after 1960.Consequently, we have witnessed several high magnitudes frequent floodings in the upper Ganga (Chawla et al 2018, Swarnkar et al 2021a, 2021b, Swarnkar and Mujumdar 2023), upper Yamuna (Kumar et al 2019, Tomar et al 2021), and upper Ghaghra (Houze et al 2017, Sudeepkumar et al 2023).Recently, high intensity rainfall in the upper Yamuna area caused the 2023 Yamuna flood.It is further noted that these regions (especially upper Ganga and upper Ghaghra; figure 1(a)) are surrounded by relatively steep topography.Hence, the quick basin response due to the steep terrain of these upper sub-basins might have further exaggerated the flood hazards to the downstream populated regions.For example, Swarnkar and Mujumdar (2023) showed how topographic steepness has influenced the basin response of the Alaknanda and the Bhagirathi sub-basin of the upper Ganga basin.Many populated cities, including the capital of India (New Delhi; population over 11 million), are situated close to these flood-prone upper mountainous regions.Overall, the results obtained from this study reveal that the probability of the occurrences of the hydrological extremes has been increased in the IGB.On the one hand, we observed a significant increment in short-duration high-magnitude monsoon rain spells, suggesting increased occurrences of flash flooding events after 1960.However, on the other hand, we also detected depletion in monsoon rain spell volume across the IGB sub-basins, implying significant pressure due to frequent droughts on freshwater availability after 1960.

Conclusion
More than half a billion people reside in the IGB, making it one of the most densely populated regions in the world.The Ganga River is a critical source of freshwater for agriculture, industrial, household demands, and other economic activities, provided by the ISM for four months (June to September).Any alterations in ISM rainfall characteristics in the IGB significantly impact freshwater availability and water-related disasters.Therefore, understanding the ISM rainfall characteristics using century-long ground-based observations is critical and, at the same time, valuable for long-term hydrological planning.In the present work, we have analyzed the monsoon rain spell characteristics, i.e. peak, volume, and duration in the IGB, and their possible hydrological implications.The extreme value distribution, i.e.GEV, is used to analyze the maximum monsoon rain spell characteristics over the past 119 years  at the 1157 grids in the IGB.The statistical analysis focuses on the maximum values of rainfall attributes (such as volume, peak, and duration) during particular monsoon seasons from 1901 to 2019 for each grid in the IGB region.The following key findings are concluded from this study: (1) There has been notable variation in the maximum rainfall (R max ) in most of the IGB zones since 1990.
Chambal, Betwa, Tons, Sunderban, northern Himalayan and foothill sub-basins have seen an increase in R max (maximum rainfall) since 1960.Consequently, our investigations indicate higher chances of flash flood-related hazards across these particular sub-basins of the IGB.(2) In contrast, the monsoon rain spell volumes have decreased in major regions of the IGB.Specifically, the southern IGB sub-basins and Indo-Gangetic plains have shown a relatively higher depletion of monsoon rain spell volume than the northern IGB sub-basins.We identified Betwa, Tons, and Sone are the worst affected sub-basins in the southern IGB regions in terms of magnitudes and total area coverage of depleted monsoon rain spell volume regions.Thus, these results reveal the increased possibility of drought events in the southern IGB sub-basins.
(3) The upper Yamuna, upper Ganga, and upper Ghaghra sub-basins have experienced increased monsoon rain spell volume, peak rainfall, and duration.The increase in large-duration & high-magnitude and steep topography exaggerated the occurrences of high-magnitude flood events, including the Kedarnath flood disaster in 2013.Therefore, these flood-prone mountainous regions should be dealt with more precise and scientifically sound flood management approaches.
These conclusions strongly indicate increased occurrence of extreme hydrological hazards in most of the IGB regions that affected the overall freshwater availability and water-related disasters.This research provides an in-depth understanding of the evolution of extreme monsoon rain spell characteristics in the Ganga basin over the past century.Additionally, this study also encourages the use of recent scenario based projected future rainfall data to improve our comprehension of the development of future monsoon rain patterns in the Ganga basin.Some additional sub-basins, including Sone, Tons, and Ghaghra, have not been well examined in terms of their historical and projected rainfall patterns and how they affect the basin's hydrology.The findings of our study strongly indicate that there has been a significant increase in the intensity of hydrological events in these less-studied river basins.Hence, conducting thorough investigations on extreme hydrological phenomena such as floods and droughts is of utmost importance for the under-researched areas within the Ganga basin.In a nutshell, this research indicates that it is crucial to design and effectively execute sustainable basin management strategies in order to alleviate the impacts of extreme hydrological events in the IGB.

Figure 1 .
Figure 1.(a) Study area map shows the topography using shuttle radar topography mission (SRTM), major rivers & tributaries, and 13 sub-basin divisions of IGB using HydroBASIN dataset, (b) total rainfall during the monsoon period over IGB on each grid point (in mm) (c) percentage contribution of monsoonal rainfall in total rainfall estimated from the daily rainfall data provide by Indian meteorology department (IMD) at 0.25 degree.

Figure 2 .
Figure 2. (a)-(c) show normalized 30-year moving mean of monsoon rain spell characteristics time series with respect to their mean and standard deviations to visualize the significant differences in the trends in pre-1960 and post-1960.(d)-(f) Anderson-Darling (AD) test results of 30-year normalized Vmax, Rmax, and Dmax dataset for the confirmation of the statistical significance (p < 0.05) difference between pre-1960 and post-1960 pdfs, (g)

Figure 3 .
Figure 3. (a) Spearman correlation analysis between the fitted GEV (theoretical) and empirical quantiles of Vmax, Rmax, and Dmax for pre-and post-1960 dataset, (b) representative plot for frequency analysis an example grid

Figure 5 .
Figure 5. Box plots for the (a) northern and (b) southern IGB sub-basins post-and pre-1960 Vmax percentage difference and (c) and (d) positive difference percentage area analysis for ∆Vmax at 2-, 10-and 100-year return period frequencies over IGB.The percentage difference is calculated by subtracting the pre-1960 values from the post-1960 values and then dividing the result by the pre-1960 values.

Figure 7 .
Figure 7. Box plots for the (a) northern and (b) southern IGB sub-basins post-and pre-1960 Rmax percentage difference and (c) and (d) positive difference percentage area analysis for ∆Rmax at 2-, 10-and 100-year return period frequencies over IGB.The percentage difference is calculated by subtracting the pre-1960 values from the post-1960 values and then dividing the result by the pre-1960 values.

Figure 9 .
Figure 9. Box plots for the (a) northern and (b) southern IGB sub-basins post-and pre-1960 Dmax percentage difference and (c) and (d) positive difference percentage area analysis for ∆Dmax at 2-, 10-and 100-year return period frequencies over IGB.The percentage difference is calculated by subtracting the pre-1960 values from the post-1960 values and then dividing the result by the pre-1960 values.

Figure 10 .
Figure 10.Combined analysis of the difference percentage of the Vmax, Rmax, and Dmax of the northern and southern region sub-basins at (a) and (d) 2-, (b) and (e) 10-and (c) and (f) 100-year return period frequencies.The percentage difference is calculated by subtracting the pre-1960 values from the post-1960 values and then dividing the result by the pre-1960 values.
(Lin and Huybers (2019)the India meteorological department (IMD; Pai et al 2014) for 1901-2019.This high-resolution rainfall dataset was developed from 6995 observed rain gauges nationwide.However, it was observed that the number of stations used to construct this dataset varied over time.For example, around 1500 and over 6000 stations were used to prepare this dataset in 1901 and 2019, respectively.This varying station density led to artificial multi-decadal variability of mean rainfall and artificial variability of extreme rainfall(Lin and Huybers (2019), Singh et al (2019), Zahan et al (2021)).Further, IMD also provides a century-long (1901-2022) daily rainfall dataset developed by Rajeevan et al (2008) based on a fixed station network (1380 stations) gridded at 1 • × 1 • resolution.Rajeevan et al (2008) dataset used a fixed number of rain gauges; therefore, station density is relatively lower compared to the latest high resolution (0.25 • × 0.25 • ) proposed by Pai et al (2014).Pal et al (2021) further compared the 0.25 • (Pai et al 2014) and 1 • (Rajeevan et al 2008) rainfall datasets and concluded that rainfall trends present in both datasets are reasonably similar for the entire Indian region, except Jammu and Kashmir regions.The observed inconsistency between the two gridded datasets over the Jammu and Kashmir regions was primarily induced by very low rain gauge station density in the Rajeevan et al (2008) dataset (Pal et al 2021).Furthermore, Pai et al (2014) compared the gridded IMD rainfall dataset proposed by