Flood frequency analysis and inundation mapping for lower Narmada basin, India

A flood is defined as any moderately high-water flow that overtops the artificial or natural banks in any portion of a waterway or stream (Şen 2018). Depending on the cause of the flood, it can be distinguished in different types such as flash floods, river floods, urban floods, coastal floods, and flooding due to dam breaking (Şen 2018; Kundzewicz et al. 2019). Floods are one of the most frequent hazards globally, which is further aggravated due to climate change impacts and human-induced activities (Arnell & Gosling 2016). The rising trend of global average temperature has led to a disturbing rainfall pattern of no precipitation for longer periods, followed by a sudden high-intensity excessive precipitation, resulting in extreme events that have a devastating impact on human life and livelihood (Tabari 2020). India is highly vulnerable to flood (Mohapatra & Singh 2003). Climate change has had a significant influence on India, making it one of the most flood-prone countries globally (Guhathakurta et al. 2011). Out of 329 million hectares of the total geographical area, more than 40 million hectares are prone to flood; that is, nearly 12.16% of the total area is susceptible to flooding (National Disaster Management Authority 2008). According to the CWC's (Central Water Commission, Govt. of India) report, 107,535 casualties occurred due to heavy rains and floods, and 53,425.50 million USD (INR 378,247.047 crores) worth of public utilities, houses, and crops have been damaged across India over 64 years, from 1953 to 2017 (Central Water Commission 2018).

In recent years, climate change induced by human activities has affected rainfall patterns and distribution, causing floods in areas that were not flood-prone earlier (Yadav & Mangukiya 2021). Moreover, the construction of infrastructures in a floodplain near the bank of a river encroaches the river's width and causes a reduction in flow carrying capacity (Patel et al. 2017). Thus, it is essential to evaluate the flow carrying capacity of the river for the prediction of a flood. Different techniques and analysis systems are available in the literature for estimating the flow carrying capacity of the river (Ng et al. 2018; Yan et al. 2018; Mehta & Yadav 2020; Mangukiya & Yadav 2021); among them, the Hydrological Engineering Centre – River Analysis System (HEC-RAS) is most widely used (US Army Corps of Engineers 2010). HEC-RAS is free open-source computer software that simulates one-dimensional (1D) and two-dimensional (2D) water flow in natural rivers and other channels (Brunner 2010). A 1D hydraulic model simulates water level and flow along a river (say, in the x-direction), whereas a 2D hydraulic model simulates water level and flow along a flood plain (i.e., x- and y-directions). Various studies have demonstrated the capabilities of HEC-RAS for analyzing river carrying capacity and flood inundation mapping using 1D and 2D hydraulic models (Cook & Merwade 2009; Timbadiya et al. 2014; Derdous et al. 2015; Papaioannou et al. 2016; Patel et al. 2017; Vora et al. 2018).

One of the most significant difficulties in hydrology is better understanding flood regimes. For this purpose, civil engineers and hydrologists most commonly use flood frequency analysis (FFA), which estimates flood peak values based on non-exceedance probabilities (Tanaka et al. 2017). The findings in the FFA application are theoretically valid if the series are independent and identically distributed (Kidson & Richards 2005). Various statistical distribution methods can do flood frequency analysis (Cunnane 1988). The most commonly used methods are Gumbel's distribution (Yue et al. 1999), Log-normal distribution (Hoshi et al. 1984), and Log-Pearson Type-III distribution (Phien & Ajirajah 1984). The Gumbel's distribution represents an enormous value from a reasonably large group of independent values from distributions with relatively rapidly decaying tails, such as exponential or normal distribution. The Log-Pearson Type-III distribution is a statistical approach for predicting the design flood at a particular site by fitting frequency distribution data. The Gumbel's Extreme Value Type-I (EV-I) and Log-Pearson Type-III distribution are commonly used by various federal agencies like US Geological Survey (USGS) and Federal Emergency Management Agency (FEMA) for FFA (Rahman et al. 2014). Several studies have been conducted to investigate the application of these statistical distribution approaches in FFA (Odry & Arnaud 2017; Onen & Bagatur 2017; Thorarinsdottir et al. 2018; Baidya et al. 2020).

The studies on FFA of Indian river basins are limited due to observed data scarcity in river discharge and water level (Yadav & Mangukiya 2021). India is one of the most flood-prone countries globally; flood studies on a regional level can help reduce its impact (Mohapatra & Singh 2003). The broad objective of the present study is to find out the potential area at risk of flood in the Lower Narmada Basin, India. The FFA was carried out using Gumbel's EV-I and Log-Pearson Type-III method. The percentage error between observed peak flood and calculated peak flood from FFA were computed to find out the best approximation of extreme flood with a return period of 10, 25, 50, and 100 years. This approximate value of extreme flood was simulated using HEC-RAS to prepare a flood inundation map and identify the potential area under flood risk. The present study demonstrates the framework for projections of the FFA in a concise way that can be adopted in various river basins around the globe. The projections of FFA and associated flood inundation maps can be directly helpful for preventing future floods and increasing the city's resilience for flood events. The obtained results from the present study can help the authorities and policymakers prioritize the flood mitigation measures and decide the region's development policy.

The Narmada River is the fifth-largest and sixth-longest river of India. The Narmada Basin is bounded by the Vindya and Satpura ranges, covers an area of 98,796 km², and is located between east longitudes 72°38′ to 81°43′ and north latitudes 21°27′ to 23°37′. The basin is divided into five distinct physiographic areas. They are: (1) the upper hilly areas, which encompass the districts of Balaghat, Durg, Mandla, Seoni, and Shahdol; (2) the upper plains, which encompass the districts of Betul, Chhindwara, Damoh, Hosangabad, Jabalpur, Narsinghpur, Raisen, and Sehore; (3) the middle plains, which encompass the districts of Dewas, Dhar, Indore, Khandwa, and part of Khargone; (4) the lower hills areas, which includes parts of Dhulia, Jhabua, Narmada, Vadodara, and west Nimar; and (5) the lower plains in the coastal region, which include mostly the Bharuch, Narmada, and Vadodara districts. The lower plain of the coastal region is selected as a study area for the present study (Figure 1). The lower Narmada basin has experienced numerous floods, including severe flooding in 1970, 1973, 1984, 1990, 1994, 2006 and 2013. The upper hilly areas of the basin receive higher annual rainfall (1,400 to 1,650 mm), which causes floods in the downstream region even though it is a semi-arid zone. The temperature in the lower part is influenced by the sea and varies from 10 to 40 °C in different seasons. The primary land-use land-cover class of the Narmada basin is agriculture cropland (60%), barren land (12%), and urban land (3%).

The required data for the present study were obtained from open-source data portals and government reports. The Shuttle Radar Topography Mission (SRTM) Digital Elevation Model (DEM) of 1-arc second resolution for the study area was downloaded from the United States Geological Survey (USGS) Earth Explorer portal (https://earthexplorer.usgs.gov/). The Moderate Resolution Imaging Spectroradiometer (MODIS) data of land-use land-cover is used for the present study (https://modis.gsfc.nasa.gov/data/). The stage-discharge data for the Garudeshwar site was collected from the Sardar Sarovar Narmada Nigam Limited (SSNL) authorities. The historical flood information of the study area was compiled from various government reports published by SANDRP (South Asia Network on Dams, Rivers and People), CWC (Central Water Commission, India), and SSNL authorities (SSNL, Gandhinagar).

The methodology includes the data collection and pre-processing, flood frequency analysis for Garudeshwar weir, and development of the hydraulic model downstream of the Garudeshwar weir to identify the area at flood risk. Figure 2 shows the flowchart of the methodology adopted for the present study.

Flood frequency analysis (FFA) is a technique used by hydrologists to predict flow levels that correspond to specified return intervals or probabilities along a river. FFA is used to derive statistical information such as standard deviation, mean, and skewness, which is then utilized to prepare frequency distribution plots using annual maximum flow data that has been available for many years. FFA is measured using a variety of ways. Gumbel's EV-I and log-Pearson Type-III techniques for FFA are utilized in this study for return periods of 10, 25, 50, and 100 years. Annual Peak discharge from the year 1948 to 2016 is considered as parameter for FFA.

The 2D HD model was prepared using HEC-RAS v6.0 for simulating the extreme flood in the Lower Narmada Basin. In the absence of the surveyed terrain data, the SRTM DEM of 1-arc second resolution was used as digital terrain data for the region. The pre-processing of the DEM was carried out using QGIS desktop v3.20 for filling the sinks and removing peaks. The sinks and peaks are relatively low and high values of the cell than the surrounding cells, which arise due to the resolution of the data. The filled SRTM DEM was then given as input in HEC-RAS Mapper as digital terrain of river and floodplain. The computational mesh of cell size 50 m in the X and Y direction was generated to simulate hydraulic parameters in a 2D floodplain area. The observed discharge data of the year 2013 at the Garudeshwar weir was given as the upstream boundary condition, and the normal slope of the river was given as the downstream boundary condition. The roughness coefficient (Manning's n) is one of the key inputs for calculating the friction applied to flow by the land surface. In absence of the surveyed flood plain roughness values, the land-use land-cover map of MODIS data was given as input for the roughness coefficient of the floodplain. For the roughness coefficient of river bed, the 1D HD model was calibrated and validated by comparing the simulated and observed data at the downstream gauging station (Bharuch). The calibrated and validated HD model was then used to simulate extreme flood events of return periods 10, 25, 50, and 100 years. The simulated inundation maps were then exported and overlaid with land-use maps and Google Earth maps in QGIS desktop v3.20 to identify and calculate the area at flood risk. The calibrated and validated HD model can also be used for the region's flood early warning system, providing real-time or forecasted data as input to the upstream boundary condition.

The observed annual maximum discharge data of 69 years (1948–2016) was used for fitting the statistical distribution of Gumbel's EV-I and Log-Pearson Type-III. The approximated extreme flood of return periods 10, 25, 50, and 100 years were simulated using HEC-RAS to identify potential flood risk areas. The results of the obtained extreme flood and inundation map are described in sections 4.1 and 4.2, respectively.

The mean and standard deviation of the annual peak discharge data series were computed first. Then the reduced variant was calculated for return periods 10, 25, 50, and 100 years. The reduced mean and standard deviation values were selected from Table S1 based on the N value. Then the magnitude of flood for return periods 10, 25, 50, and 100 years was calculated using Equation (1). The annual maximum discharge data series was converted to a logarithmic (base 10) form to fit the Log-Pearson Type-III distribution. Then the coefficient of skew was calculated. The value frequency factor was obtained from Table S2 based on the skew coefficient and return period values. Then the magnitude of flood for return periods 10, 25, 50, and 100 years was calculated using Equation (2) and transformed in antilog form. The calculated parameters of the Gumbel's EV-I and Log-Pearson Type-III distribution are shown in Table 1. The estimated extreme values of peak flood for return periods 10, 25, 50, and 100 years using Gumbel's EV-I and Log-Pearson Type-III methods are shown in Tables 2 and 3, respectively.

Based on the obtained results from both methods, the flood frequency curve was fitted using a logarithmic relationship (Figure 3). The result shows that both approaches have high accuracy to be used for FFA. To evaluate further, the obtained result from both methods was compared with the fitted curve line of observed annual peak flood data, and percentage deviation was calculated (Figure 4). This comparative analysis shows that the Log-Pearson Type-III method gives a more realistic estimation of peak floods for less return period (T < 60 years). In comparison, Gumbel's EV-I method provides a more realistic estimation of peak floods for a high return period (T > 60 years). Keeping this in view, the estimated flood peak from the Log-Pearson Type-III method corresponds to 10-years, 25-years, and 50-years, and the estimated flood peak from Gumbel's EV-I method corresponds to 100-years considered for generating flood inundation maps (Table 4).

The 2D HD model was prepared in HEC-RAS v6.0 using SRTM DEM, MODIS land-use land-cover data, and observed discharge at the Garudeshwar weir. The observed discharge data of the year 2013 at the Garudeshwar weir (peak discharge 32056 cumec) was given as upstream boundary condition, and the normal slope of the river was given as downstream boundary condition. For roughness coefficient of river bed, the 1D HD model was calibrated by changing the roughness coefficient of the river within the range of 0.02 to 0.03 (Chow 1959) and comparing the simulated and observed data at the Bharuch gauging station. The value of roughness coefficient as 0.022 was selected with the minimum absolute error of 0.3 m between maximum observed and simulated water level. The validation of the 1D HD model was carried out by comparing the observed and simulated water level of independent flood event of 2006. The calibrated and validated HD model was then used to generate the flood inundation maps corresponding to return periods of 10, 25, 50, and 100 years (Figure 5). The result indicates that the inundation area with a high flood depth that significantly increases with the return period. The depth-wise inundation area analysis shows that the inundation area with less than 1 m water depth remains almost constant (50 km²) for different T; the inundation area with a water depth of 1 to 3 m decreases with an increase in T; the inundation area with water depth 3 to 6 m increases initially and then decreases with increase in T; the inundation area with water depth more than 6 m significantly increases with increase in T (Figure 6).

Further, the inundation area of different land-use classifications was obtained from the flood inundation maps. It shows a significant increment of inundated cropland and urban land with an increase in T (Table 5), which indicates that the flood vulnerability in the area is significantly increasing with an increase in T. Based on the obtained flood inundation map, the villages under risk were identified (Figure 7). The level of risk was decided based on the flood depth and corresponding return period T. The high flood depth with less return period was considered a high-risk area, while the low flood depth with a high return period was considered a low-risk area. The identified risk of different villages (Table 6) will be helpful for authorities and policymakers to decide the flood prevention and mitigation measures and policy.

The present study demonstrates the framework for projections of the FFA in a concise way that can be adopted in various river basins around the globe. In the present study, the flood frequency analysis (FFA) was carried out using Gumbel's EV-I and Log-Pearson Type-III method for the Lower Narmada Basin, India. The percentage error between calculated peak flood from FFA and observed peak flood was computed to find the best approximation of extreme flood with a return period of 10, 25, 50, and 100 years. The result shows that the Log-Pearson Type-III distribution method best approximates the extreme event of the lower return period (T < 60 years). In comparison, Gumbel's method gives the best approximation of the extreme event for the higher return period (T > 60 years). The two-dimensional (2D) hydrodynamic (HD) model was prepared using HEC-RAS v6.0 for the Lower Narmada Basin. The roughness coefficient value of 0.022 for the river bed was calibrated and validated by 1D HD model with a minimum absolute error of 0.3 m. The calibrated and validated HD model was used to simulate estimated peak floods of 47,205.65 cumec for 10-year, 59,773.21 cumec for 25-year, 68,839.94 cumec for 50-year, and 77,859.72 cumec for the 100-year return period. The result indicates that the inundation area with a high flood depth significantly increases with the return period. The depth-wise inundation area analysis indicated that the inundation area with depth 3 to 6 m decreases with an increase in return period while the inundation area with depth more than 6 m increases with an increase in the return period. The inundation area of different land-use types shows a significant increment of inundated cropland and urban land with an increase in return period, which indicates that the flood vulnerability in the area is significantly increasing with an increase in the return period. Based on the obtained flood inundation map, the villages under risk were identified, which shows that Bharuch and Ankleshwar (two major cities in the study area) are at high and medium flood risk, respectively. The projections of FFA and associated flood inundation maps can be directly helpful for preventing future floods and increasing the city's resilience for flood events. The flood inundation projections from the FFA and associated inundated land-use land-cover class information can assist the authorities and policymakers in prioritising the flood mitigation measures and deciding the region's development policy.

The authors are not affiliated with or involved with any organisation or entity with any financial interest or non-financial interest in the subject matter or materials discussed in this paper.

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

The authors are thankful to the Sardar Sarovar Narmada Nigam Limited (SSNNL), Gandhinagar, and Central Water Commission (CWC), Surat, for providing necessary data for the study reported in the paper.

All relevant data are included in the paper or its Supplementary Information.