Trends and Non-Stationarity in Groundwater Level Changes in Rapidly Developing Indian Cities

: In most of the Indian cities, around half of the urban water requirement is fulﬁlled by groundwater. Recently, seasonal urban droughts have been frequently witnessed globally, which adds more stress to groundwater systems. Excessive pumping and increasing demands in several Indian cities impose a high risk of running out of groundwater storage, which could potentially a ﬀ ect millions of lives in the future. In this paper, groundwater level changes have been comprehensively assessed for seven densely populated and rapidly growing secondary cities across India. Several statistical analyses were performed to detect the trends and non-stationarity in the groundwater level (GWL). Also, the inﬂuence of rainfall and land use / land cover changes (LULC) on the GWL was explored. The results suggest that overall, the groundwater level was found to vary between ± 10 cm / year in the majority of the wells. Further, the non-stationarity analysis revealed a high impact of rainfall and LULC due to climate variability and anthropogenic activities respectively on the GWL change dynamics. Statistical correlation analysis showed evidence supporting that climate variability could potentially be a major component a ﬀ ecting the rainfall and groundwater recharge relationship. Additionally, from the LULC analysis, a decrease in the green cover area (R = 0.93) was found to have a higher correlation with decreasing groundwater level than that of urban area growth across seven rapidly developing cities.


Introduction
Groundwater, one of the Earth's largest available freshwater resources, is being extracted at a rate of~982 km 3 /year [1,2]. Globally, groundwater meets about half of domestic water needs and 38% of irrigation water demand for sustaining food security [3,4]. The increasing global population is a major concern since presently, about 80% of the world's population is at an incident threat to water security [5]. Unlike developed countries, developing countries like India are at high risk of water activities affect water availability [8], the high dependency of groundwater in larger urban centers of India has created an acute supply-demand deficit [29]. Recent studies also suggest that the frequency of droughts over the massive urban environments is likely to increase throughout the Indian cities in the near future [30,31].
The migration to cities from rural areas is exorbitant and likely to increase in the future. Thus, the projected rapid urban growth is concentrated particularly in secondary rapidly developing cities. This would consequently lead to unprecedented challenges and uncertainties in meeting water needs in the future. While there is a little scope for improving the water infrastructure in developed cities conceivably due to feasibility and implementation issues, proper planning and management of water systems are crucial in secondary developing cities.
In the Indian context, no study reported the trends in groundwater level (GWL) changes relating to non-stationarity neither for developed cities nor secondary cities. However, in other parts of the world, the spatio-temporal changes and non-stationarity in groundwater levels have been analyzed using several statistical methods [32][33][34][35]. The major findings from these studies are (i) majority of the wells exhibited the non-stationarity with significant trends in GWL changes, and (ii) precipitation and recharge of groundwater levels are affected by climate variability (both on an Annual and Decadal scale), which has more impact than the increasing temporal patterns in pumping [32][33][34][35].
While previous studies focused on assessing the GWL trends influenced by the land use/land cover (LULC) changes at multiple spatial scales [36,37], this paper focused on a comprehensive analysis of GWL change dynamics considering the climate variability (rainfall) as well as LULC change (anthropogenic activity), especially in rapidly developing cities, and discussed the potential driving forces responsible for non-stationarity. For this endeavor, seven densely populated and rapidly growing secondary cities that share similar socio-economic developmental goals were selected from different geographic regions across India for assessing the condition of groundwater systems. Various statistical analyses, widely employed in the literature, such as Theil-Sen estimator [38], modified Mann-Kendall test [39], Augmented Dicky Fuller test [40], Phillip-Perron test [41], and Kwiatkowski-Phillips-Schmidt-Shin test [42], were used to determine the trends and non-stationarity in the groundwater level (GWL) changes [32,33,35,43,44]. The specific objectives of the paper are to (i) provide an overall outlook of trends in GWL change dynamics, (ii) compare the GWL changes with the Total Water Storage (TWS) anomalies, (iii) determine the non-stationarity in the GWL changes, and (iv) assess the contribution of climate variability (rainfall) and LULC change to the non-stationarity.

Study Area and Data Collection
Seven cities (Table 1) were selected based on the following criteria: (i) population~1 Million (2011 Census) [45], (ii) decadal population growth rate greater than~15% (2001 to 2011 Census) [45], (iii) high Gross Domestic Product (GDP) and economic growth rate [46], (iv) population density, and (v) shortlisted as a smart city under the urban renewal program, National Smart City Mission, by the government of India [47]. The groundwater level data for the years 1996-2018 has been collected for the selected seven cities from the Central Groundwater Board (CGWB) through the India Water Resource Information system [48]. The CGWB monitors seasonal groundwater levels quarterly (i.e., January-Post-monsoon Rabi, May-Pre-monsoon, August-Monsoon, and November-Post-monsoon) [49]. Grubbs test was used to remove the outliers in the dataset.

The Grubbs test statistics are defined as
, where Y i is the ith data point and Y and s denote the sample mean and standard deviation, respectively. The wells located in the selected cities having at least 18 years of observed data out of a total of 23 years (~80% observations) were selected for the analysis (refer to Figure 1 for the map of the study area and the wells selected for the analyses) [50]. The monthly Terrestrial Water Storage anomalies were obtained from Gravity Recovery and Climate Experiment (GRACE) satellite observations for the months of January, May, August, and November within the study period (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017) to analyze the total water storage in the selected cities [51]. Daily gridded rainfall data (0.25 ø × 0.25 ø ) obtained from the Indian Meteorological Department (IMD) were used to estimate the monthly rainfall values for the selected cities [52]. Moderate Resolution Imaging Spectroradiometer (MODIS) classified LULC maps (MCD12Q1v006: LC_type 4 layer) obtained for the years 2001 and 2018 were used to perform the LULC analysis [53].

Trend Analysis
The Theil-Sen estimator, a widely used method in detecting trends in a hydrological time series, was employed to determine the nonparametric linear GWL trend (in cm year −1 ) [38,54], and the significance of the determined trend was estimated by applying the modified Mann-Kendall (MMK) test [39,55]. The Theil-Sen estimator determines the slope between all possible data pairs as, where Q i is the slope between data points y j and y k at time j and k, such that j > k and i = 1, 2 . . . n.
if N is even [38]. The MMK test uses a modified variance to reduce the influence of autocorrelation on the Mann-Kendall test results [56,57]. The test statistics of MMK are the same as the Mann-Kendall test and are defined as S = n−1 i=1 n j=k+1 Sign x j − x k , where S is the test statistics, x j and x k are the data points observed at time j and k respectively ( j > k), and n is the length of the dataset [39,56]. The trends in the GWL were estimated in terms of changes in the depth to groundwater table below ground level in centimeters per year. Although several methods exist to recover the GWL from Satellite products (e.g., GRACE [51]), these methods provide only gross estimates with higher uncertainty, bias, and errors [58,59]. To overcome the above issues, GWL observations from monitoring wells by CGWB were used. Though the trend analysis was focused on estimating the temporal trends at the point level (individual wells), the selected wells are uniformly distributed across the city ( Figure 1), with each well representing a block (or part) of a city, which subsequently helps to understand the spatial trends as well.

TWS Anomalies and Groundwater Levels
To understand the extent of changes in the overall water storage (includes surface water, soil moisture, groundwater, etc.), the Terrestrial Water Storage (TWS) anomalies were analyzed. The TWS anomalies data were obtained from the Gravity Recovery Climate Experiment (GRACE) interpolated at 1 ø × 1 ø grid cells for the selected seven cities for the period 2002 to 2017 [51]. The trends and non-stationarity in the TWS changes were determined, and the results of TWS were correlated with GWL changes in the selected cities.

Non-Stationarity Analysis
Non-stationarity analyses were performed to assess the impacts of climate change caused by low-frequency climate variability and LULC changes on GWL changes [60]. The most frequently used Augmented Dickey-Fuller (ADF), Phillip-Perron (PP), and Kwiatkowski-Phillips-Schmidt-Shin (KPSS) tests were applied to estimate the non-stationarity in GWL changes [40][41][42].
The ADF test assesses the null hypothesis of a unit root using the model: y t = y t − 1 + β 1 ∆y t − 1 + β 2 ∆y t − 2 + . . . + β p ∆y t − p + ε t against the alternative model: where y t is the time series with time t, ∆ is the differencing operator, β is constant, p is the number of lagged difference terms, and ε t is a mean zero innovation process [40,61].
The alternate hypothesis (φ < 1) in the ADF test is generally stationary or trend stationary [61].
Similar to the ADF test, the PP test assesses the null hypothesis of a unit root in a univariate time series; however, it makes a nonparametric correction to t-statistics, which makes it more robust with unspecified autocorrelation [62,63]. The PP test assesses the null hypothesis using the model: y t = y t−1 + e(t), against the alternative model: y t = a y t−1 + e(t), where y t is the time series with time t, e(t) is the innovations process and a is the auto regression coefficient, such that a < 1 [41,63].
The KPSS test, in contrast to ADF and PP tests (which tests for non-stationary as a null hypothesis), assesses the null hypothesis that a univariate time series is trend stationary against the alternate hypothesis of non-stationarity [42]. The KPSS test assesses the null hypothesis of trend stationarity using the models: y t = c t + δ t + u 1t and c t = c t − 1 + u 2t against the alternate hypothesis of non-stationarity, where c t is the random walk term, δ is the trend coefficient, u 1t is a stationary process, and u 2t is an independent and identically distributed process with mean 0 and variance σ 2 [64]. The null hypothesis σ 2 = 0 implies that the random walk term (c t ) is constant, and the alternative hypothesis σ 2 > 0, which introduces the unit root in c t [42,64].
Although the power of non-stationarity tests varies, all the non-stationarity tests are known to have low power (power is determined by auto-regression parameter φ 1 in case of ADF and PP tests) when the length of the time series is short [65]. However, the KPSS test is an exception and performs well in shorter time series, which made it suitable for our analysis (although the KPSS test has a high rate of Type I Errors (rejection of null hypothesis)) [65,66]. As the availability of the observed GWL data is very limited, the non-stationarity in the GWL changes was confirmed if the majority of the three tests (ADF, PP, and KPSS) suggest non-stationarity. Thus, the bias due to data length was minimized.

Correlation Analysis: Rainfall and Groundwater Levels
The impacts of climate variability on GWL changes were grossly analyzed by correlating rainfall with GWL changes in the wells having non-stationarity [67]. The relationship between rainfall-groundwater recharge is complex; however, few studies show a high correlation between them [68][69][70]. Since the groundwater system's recharge characteristics are complex, there is always a time delay between rainfall and the water to reach the groundwater table [68,71]. The cross-correlation was carried out for two conditions such as: (i) no lag and (ii) 3-month lag (or 1 season lag) with the groundwater level data to account for the delay in the recharge of groundwater after the rainfall. For this, the gridded rainfall data and GWL changes from the wells exhibiting non-stationarity behavior and falling within the grid were used. It may be noted that about 94% of the wells used in the analysis are dug wells (lies predominantly in the unconfined aquifers) with depth generally varying between 10 and 30 m, and therefore directly relating the monthly rainfall to GWL changes is justified (not accounting for the recharge from the confined aquifer) [49].

LULC Change Analysis
The impacts of LULC change on GWL were also quantified. The Moderate Resolution Imaging Spectroradiometer (MODIS) classified LULC maps (MCD12Q1v006: LC_type 4 layer) for the years 2001 and 2018 were used to analyze the LULC changes in the respective cities. This MODIS-MCD12Q1v006 data product for global land cover is available at 500 m resolution [53]. A 15 km radial buffer around the center point of the cities was created since the overlay analysis revealed that the urban expansion between 2001 and 2018 at these selected cities was very much contained within this buffer zone. This 15 km buffer zone shapefile was used to clip out the LULC maps of the respective cities. Google Earth Engine (GEE) algorithms were used to analyze the entire spatial data, starting from accessing the data to the creation of the LULC change map. The impact of change in LULC components on GWLs was assessed through performing correlation analysis between the percentage of wells having significant increasing/decreasing GWL trends (>2 cm/year) and LULC change between 2001 and 2018 at the city level. The statistical analyses were also performed for the wells located outside the urban boundary as these cities have high scope for expansion in terms of infrastructural growth and socio-economic activities.

Trends in Groundwater Level Changes
Trend analysis for GWL changes revealed that the GWL majorly varied between ±~10 cm/year ( Figure 2) across the selected cities. Spatially, GWL trends ( Figure 3) were found to be decreasing in the majority of the wells, except in Rajkot and Guntur. In general, consistent high decreasing and low increasing GWL trend patterns were observed in January relatively, compared to other months of observation ( Figure 2). A high increasing GWL trend was found in May (Figures 2-4), though May is considered to have peak summer with high evapotranspiration losses and very few rain spells in India [72,73]. This anomaly, however, can be attributed to changes in the characteristics and spatial-temporal variability of rainfall events in India [74,75]. Very similar GWL trends across the seasonal data collected in January, May, August, and November ( Figure 2) were observed, which suggests a less interdependence between seasonality and GWL trends. However, the lower tail of the box plot ( Figure 2) indicates a sharp increase in GWL, and this behavior could be attributed to groundwater recharge from the monsoon rainfall [15,16]. Similarly, the upper tail of the box plot indicates a sharp decline in GWL that might be due to excessive pumping [16]. Overall, a higher magnitude GWL trend (increasing as well as decreasing) greater than~25-100 cm/year throughout the seasons from January to November was observed, which might be due to the effect of climate variability and/or anthropogenic activities [76][77][78]. Figure 2 also illustrates the magnitude of GWL variability among the selected cities. A higher magnitude of decreasing GWL trends was noticed in Allahabad (especially in August and November), followed by Jodhpur. However, in Jodhpur, a high-magnitude increasing trend was also found in some proportion of the wells throughout the seasons (January to November). A high magnitude of increasing GWL trend was detected in Rajkot, which confirms a very high rate of increasing GWL in the western arid regions of India, as reported in previous studies [16,79]. A similar magnitude of increasing and decreasing GWL was observed in Guntur from January to November. However, in May, a relatively higher decreasing GWL trend in Guntur was found, with more than 20% of wells declining at a rate greater than 10 cm/year. This might be due to a high surge in pumping during the peak of a seasonal dry spell since Guntur's economy is highly dependent on agricultural and allied activities [80]. Tiruchirappalli and Solapur exhibited similar behavior in terms of resulting slightly higher magnitude of declining GWL. Except for Rajkot, where the GWL was increasing in the majority of wells betweeñ 2 to 20 cm/year, the overall groundwater level was found to be increasing or decreasing at~10 cm/year in all other selected cities.   . City-wise seasonal split-up of the wells with groundwater level trends. In majority of the selected cities, the decreasing GWL trend was found to be higher in majority of the wells throughout the season. In about~9-47% of the wells, GWL changes occurred only between ±2 cm/year throughout the seasons. At the city level, Rajkot has the highest percentage of wells with increasing GWL and Allahabad has the highest percentage of wells with decreasing GWL with a significant trend (>2 cm/year).
There were no significant GWL changes ( Figure 4) in about~7-29% (with exception to Rajkot during January (4%), August (6%), and November (4%), and Raipur in May (4%) and August (47%)) of the wells throughout the study period (we assumed that there is no significant GWL trend if the GWL is varying at ±2 cm/year) in all these selected cities. In the selected cities, increasing and decreasing GWL trend was found in~8-79% and~17-64% of wells respectively, from January to November. The city-level results suggested that Rajkot has the highest (60-79%) percentage of wells with a significant increasing GWL trend (greater than 2 cm/year). In contrast, Allahabad has the highest percentage of wells (55-64%) with a significant decreasing GWL trend (greater than 2 cm/year). Raipur was determined to have an almost equal percentage of wells with significant increasing and decreasing GWL in May and November. Overall, the percentage of wells with a decreasing trend was found to be higher than the increasing trend. This reveals groundwater depletion is evident in most of the developing cities across India, which may threaten urban water security [8].
Further, analyses were conducted to find the percentage of wells having extreme increasing or decreasing GWL trends ( Figure 5). Since the Theil-Sen slope estimates the median trend from a combination of slopes, it is important to note that the magnitude of the actual trend (increasing or decreasing) might be much higher than the trend estimated using the Theil-Sen slope method. Based on the thresholds of~25 and~50 cm/year, the percentage of wells having an extreme GWL trend during the study period was determined. In the first half of the year (i.e., January and May), the percentage of wells with increasing GWL was found to be higher than that of decreasing GWLs in the selected cities. However, at the end (November), opposite behavior was observed for which the percentage of wells having a decreasing GWL was high. This might be due to the excessive rate of pumping in November, mainly to manage the water shortage caused by the lack of rainfall during the post-monsoon season across India [74]. Figure 5. Percentage of wells in the seven cities with extreme (high magnitude) groundwater level trends for January, May, August, and November. About more than 5% and 2% of the wells in the selected cities have extreme increasing and decreasing GWL trends respectively, from January to November.

Non-Stationarity and Significance of Groundwater Level Trends
The significance of GWL trends was estimated using the MMK test. Based on the results (Figure 6), we found that~40-93% of the wells in the selected cities had a significant trend during the different months of study (i.e., a monotonic trend is present, and the trend is statistically significant), with the highest MMK significance observed consistently in November (60-91% of wells showed a significant trend). Similarly, the significance of the trend at the city level was found to be in higher proportions in Solapur (79-95% of wells), followed by Tiruchirappalli (80-93% of wells, except in August). The significance of trend was found to be relatively minimal in Jodhpur and Guntur.
The combined results ( Figure 6) of ADFT, KPSS, and PP tests showed that~70-100% of the wells exhibited non-stationarity across different months (i.e., January to November). Notably, in November, about 90-100% of the wells ( Figure 6) exhibited non-stationarity. At the city level, non-stationarity was observed in almost all the wells in Rajkot (98-100%) from January to November. On the other side, non-stationarity is least in Solapur, with 33% of wells (observed in August), followed by Jodhpur, which has non-stationarity in only 54% of wells during May and August. It can be understood from these results that the majority of the wells in these cities (at least more than 50%) may be significantly affected by climate variability and LULC change. The presence of such high non-stationarity would create huge uncertainty in the behavior of groundwater systems, which makes it difficult to predict future GWLs. Subsequently, this might lead to complex problems in different spheres, especially in development, public health, and irrigation, in the selected rapidly growing secondary cities [81]. About~27-33% of wells in the selected cities were determined to have a decreasing GWL trend along with MMK trend significance and non-stationarity in GWL changes. At the city level, this behavior was found to be highest in Tiruchirappalli (64% of wells) and least in Rajkot (12% of wells). Notably, Jodhpur, the city with high decreasing GWL, has only~20% of wells with this behavior. Interestingly, we determined that there was no non-stationarity as well as the MMK trend significance in wells present in the regions of high decreasing GWL trend, i.e., Rajkot and Allahabad. Although these cities (i.e., Rajkot and Allahabad) have high increasing and decreasing GWL (Figures 4 and 5, respectively), the percentage of wells having non-stationarity and MMK trend significance were relatively lower among the selected cities.

Terrestrial Water Storage Anomalies
The Theil-Sen trend in TWS anomalies (Figure 7) was estimated, and the highest declining TWS trend was observed in Jodhpur (−1.79 cm/year) ( Table 2). On the contrary, TWS is highly increasing in Rajkot at 0.67 cm/year among the selected cities. While high non-stationarity was found in the GWL changes, non-stationarity in TWS anomalies was not observed in any of the selected cities. Similarly, MMK trend significance was also found to be absent in cities with a higher magnitude of increasing or decreasing TWS trends (Rajkot and Jodhpur). Unlike the GWL changes, the TWS trends indicate that the overall water storage or the TWS cycle (including surface water, soil moisture, etc., including groundwater) in the selected secondary cities have not majorly been affected by climate variability and LULC changes [51]. However, comparing the changes in TWS anomalies with GWL, a high positive correlation was observed both in the case of increasing (R = 0.85) and decreasing (R = 0.87) trends. Refer to Appendix A Figure A1 for TWS trends of all the selected cities.

Rainfall and Groundwater Level Change Relationship
The response of GWL to rainfall was investigated through correlation analysis for the period 1996-2018. The correlation analysis on the wells having non-stationarity in the selected cities revealed that there is no strong dependency between rainfall and GWL changes. Figure 8 shows the box plot of the Correlation Coefficient (R) values for all the wells with non-stationarity. The positive R-value indicates that the GWL increases with rainfall and vice versa. The R-value varying between ±0.2 was found in~95% of the wells for both no lag and three months lag (with the exception of Allahabad in January (under 3-month lag scenario) and Solapur and Rajkot in August (under no lag scenario)), as illustrated in Figure 8. In the 3-month lag scenario (Figure 8a), negative R values were observed in the majority of the wells, notably in Rajkot (both in January and November). A similar observation was noted for the no lag scenario (Figure 8b); however, the tail ends were shorter in the box plot in August and November. Since the majority of the wells (dug wells with depth~10-30m) considered in this paper are shallow wells, the recharge time delay is comparatively small. Therefore, no lag scenario might be considered as a realistic scenario [71]. Although from the results of correlation analysis, no good rainfall-GWL response relationship could be drawn-out in the selected cities for January, May, and November, the R values were relatively high in Solapur, Rajkot, and Raipur in the month of August in the no lag scenario. While there is comparatively a good correlation between rainfall and GWL observed in August in the above cities, the percentage of wells with significant declining GWL trends was observed only in 20-36% of wells (Figure 4). In contrast to August, poor rainfall-GWL response was observed in January in which almost all the wells in Rajkot, Guntur, and Allahabad showed a negative correlation. January, which has a poor R-value in the no lag scenario (Figure 7b), has a relatively higher proportion of wells with decreasing GWL in Allahabad, Jodhpur, and Solapur (Figures 3 and 4), which suggests that declining GWL have a strong correlation with decreasing recharge from rainfall. The correlation analysis suggests that the rainfall and groundwater recharge relationship has been seriously affected and also ascertains the lack of seasonality in the GWL trends in the selected cities. This confirms the high variability in the rainfall events and the delay (and sometimes even failure) in the monsoon rainfall, as reported in previous studies across India [74,75,82]. The results also strongly suggest that climate variability (rainfall) has affected the GWL trends in the selected cities.

Impact of Land Use/Land Cover Changes on Groundwater Level
LULC analysis (Figure 9) was performed using MODIS images collected between the years 2001 and 2018 (refer to Table 3 for detailed LULC area changes in the selected cities). Overall, the conversion of vegetation into grassland and croplands was majorly observed in all the selected cities. A relatively high growth in the urban area was seen in Raipur and Rajkot (1.6% and 2% increase, respectively), followed by Jodhpur (1.4% increase). However, in other cities, the urban area growth was less thañ 0.5%, with almost no growth in the urban area of Tiruchirappalli. The higher rate of population growth (12-45% based on the census), burgeoning real-estate activities, mainly to convert agricultural fields and uncultivated land into barren land (plots) for domestic housing and industrial expansion, clearing land for road construction and highways projects, increasing the urban built area, etc. [83,84], could be considered as the main factors for the LULC changes during the study period. Notably, the high growth of the urban area in Raipur might be attributed to the city being declared as the capital of the state Chhattisgarh after its partition from Madhya Pradesh in 2001 [85]. Exponential industrial and infrastructural growth in Raipur in the past two decades has led the city to transform into a major economic hub in the region [85][86][87]. The LULC change trends observed in Raipur during the study could be expected in other secondary developing cities in the future, where rapid migration, infrastructural and industrial growth, etc., could lead to major LULC changes (especially the conversion of the green cover into barren land) and the sudden change might create acute stress on groundwater systems [22,87]. Unlike the metropolitans, where the growth of the urban area is very high in the outskirts of the city, the seven cities in this study show no significant expansion of urban areas outside the pre-existing urban boundary. This implies that major economic activities in the secondary Indian cities primarily take place within the existing urban boundaries, as reported in previous studies [88,89]. This subsequently indicates that most of the infrastructural developments are progressing without significantly reducing green cover in these selected cities. Within the seven cities, the highest LULC changes ( Figure 9) were observed in Jodhpur, followed by Solapur and Rajkot. In Jodhpur, the cropping practice has changed, and the cropping area has increased significantly (increased by 8%) during the study period, which might be the reason for the high decreasing GWL in the region due to increased water requirement for irrigation activities [90].
On the other hand, in Rajkot, although a significant decline in the Grasslands area (green cover) occurred, the GWL trends in Rajkot are better than other cities. This trend, however, is recently observed mainly in the western arid states, which might be due to improved rainfall-GWL relationships. The positive (R = 0.14) and negative (R = −0.22) correlation coefficients were observed when cross-correlating decreasing and increasing GWL respectively, with the increase in urban areas. This strongly suggests that GWL decreases with an increase in the urban area (Table 4). These results are also in good agreement with other recent studies and suggest that LULC has a significant impact on GWL in developing cities [36,91]. Also, a very high correlation between the increase in the bare soil (decreasing green cover) area and decreasing GWL (R = 0.93) was observed. Thus, the major influence on decreasing GWL could be attributed to the conversion of the green cover into barren land for various anthropogenic activities such as infrastructural development, real-estate, etc. It may be inferred from these results that the irrigation and agricultural water requirements, although expected to have a very high impact on decreasing GWL in the selected developing cites (since most of the selected cities are still having major agricultural activities outside the urban boundary), are not the major driver; however, the conversion of the green cover (majorly the cropping area) into barren land has more influence on decreasing GWL. LULC change analyses also reveal that an increase in the urban area has a significant impact on decreasing GWL. However, the extent of the impact is comparatively lower (R = 0.14), which could be due to the proper functioning of the public water supply schemes and lower domestic water demands. However, as the urban area increases in the future, the stress on GWL might also increase in the selected cities, as previously observed in other already developed Indian cities [92,93].

Conclusions
This paper focused on understanding the GWL trends in the secondary developing cities in India. Seven cities were selected for the analysis based on various aspects, including the socio-economic and population growth trajectories. The key findings were derived based on statistical and LULC change analysis performed on the GWL data. Overall, it was found that the GWL varies between ±10 cm/year, with the trend to be significant (MMK test) in the majority of the wells. Non-stationarity was observed in the majority of the wells throughout the seasons, which suggests GWL has been affected by both climate variability (in monsoon rainfall) and LULC change. However, wells having a high decreasing GWL trend mostly did not exhibit non-stationarity in GWL changes. This indicates that the declining GWL of higher magnitudes has no substantial variation, possibly due to irreversible aquifer depletion caused by a high rate of pumping over a long time, where the effects of climate variability and LULC changes might be negligible. Less interdependence between seasonality (based on monsoon) and GWL trends was observed, which might have been caused by the impact of climate variability on groundwater recharge. The cause of non-stationarity in GWL was analyzed by linking climate variability (rainfall) and LULC with the GWL changes. Correlation analysis between GWL and rainfall indicates that the rainfall-groundwater level relationship has been seriously affected across all the selected cities, plausibly due to the failure of monsoon rainfall. Based on LULC change analysis, a high correlation between decreasing green cover (increasing barren land) (R = 0.93) and decreasing GWL indicates a high anthropogenic activity. Although agriculture-based activities are still actively practiced in secondary developing cities, its impact on GWL was surprisingly relatively low. The change in the urban area was majorly within the pre-existing urban boundary, and the domestic water requirement presently might not contribute much to the decreasing GWL. However, with more urbanization, this might change in the future. The major findings of this study can very well be adopted in other secondary developing cities with similar socio-economic growth patterns, especially in developing countries.
Instead of using the derived GWLs from satellite products, which are considered to contain more error and uncertainty, this paper used point scale well observations to assess the trends and non-stationarity in GWL changes. However, the bias in the trend and non-stationarity estimates are still possible due to the limitations in the data length and the number of available observation wells. Therefore, future research works, including the more observed data to conduct a similar study, are recommended.