Land subsidence caused by groundwater exploitation in Quetta and surrounding region, Pakistan

Land subsidence is effecting several metropolis in the developing as well as developed countries around the world such as Nagoya (Japan), Shanghai (China), Venice (Italy) and San Joaquin valley (United States). This phenomenon is attributed to natural as well as anthropogenic activities that include extensive groundwater withdrawals. Quetta which is facing similar subsidence phenomenon is the largest city of Balochistan province in Pakistan. This valley is mostly dry and ground water is the major source for domestic and agricultural consumption. The unplanned use of ground water resources has led to the deterioration of water quality and quantity in the Quetta valley. Water shortages in the region was further aggravated by the drought during (1998–2004) that affected the area forcing people to migrate from rural to urban areas. Refugees from the war torn neighboring Afghanistan also contributed to rapid increase in population of Quetta valley that has increased from 0.26 million in 1975 to 3.0 million in 2016. The objective of this study was to measure the land subsidence in Quetta valley and identify the effects of groundwater withdrawals on land subsidence. To achieve this goal, data from five Global Positioning System (GPS) stations in Quetta were acquired and processed. Furthermore the groundwater decline data from 41 observation wells during 2010 to 2015 were calculated and compared with the land deformation. The results of the GPS readings revealed that the land of Quetta valley is subsiding from 30 mm yr⁻¹ on the flanks to 120 mm yr⁻¹ in the central part. 1.5–5.0 m yr⁻¹ of groundwater level drop was recorded in the area where the rate of subsidence is highest. Whereas 9–10 cm of subsidence was recorded in the surrounding areas of Quetta where agriculture and settlements are high. The surrounding areas include Kuchlak, Mastung, Pishin, Gulistan and Hurumzai districts. These results were acquired using InSAR imagery collected from October 2014 to march 2019. So the extensive groundwater withdrawals in Quetta valley and surrounding areas is considered to be the driving force behind land subsidence.

1 Introduction

Quetta is the capital city of Balochistan located in the Quetta valley. The total watershed area of the valley is 1757 km² out of which 792 km² is covered mainly by alluvium (Fig. 1). Quetta valley falls in the Southern Basin watershed also known as Quetta Pishin Sub basin which is part of North East Pishin Lora Bain (NEPL).The hills of Mian Ghundi and Landi divide the valley into two basins namely Quetta valley to the north and Dasht plain to the south (WAPDA, 1988; Kazmi, 1973). Quetta valley is formed by three kinds of distinguishing landforms. The first landform is valley floor which is present in the central part of the valley. The second landform is Piedmont areas which is the zone between valley floor and mountains. Here the groundwater is mainly recharged due to high hydraulic conductivity (WAPDA, 1986). The third landform is the elevated mountains with steep inclinations that include Zarghoon, Takatu, Daghari, Chiltan and Murdar mountain ranges (Haque, 1986).

The aquifer system of Quetta valley is divided into four types. The uppermost aquifer is called the alluvial fan, this zone has high porosity and high permeability due to coarse grained sediments. The hydraulic conductivity is very high in this zone and groundwater is mainly recharged from this zone. Currently this zone is dry due to water depletion of the valley. The second aquifer of Quetta valley is called alluvium. This zone has intermediate porosity and permeability and is composed of gravels, sand and clay particles. The third zone is called the Boston formation which has very low hydraulic conductivity due to the presence of silt and clay particles. The forth zone of Quetta aquifer system is the Chilton Limestone. The primary porosity of this zone is low but due to high tectonic activity the fractures are developed making the secondary porosity very high. This zone is composed mainly of limestones. The depth of Chilton Limestone varies in the valley. In the northern end of Quetta valley the depth of this bed rock is about 1500 m and about 1058 m to the southern end (Alam and Ahmed, 2014). The unconsolidated alluvial aquifer and the bed rock are hydraulically connected to each other (Halcrow and Cameos Consultant Companies, 2008; TCI, 2004; Kazmi et al., 2005).

Figure 2 (a) shows the fissures in Kechi Baig area whereas (b) shows the damage that is caused to Dr. Sarfaraz hospital present on Gai khan Chowk.

The active Left-Lateral Chaman fault traverse on the western side of Quetta valley that splits the Eurasian and Indian plates making the geology of the area very complex. The stratigraphic sequence starts from Early Jurassic to Quaternary deposits in Quetta valley (UNDP, 1982; Kazmi, 1973; Hunting, 1961; Kazmi and Reza, 1970)

The population of Quetta valley has increased several folds in the past 30 years. Groundwater is the main source for domestic and agricultural consumption. The Quetta valley alluvium aquifer is under substantial stress due to the decreasing groundwater resources (WAPDA, 1989). Over exploitation of groundwater has created several problems in the valley that include significant groundwater decline in several parts of the valley, deterioration of water quality and subsidence of land (TCI, 2004; Khan et al., 2013, 2010). Land subsidence due to urbanization is reported in several metropolis around the world that include Shanghai in China (Yin el al., 2006), Tokyo in Japan (Hayashi el al., 2009) and Santa Clara valley in United States (Galloway et al., 1999).

Fissures were being developed in the central part of the Quetta valley and surrounding areas such as Kuchlak, Mastung, Pishin, Gulistan and Hurumzai areas. The unconsolidated material is usually thick in these areas. These fissures cover an area of about 1000 m long and 200 m wide. They were 3–5 m deep and 1–2 m wide at some places and have damaged several buildings in its path (Fig. 2). These fissures may be developed due to the compaction of sediments after the withdrawal of the groundwater. GPS data also suggest increased subsidence in the central part as compared to the flanks of the Quetta valley where the rate of subsidence is minimum, here the unconsolidated material is thin.

The objective of this study was to evaluate the groundwater decline in Quetta valley and surrounding areas, to identify land subsidence in Specific zones, to find out the causes of groundwater decline and to quantify the land cover over time using remote sensing in ArcGIS.

2 Material and methods

Groundwater decline data of the year (2010–2015) was obtained from Irrigation Department and Water and Sanitation Authority (WASA). This data was collected from 41 observation tube wells installed throughout the city. Land subsidence was measured using GPS units installed at five locations throughout the city. The GPS data were collected using Trimble NetRs 5700 and R7 receivers that download data after every 30 s. The accuracy of GPS unit was 1 mm in x, y and z axis. Its antenna is placed on a concrete flat floor which needs at least 4 satellite connections to work properly. GPS data was processed in GAMIT software that was referenced in ITRF2008 frame then linear trend was fit into the result to get the rate, this data was processed in Department of Earth and Atmospheric sciences, University of Houston. InSAR velocities were derived from C-band Sentinel-1 SAR data collected from October 2014 to March 2019. To ensure high phase coherence, the geometric and temporal baselines are limited to be less than 150 m and 90 days, respectively. InSAR data was processed in Department of geological sciences, university of Boulder, Colorado. Field survey was done to assess the damages caused by land subsidence. Water Decline data was processed in ArcGIS 10.2 and contour maps were made using (IDW) inverse distance weighted interpolation technique. Agricultural Remote sensing images were processed in ArcGIS 10.2 to assess the increase in vegetation over time. The Precipitation data was obtained from the surface water monitoring wing of the irrigation department.

Figure 3 (a) Graph of GPS data for QTIT station, 120 mm yr⁻¹ subsidence recorded on this site. (b) Graph of GPS data for QTAG station, 106 mm yr⁻¹ subsidence recorded on this site. (c) Graph of GPS data for JHAK station, 88.3 mm yr⁻¹ subsidence recorded on this site. (d) Graph of GPS data for SBKW station, 29 mm yr⁻¹ subsidence recorded on this site. (e) Graph of GPS data for KHAL station, 35.1 mm yr⁻¹ subsidence recorded on this site.

Figure 4 Map showing subsidence in the Quetta valley. The dark brown color is showing subsidence of 100–120 mm yr⁻¹, the light brown color is showing subsidence of 80–90 mm yr⁻¹, the yellow color is showing subsidence of 30–40 mm yr⁻¹, the green color is showing subsidence of 20–30 mm yr⁻¹. The fissures can also be seen on the map whereas the Kotwal Weather station is shown by a black arrow. 30 m Digital Elevation Model (DEM) of Shuttle Radar Topographic Mission (SRTM) were downloaded from ©NASA's Earth Explorer website to identify the topography and profile of the study area.

3 Land subsidence

Land Subsidence is the sinking or lowering of land which is mostly caused by fluid withdrawal. The underground microscopic pore spaces are filled with water molecules which makes up the hydrostatic pressure in a pore. When the water molecules are extracted through a tube well, a vacuum is formed and the hydrostatic pressure diminishes in the void as a result the pore space collapses and the whole ground become subsided.

Land subsidence in Quetta valley was first reported by (Khan et al., 2013) when a drop of about 10 cm yr⁻¹ was recorded by GPS data in the valley. The network of these GPS points were extended to five locations in the valley to get more accurate results. Two of these sites were located in University of Balochistan that were codenamed as QTIT and QTAG. The other point was located at Arbab Karam Road that is codenamed as JHAK. The next site was located at Sardar Bahadur Khan Woman University that is codenamed as SBKW. The last site is located on a gas station that is adjacent to the Frontier Corps Camp and is codenamed as KHAL. QTIT is a permanent station that continuously download data from the GPS device whereas the rest of the stations were campaign sites where GPS device was installed for one week during the field.

Figure 5 InSAR Map© showing subsidence in Blue which is 8–10 cm yr⁻¹, red color (saturated) is showing uplifting due to tectonic activity. Black dash line is showing Chaman transform plate boundary, pink dots shows recent earthquakes. InSAR velocities were derived from C-band Sentinel-1 SAR data collected from October 2014 to March 2019.

Figure 6 (a) Graph of water decline on eastern flank of Quetta valley. The average decline recorded in this zone is 1.5–5.0 m yr⁻¹. (b) Graph of water decline on western flank of Quetta valley. The average decline recorded in this zone is 0.5–1.5 m yr⁻¹. The observation well in Chilton housing Scheme shows rise in water table.

4 Results

4.1 QTIT

The first reading on QTIT station was recorded on 25 August 2008 whereas, the last reading was taken on 22 January 2016. Its Latitude is 30.16291436 and Longitude is 66.9888975. The total subsidence recorded during this duration was 890 mm (≈ 3 ft) whereas the annual subsidence recorded on this station was 120.2 mm yr⁻¹ (Fig. 3a). This is the central part of the valley where the unconsolidated material is thick.

4.2 QTAG

The first reading on QTAG station was recorded on 17 July 2006 whereas the last reading was taken on 24 January 2016. Its Latitude is 30.16624095 and Longitude is 66.9908418. The total subsidence recorded during this time period was 1036 mm (3.3 ft) whereas the annual subsidence recorded on this station was 106.8 mm yr⁻¹ (Fig. 3b).

4.3 JHAK

The first reading on JHAK station was recorded on 6 June 2011 whereas the last reading was taken on 17 January 2016. Its Latitude is 30.18506981 and Longitude is 66.99151868. The total subsidence recorded during this duration was 430 mm. whereas the annual subsidence recorded on this station was 88.3 mm yr⁻¹ (Fig. 3c).

4.4 SBKW

The first reading on SBKW station was recorded on 9 January 2012 whereas the last reading was taken on 12 January 2016. Its Latitude is 30.19146625 and Longitude is 66.9591966. The total subsidence recorded during this duration was 116 mm. whereas the annual subsidence recorded on this station was 29 mm yr⁻¹ (Fig. 3d). This is the area located on the western flank of the valley where unconsolidated material is thin.

4.5 KHAL

The first reading on KHAL station was recorded on 10 April 2009 whereas the last reading was taken on 9 January 2016. Its Latitude is 30.30112026 and Longitude is 66.93319597. The total subsidence recorded during this duration was 218 mm. whereas the annual subsidence recorded on this station was 35.1 mm yr⁻¹ (Fig. 3e).

Figure 7 Map of water decline in Quetta valley. The red and dark red color show the water decline of 1.5–15 m yr⁻¹ whereas, the green and dark green color show the water decline of 0–1.5 m yr⁻¹. 30 m Digital Elevation Model (DEM) of Shuttle Radar Topographic Mission (SRTM) were downloaded from ©NASA's Earth Explorer website to identify the topography and profile of the study area.

The annual subsidence rate determined from the GPS data is shown in Fig. 4. The dark brown color shows the highest rate of subsidence i.e., 100–120 mm yr⁻¹, the light brown color shows the subsidence of 80–90 mm yr⁻¹, the yellow color shows subsidence of 30–40 mm yr⁻¹ and the green color displays subsidence of 20–30 mm yr⁻¹. Location of the fissures can also be seen on this figure.

InSAR images collected from October 2015 to March 2019 revealed that the subsidence phenomenon is not confined to Quetta district only, in fact the surrounding districts where agriculture is utilizing significant ground water are also witnessing significant land subsidence due to ground water decline. The dark blue color in Fig. 5 shows 8–10 cm yr⁻¹ of subsidence in these area. These include Kuchlak, Mastung, Pishin, Gulistan and Hurumzai areas.

Figure 8 Graph of population in Quetta valley which increased from 0.26 million in 1975 to 3.0 million in 2016.

Figure 9 ©Landsat Thematic Mapper bands 7, 4 and 2 displayed as red, green and blue showing the vegetation increase overtime in Quetta valley. (a) is from 21 June 1987 whereas (b) is from 16 August 2011. Increase in vegetation can be seen in bright green colors in the northern and southern end of Quetta valley.

Figure 10 Location of drilled well sites at Sheikhmanda area. The green dots are the bore holes drilled in area of 1.3 km². Base map taken from ©Google earth.

5 Groundwater depletion

The pioneer artesian well in Quetta valley was dug in 1889 at railway station using steam engine (Oldham, 1892). The depth of the well was 35 m while the hourly discharge was 20 000 gallons (Oldham, 1892). But over time water started to decline in the Quetta valley and WAPDA started monitoring of groundwater in Quetta valley in 1988. Around 177 observation tube wells were installed in the Pishin Lora Basin to observe the groundwater conditions. The depth of those wells were ranged from 3.5 m (11 ft) to 150 m (500 ft). Due to overexploitation of groundwater, around 123 observation wells were completely dried up and now 54 observation wells are operational. A drop of 0.25 m yr⁻¹ was first observed by WAPDA in 1989 (WAPDA, 2001). Ahmad (2007) recorded decline of 1.09 m in the 1990s. In the 2000s Water and Sanitation Authority (WASA) installed several new observation wells in Quetta valley to monitor the groundwater decline. The depletion of Quetta valley water table was calculated from the year 1987 to 2010 from 55 observation wells. The average decline during this period was 1.5 m yr⁻¹ in several parts of the valley (Khan et al. 2013).

To get the latest groundwater level, data from 41 observation wells were collected from WASA and Irrigation Department for the year 2010–2015. The average decline from 2010 to 2015 were calculated and surface was generated using Inverse Distance Weighting (IDW) interpolation technique in ArcGIS. It was observed that the average decline on the eastern side of the valley is high (1.5–5.0 m) (Fig. 6a) where the population is congested while on the western side the depletion is low (0.5–1.5 m) where the population is scarce (Fig. 6b).

Figure 11 Graph showing the precipitation data, precipitation gauge is located in Killi Kotwal Station in Quetta. Drought can be seen from the year 1998–2004 whereas heavy precipitation in 2011.

In Fig. 7 the dark red color shows the average decline of about 5–15 m on the location of Ali Bahadur and Al-Gilani Road. The light red color represents the annual decline of 1.5–5.0 m on the location of Killi Shabak road and Killi Nasaran Malezai Kuch Road in the north of the valley to Pashtoonabad and office of the Geological survey of Pakistan to its south. This zone covers the whole eastern part of Quetta valley. The yellow color is showing the annual decline of 0.5–1.5 m on the location of Chasma Achozai to the north and MNW-QA-10 to the south of the valley. The green color represent the zone where the groundwater level is rising, this zone is located on Killi Kirani and Samungali Housing Scheme.

There was no observation well installed nor the water table data was available outside the Quetta district where field observations took place.

6 Causes of water depletion and land subsidence

6.1 Population

According to the census of 1975, the population of Quetta valley was just 260 000 and it steadily increased and reached 430 000 in 1981. Due to urbanization and migration of refugees from Afghanistan's civil war, the population increased to 1.1 million in 2000. Due to drought in the last decade a lot of people were forced to migrate to Quetta to earn their livelihood. According to recent estimates the population of Quetta valley is around 3.0 million in 2016 (Fig. 8). This huge population is a tremendous burden on the aquifer of Quetta valley and it is creating the problems of supply and demand in the valley due to which the aquifer is depleting at an alarming pace.

6.2 Agriculture

Different crops are grown in Quetta valley and surround areas that include wheat, apricot and apples. The agriculture has increased several folds as more people earn their livelihood from it. Groundwater is the only source used to grow these crops. Due to this phenomena the groundwater is depleting at a faster rate. Two satellite images were processed in ArcGIS that were downloaded from USGS Global Visualization Viewer (GLOVIS) to see the land use overtime. These images were Landsat Thematic mapper at 30 m resolution from the year 1987 and 2011. Dramatic increase of vegetation can be seen in northern and southern end of Quetta valley (Fig. 9). The bands that were used to enhance the vegetation were 7, 4 and 2 displayed as red, green and blue. Similar situation was field verified in Kuchlak, Mastung, Pishin, Gulistan and Hurumzai areas.

6.3 Illegal drillings

Illegal drillers are drilling tube wells to extract groundwater resources without any authorization from the government. The data of illegal drilling sites is not available. To demonstrate the congestion of wells in the Quetta valley, an area of 1.3 km² was observed (Fig. 10). Out of the 20 wells only one of the well was government owned. The illegal tube wells are drafting groundwater at unprecedented rate due to no check and balance.

6.4 Precipitation

The precipitation data of the year 1976–2015 was obtained from the Irrigation Department and plotted in Fig. 11. Several fluctuation can be seen with the substantial drought which occurred during 1998–2004. During this drought the reservoirs of Spin Karez and Hana Urak which are situated to the eastern side of Quetta valley were completely dried up which resulted in significant water decline in Quetta valley aquifer.

7 Discussion

GPS results of Quetta valley illustrates that the subsidence is intense in the central part of the valley where the thickness of unconsolidated material is high. The average subsidence of QTAG and QTIT station has slightly increased over time. Khan et al. (2013) observed the average subsidence at QTAG station (2006–2008) to be 81 mm yr⁻¹ and the average subsidence at QTIT station (2008–2009) to be 116 mm yr⁻¹ – whereas the average subsidence at QTAG station (2006–2016) has increase to 106 mm yr⁻¹ and the average subsidence at QTIT station (2008–2016) has increased to 120 mm yr⁻¹. The subsidence at the flanks of the Quetta valley is far less than the central area e.g. the subsidence at SBKW is 27 mm yr⁻¹ and KHAL is 35 mm yr⁻¹. In this zone the thickness of unconsolidated material is far less than the central part.

The InSAR imageries collected from 2014 to 2019 shows that the subsidence due to groundwater decline is not limited to Quetta district alone. Other districts in Balochistan province such as Kuchlak, Mastung, Pishin, Gulistan and Hurumzai districts are also facing about 8–10 cm of subsidence per year.

Increase in water depletion was recorded during the study. 1–1.5 m yr⁻¹ water decline was recorded at several places in Quetta valley during 1987–2010 (Khan et al., 2013; Ahmad, 2007) whereas the water depletion during 2010–2015 had reached 1.5–5.0 m yr⁻¹.

8 Conclusion

GPS data show that the central part of the Quetta valley has subsided more than 1 m during the last decade and continue to subside at about 120 mm yr⁻¹. Whereas the subsidence at the flanks of the valley is about 30–40 mm yr⁻¹. Due to land subsidence, fissures are being developed in several districts in Balochistan that has damaged numerous buildings and houses. The groundwater is also depleting at an alarming rate that is 1.5–5.0 m in several parts of the city. The increase in agriculture and population growth are the main causes of water depletion and subsidence in Quetta valley. The land subsidence and its associated fissures will cause billions of rupees of damages to the buildings in the future if the trend continues. The process may also mix fresh and sewerage water which will create health hazards in the city.