Use of InSAR Data for the Measurement of Land Subsidence in Response to Groundwater Fluctuations and Climate Change (Case Study: Ardabil Plain, Iran)

doi:10.21203/rs.3.rs-1410400/v1

The Ardabil plain, with an approximate area of 1097.2 km² in northwestern Iran, has experienced land subsidence due to intensive groundwater withdrawal and long seasons of drought in recent years. Different techniques have been used to investigate and evaluate subsidence in this region including: Global Positioning Systems (GPS), Levelling, and Geotechnical methods. These methods are typically expensive, time-consuming, and identify only a small fraction of the areas prone to subsidence. This study employs an Interferometric Synthetic Aperture Radar (InSAR) technique to measure the long-term subsidence of the plain. An open-source SAR interferometry time series analysis package, LiCSBAS, that integrates with the automated Sentinel-1 InSAR processor (COMET-LiCSAR) is used to analyze Sentinel-1 satellite images from October 2014 to January 2021. Processing of Sentinel-1 images shows that the Ardabil plain has been facing rapid subsidence due to groundwater pumping and reduced rainfall, especially between May 2018 to January 2019. The maximum subsidence rate was 45 mm/yr, measured at the southeastern part of the plain. While providing significant advantages (less processing time and disk space) over other InSAR processing packages, implementation of the LiCSBAS processing package and its accuracy for land subsidence measurements at different scales needs further evaluation. This study provides a procedure for evaluating its efficiency and accuracy for land subsidence measurements by comparing its measurements with the results of the GMTSAR and geotechnical numerical modeling. The results of geotechnical numerical modeling showed land subsidence with an average annual rate of 3.8 cm between 2006 and 2020, which was close to measurements using the InSAR technique. Comparison of the subsidence measurements of the Ardabil plain using the LiCSBAS package with results obtained from other techniques shows that LiCSBAS is able to accurately detect land deformation at large scales (~ km). However, they may not be optimized for more local deformations such as infrastructure monitoring.

LiCSAR

GMTSAR

SBAS-InSAR

Coherence

Geotechnical

Drought Index

Land subsidence caused by compaction of over-drafted aquifer systems and climate change has caused severe damage to agricultural lands, industrial activities, and urban infrastructure [1]. Land subsidence has been documented in several areas in the world including: Tehran, Iran [2], Los Angeles [3], Xi'an, China [4], Calcutta [5] and Las Vegas Valley and vicinity, Clark County, Nevada [6]. This phenomenon is a frequent environmental problem that needs to be examined in long-term, systematic research [7, 8]. Therefore, a permanent monitoring system with a correct estimate of surface deformation is essential to show the existing threats posed by the land subsidence.

Over the past two decades, different techniques have been proposed and implemented to detect and evaluate land subsidence including Leveling [9, 10], Global Positioning Systems (GPS) [11, 12], and geotechnical techniques [13, 14, 15]. These methods are typically expensive, time-consuming, and provide continuous monitoring of land deformation over only a small fraction of the areas prone to subsidence. Compared with traditional measuring alternatives, the InSAR technique can detect a wider spatial range of surface deformation due to its broad spatial coverage and high accuracy [16, 17]. The results obtained from this technique are capable of measuring small to large-scale displacements on the land surface. The most common families of time-series InSAR techniques used for land subsidence and land surface deformation measurements are multi-temporal interferometry methods including Persistent Scatterer – PS-InSAR (PSI) [18, 19] and Small Baseline Subset (SBAS) [20, 21]. The PSI technique generates differential interferograms with one common master: identifying persistent point-wise reflectors, such as urban areas, with a high PS density [22]. The SBAS-InSAR technique, on the other hand, relies on a redundant network of image pairs with a short spatial and moderate temporal baseline. It detects the temporal evolution of surface deformations and increases the spatial coverage, especially over non-urban areas where the PS density may be low [21, 22]. Compared to PSI, SBAS requires smaller Synthetic Aperture Radar (SAR) images and has a stronger ability to obtain nonlinear deformation information, particularly in regions with higher signal decorrelation [23]. Due to the advantages of SBAS over PSI, particularly in agricultural regions, the SBAS-InSAR technique was considered for this study.

The SBAS-InSAR technique with Sentinel-1 data has been widely used to estimate surface dsiplacement [7, 24, 25, 26, 27, 28, 29, 30, 31, 32]. For example, Zhou et al. (2017) analyzed the InSAR time-series of land subsidence in the Eastern Beijing plain, China. Different types of land use (e.g., expressway, subway, and urban residential land; railway, agricultural, and urban residential land; and agricultural, villages and towns residential land) and water resource utilization were considered in this study. Their analysis showed a maximum land subsidence rate of 146 mm/yr, mainly in inhabited lands during the years of 2011 to 2015 [7]. Additionally, a comparison of the land subsidence measurements and changes in groundwater level indicated that the subsidence is directly related to changes in groundwater level. This study also compared the SBAS-InSAR land subsidence results with ground leveling measurements, indicating a precision of 2 mm for the SBAS-InSAR technique. Similar conclusions were made by [27, 33] when using Sentinel-1 SAR and SBAS for land subsidence monitoring in Hanoi, Vietnam, and Gorgan Plain, Iran, respectively. The results of [27] also indicated that InSAR-derived land deformation data on a monthly scale could be used for characterization of the hydraulic head and aquifer properties in confined aquifer systems, even those experiencing significant subsidence. Although data availability has been greatly improved by Sentinel-1, there are still barriers to making full use of this almost universal coverage and the vast amount of data involved. Also, the use of open-source software tools to process Sentinel-1A (S1A) images such as Generic Mapping Tools SAR (GMTSAR: an open-source - GNU General Public License - InSAR processing system designed with Generic Mapping Tools - GMT) [34], Sentinel Application Platform (SNAP) (http://step.esa.int/main), and time series analysis by Stanford Method for PS (StaMPS) [35] is challenging. These concerns are particularly salient in cases where long time series and large stacks of Sentinel-1 scenes are to be analyzed, often requiring more disk space for processing [36, 37].

Recently, an open-source SAR interferometry time series analysis package, LiCSBAS, has been developed and integrated with the automated Sentinel-1 InSAR processor (COMET-LiCSAR) [38]. Because LiCSBAS utilizes products published directly on the COMET-LiCS web portal (https://comet.nerc.ac.uk/COMET-LiCS-portal/), it does not need to pre-process and produce interferograms, thus saving significant time series processing time and disk space. In the LiCSBAS processing scheme, interferograms with many unwrapping errors are automatically identified by loop closure. Also, reliable time series and velocities are derived with the aid of masking using several noise indices. Furthermore, relative displacements on a large scale (> 100 km) and locally (2 km) can be measured with an accuracy of < 1 cm/epoch and ~ 2 mm/yr in velocity [36].

This study aims to measure the long-term subsidence of the Ardabil plain from 2014 to 2021 using Sentinel-1 SAR data processed in the LiCSBAS processing scheme along with a geotechnical based numerical modeling technique in response to groundwater fluctuations and seasonal weather change. Over the last decade, the Ardabil plain, one of the most important plains agriculturally and industrially in northwest Iran, has experienced significant subsidence due to groundwater level fluctuations and climate change [39]. In the past, InSAR data for monitoring subsidence in the Ardabil plain was very limited. With the advent of Sentinel-1 satellite in 2014, abundant and useful SAR data have been provided for the Ardabil plain, which have the potential to detect ground deformation in high spatial and temporal resolution. In this research, the LiCSBAS package was used to estimate spatially and temporally detailed deformation time series of the Ardabil plain from the LiCSAR products for the period of October 2014 to January 2021. Although LiCSBAS accurately detects large-scale deformation (~ km), it has not been fully tested in local contexts. Specially, LiCSAR products are not optimized for local deformations such as infrastructure monitoring. This is mainly because of its generated interferograms which are multilooked by factors of 5 in the range and 20 in the azimuth directions, geocoded onto a 100 m grid, and spatially filtered for the main purpose of large-scale deformation monitoring and tectonic strain mapping [36]. To evaluate consistency of the LiCSBAS package in measuring land subsidence at scales relevant to Ardabil plain (~ 20 km), the results of the LiCSBAS were compared with the results of GMTSAR. GMTSAR is a well-tested tool for the analysis of satellite images and has been used for different applications such as land subsidence estimation [27, 40, 41, 42], landslide estimation [43, 44], and earthquake and other deformation estimation [45, 46]. The motivation for choosing GMTSAR in this paper was due to its high processing performance, C-shell scripting susceptibility, good potential for parallel implementation, and free access compared to other commercial SBAS processing systems (such as SARscape - https://l3harrisgeospatial.com and GAMMA [47]). Also, the validity of SBAS-InSAR technique results was discussed using the time-series of land deformation obtained from permanent GPS observations, groundwater changes, and a sophisticated geotechnical based numerical modeling technique. The findings of this research provide researchers and practitioners with an accurate technique for measuring land subsidence over large areas due to the spatial coverage of available InSAR data and LiCSBAS’s efficient processing strategy. This is a critical need for groundwater management in areas that are prone to subsidence and provides a way for efficacy assessment of different strategies to limit land subsidence.

Description of Study area

The Ardabil plain, with an approximate area of 1097.23 km², is an inter-mountain plain located in northwestern Iran with geographical coordinates of 38°5' to 38°27' N and 48°9' to 48°37' E (Fig. 1).

Figure 1. The geographical location of Ardabil plain.

The Ardabil plain has been formed out of Quaternary alluvial deposits due to weathering and erosion of surrounding mountains [48]. Conglomerate, volcanic ashes, and lahars are the predominant lithological outcrops in the west of the plain. There are abundant springs in these units which infiltrate into the groundwater and affect aquifer recharge. In the past, most springs, qanats, and rivers were used for water supply, but now water supply is mostly based on the exploitation of wells. Based on the results of geophysical studies, pumping test data, and drilling logs, Ardebil aquifer has a maximum thickness of 220 m and is mainly composed of gravel, sand, and a small amount of clay. The transmissivity of the aquifer varies between 50 and 2200 m²/day, and the specific yield ranges from 0.021 to 0.14. The general groundwater flows from other directions into the north–west of the plain [48]. More than 85% of exploited groundwater is used for agricultural activities, 14% is used for drinking purposes and 1% for industrial use. The yearly maximum and minimum groundwater levels are usually recorded in May and September, respectively [49]. Groundwater levels in the plain are directly affected by human activities and infrastructure, including industrial and agricultural activities [50]. The climate of the study area is cold and semi-arid. Figure 2 shows precipitation and temperature times series of the region for the priod of study. This data was obtained from Tropical Rainfall Measurement Mission (TRMM) satellite precipitation data (https://disc.gsfc.nasa.gov/datasets/TRMM_3B43_7) using the TRMM adjusted merged microwave infrared precipitation rate (in mm/hr) and Root-Mean-Square (RMS) precipitation-error estimates, and Famine Early Warning Systems Network Land Data Assimilation System (FLDAS) datasets (https://das.gsfc.nasa.gov/fldas). Times series of precipitation and temperature show average annual precipitation and temperature of ~ 0.13 mm/hr (mostly happening during the rainy season - August to April), and 11.1˚C (with the highest temperature in August and the lowest in February), respectively.

Figure 2. Total precipitation from (2000–2020) based on the TRMM rainfall estimate month data, and annual temperature rate from (2002–2020) based on the FLDAS datasets (https://das.gsfc.nasa.gov/fldas).

Research methodology

InSAR data

In this study, Sentinel-1 satellite images from October 2014 to January 2021 in a time interval of 12 days were used (Fig. 3). The Interferometric Wideswath (IW) mode was used, which creates images with a 250 km swath at 5×20 m spatial resolution. Details of the SAR images that were used for the InSAR time series analysis are presented in Table 1. We processed the LiCSAR generated interferograms over two LiCSAR frames covering Ardabil plain in ascending and descending orbits (LiCSAR frame ID 101A_05193_131313 and 006D_05111_131313) using LiCSBAS package. The interferograms were created with the GAMMA SAR software and filtered using an adaptive phase filter [38]. GeoTIFF files of unwrapped phase and coherence are published on the COMET-LiCS web portal on a consistent geographic frame basis in 250 km². Both ascending and descending data were available for the whole or a specific period. The numbers of processed epochs for the ascending and descending frames were 137 and 136, respectively.

Figure 3. Geographical locations of SAR images used to study land subsidence in Ardabil plain (https://comet.nerc.ac.uk/COMET-LiCS-portal/).

Table 1

Details of the SAR images used for the InSAR time-series analysis.
SAR Sensor	Orbit	(λ)	Revisit cycle	(θ°)	N	Temporal coverag	Mode
S1A-101	Ascending	C(5.56)	12 days	39.46	137	20141129–20200910	IW2
S1A-6	Descending	C(5.56)	12 days	39.50	136	20141006–20210114	IW1
θ°: Incidence angle of the sensor N Number of processed epochs

LiCSBAS time-series analysis

LiCSBAS is an open-source package in Python3 to carry out InSAR time series analysis using LiCSAR products (i.e. unwrapped interferograms and coherence) that are freely available on the COMET-LiCS web portal [36]. LiCSBAS can run all steps selected for analysis with specified parameters at once using a batch script. The first step is to run the script, download LiCSAR products and then convert file format. In the next step, to reduce the impact of unwrapping errors in Ardabil plain, before the subsequent time-series analysis, pixels with a coherence of 0.1 are masked and then the data are unwrapped automatically using SNAPHU [51]. Thereafter, all the unwrapped data and corresponding coherence images are clipped to the Ardabil plain area. The unwrapped coverage and coherence threshold are set to < 0.5 and 0.06, respectively. The unwrapped interferograms are resampled and geocoded using the Shuttle Radar Topography Mission Digital Elevation Model (SRTM DEM 30 m). The COMET-LiCSAR system processes the InSAR data and generates interferograms that connect each epoch to three or four nearest acquisitions in time, backward and forward. LiCSBAS employs a modified small-baseline NSBAS approach [52, 53] for the time-series analysis. In the time-series analysis, deficient interferograms are identified and then discarded based on the coherence and coverage of the unwrapped data by checking loop closure [36, 54].

The whole interferometric networks are then inverted for incremental displacements between the acquisition dates, using the least-squares method [36], to estimate the displacement value of each pixel in the Ardabil plain. Overall, the baseline networks obtained from S1A images for the ascending (A-101) and descending (D-6) tracks include 453 and 365 interferograms, respectively (Figs. 4a and 4b). Interferograms that contain noisy data or do not pass a phase loop closure test (indicating severe unwrapping errors) are automatically excluded from the analysis [55, 56]. Finally, the maps are high-pass filtered in time and low-pass filtered in space using a Gaussian filter kernel [57], to separate noise components from the displacement time series (set to 0.50 year − 182 days).

Figure 4. Perpendicular baseline configuration and network of the 453 and 365 interferograms formed from track (a) A-101 and (b) D-6 S1A images used in the study area.

In order to compare the LiCSBAS and permanent GPS data, the following equation was used to obtain GPSLOS measurements in LOS direction [58]:

\({GPS}_{LOS}={GPS}_{up}cos\theta -\left({GPS}_{e}cos\alpha -{GPS}_{n}sin\alpha \right)sin\theta\)

(1)

where \(\alpha\) is the azimuth of the LOS vector, θ is the incidence angle, and GPS_up, GPS_e, GPS_n are the velocities in vertical directions, the eastwest and northsouth, respectively.

We estimated the surface displacements in the Ardabil plain using LiCSBAS. Then, the accuracy and consistency of the results were analyzed in both LiCSBAS and GMTSAR InSAR tools.

Interpretation of the S1A data by LiCSBAS velocities

The mean velocity maps derived from the time series analysis of S1A images by LiCSBAS are shown in Fig. 5. Several surface displacement signals including subsidence and uplift were detected in the Ardabil plain. Despite having almost all pixels masked in highly vegetated areas due to decorrelation, more than 110 thousand valid pixels (60%), mainly distributed on the Ardabil plain, still remain (Fig. 5). To obtain meaningful and accurate displacement measurements, a suitable reference area, where the displacement signal is almost zero and marked by the black circle in Figs. 5(a) and 5(b), was considered for the study area. The calculated mean LOS velocities of S1A images from A-101 (measured in the period of November 2014 to September 2020) and D-6 (measured in the period of October 2014 to January 2021) tracks are between − 45 mm/yr to + 14 mm/yr and between − 53 mm/yr to + 34 mm/yr, respectively. Since the displacement is illustrated in LOS direction, positive and negative velocities represent the motion of the ground towards and away from the satellite, respectively. The mean LOS velocity values show a maximum displacement rate of 45 mm/yr in the LOS direction in an area located at the southeast part of the plain. This area is located in a cropland area that suffers from excessive extraction of groundwater. Results presented in Fig. 5 also show signs of land subsidence with maximum displacement rate of 20 mm/yr in the cropland areas in the northwest region (Samarin). Other areas of the plain show no significant subsidence.

Figure 5. Mean LOS velocity of the results S1A images by LiCSBAS and displacement relative to the reference point. (a) From October 2014 to January 2021 for D-6 track, and (b) From November 2014 to September 2020 for A-101 track.

Figure 6 shows displacement time-series in three different sites, near the cities of Ardabil, Khalilabad, and Araloyebozorg. Khalilabad and Araloyebozorg are both located in the southeast with the highest rates of subsidence, and obtained using LiCSBAS. Results indicated a decreasing trend of displacement with time over the years of 2014 to 2021, with a maximum LOS displacement rate of 45 mm/yr near Araloyebozorg. Also, results indicate an increase in the subsidence rate in 2018. Different factors are responsible for this high rate of displacement in the southeast.

Figure 6. The displacement time-series from S1A images at three representative locations. (a) From October 2014 to January 2021 for D-6 track and, (b) From November 2014 to September 2020 for A-101 track.

Figure 7 shows some of the factors that play the most significant roles in land subsidence, including groundwater table, alluvial deposit thickness, type of soil layers and stratigraphy. According to Fig. 7, the southeast part of the plain has the maximum alluvial thickness of 190 m, especially in areas near the cities of Araloyebozorg and Khalilabad. According to studies by [48], this area has a clay layer with a thickness of more than 120 meters, starting from a depth of about 25 meters. The thickness of the sand layer is only 25 meters in this area. Also, the ratio of fine-grained to coarse-grained materials in this area is higher compared to the other parts of the plain, increasing potential for land subsidence [59]. The groundwater table measurements from the observation wells indicate an average depth of 60 to 70 meters in the southeast, which are deepest in the whole plain. The central part of the plain has the highest illegal and legal well density with significant water withdrawal, and therefore, it is expected to observe higher rates of subsidence in this area compared to other parts of the plain. However, since the groundwater level in this area is high (due to the bowl-shaped nature of the plain), the plain in this area experiences lower rates of subsidence with groundwater withdrawal.

Figure 7. Parameters involved in the issue of subsidence: groundwater table, alluvial deposits thickness, type of soil layers and stratigraphy.

Consistency assessment

To evaluate consistency of the LiCSBAS package in measuring land subsidence of the Ardabil plain, results of the package were compared with subsidence measurements obtained using the GMTSAR [60]. For this purpose, S1A images of 30 (A-101) and 28 (D-6) tracks, spanning March 2017 to January 2019, were analyzed using the GMTSAR. The interferometric process by GMTSAR based on SBAS-InSAR followed three main steps: First, a preprocessor was used for each satellite data type to convert the native format and orbital information into a generic format. Second, the InSAR processor was used to focus and align stacks of images and convert map topography into phase which is then used in radar coordinates, and form the complex interferogram. The topographic phase component was corrected using the SRTM 30 m DEM. Third, a postprocessor was applied based on GMT, to filter the interferograms and construct interferometric products of phase, coherence, phase gradient, and LOS displacement in both radar and geographic coordinates [60]. GMTSAR uses an adaptive power spectrum filter and a Gaussian filter for the unwrapping process (set to default filter parameters of alpha = 0.1 and patch size = 32). In general, the numbers of interferograms associated with a perpendicular baseline of smaller than 120 m and a temporal baseline of smaller than 130 days obtained from S1A images of (A-101) and (D-6) tracks are 101 and 97, respectively. Using the least-squares method, the displacement value of each pixel was estimated in the Ardabil plain. To make the results comparable with the LiCSBAS/LiCSAR processing chain, the same multilooking factors of LiCSAR data were employed for GMTSAR processing and the coherence threshold was set to 0.06.

Figure 8 compares the LiCSBAS velocities with those measured using GMTSAR for the period of March 2017 to January 2019 for A-101 track. Figures 8(a) shows the maximum subsidence rates measured using the GMTSAR and Fig. 8(b) shows the results obtained using LiCSBAS. As observed in this figure, in terms of the maximum subsidence rates, both packages show close values of -53 mm/yr (GMTSAR) and − 52 mm/yr (LiCSBAS). For the D-6 track, these values changed to -46 mm/yr (using GMTSAR) and − 49 mm/yr (using LiCSBAS), as observed in Figs. 8(c) and 8(d). Both GMTSAR and LiCSBAS results showed a major land subsidence concentrated in the southeast of the Ardabil plain, near the city of Araloyebozorg and Khalilabad. However, the difference between LiCSBAS and GMTSAR results is the number of masked pixels in areas of low coherence which is significantly larger in GMTSAR.

Figure 8. Comparison of the LOS displacement maps from S1A images (March 2017 to January 2019) for A-101 track using (a) GMTSAR, (b) LiCSBAS, and for D-6 track using (c) GMTSAR and (d) LiCSBAS.

In Fig. 9, the average coherence maps (March 2017 to January 2019) in Ardabil plain using GMTSAR and COMET-LiCSAR system are presented. Based on the results presented in Figs. 9(a) and 9(b), the values of average coherence for GMTSAR and COMET-LiCSAR are 0.27 and 0.37, respectively. Both results showed high coherence values (higher than 0.7) in areas around the city and part of the airport. Low coherence values (less than 0.3) were observed in areas such as croplands on the outskirts of Ardabil plain, vegetated sides of the Talesh heights in the eastern side of Ardabil plain, Shourabil Lake located in the south of Ardabil city, steep slope, rangelands, and vegetated slopes of Sabalan mountain in the western part of Ardabil plain.

Figure 9. Comparison of the average coherences for the study area: (a) GMTSAR and (b) COMET-LiCSAR system.

Figure 10 shows scatter plots (Figs. 10a) and histograms (Figs. 10b) of average coherence obtained using GMTSAR and COMET-LiCSAR system. As presented in Fig. 10(a) for the coherences of 0.15 to 0.7, the pixels are distributed along the red line, showing consistency between the GMTSAR and COMET-LiCSAR systems. This figure clearly shows a good correlation with an R² value of 0.84. Based on the results presented in Fig. 10(b), the mean coherence difference between GMTSAR and COMET-LiCSAR is 0.1, with a standard deviation (Std) of 0.01. It is observed that the coherence values produced by GMTSAR are slightly lower than the COMET-LiCSAR system.

Figure 10. Comparison of the average coherence in GMTSAR and COMET-LiCSAR system using (a) Scatter plots and (b) Histograms.

Also, for further comparison of differences in coherence values of the two systems (GMTSAR and COMET-LiCSAR), a pair of unwrapped interferogram phases (20170423–20170610) under challenging conditions were analyzed (Fig. 11). As presented in Figs. 11(a) and 11(b), despite its substantially higher coherence values, the COMET-LiCSAR unwrapped interferogram phases for the scene pair are not different from GMTSAR results for most parts of the Ardabil plain (R² = 0.79 and Std = 0.06), and the interferometric phase results have similar patterns. This behavior may be the result of using different approaches to calculate the unwrapped phase using the applied Gaussian filter.

Figure 11. Comparison of the unwrapped interferograms and the coherence values used in GMTSAR and COMET-LiCSAR system: (a) The unwrapped phase pair 20170423–20170610. (b) The coherence pair 20170423–20170610.

The major difference is in the employed SBAS algorithm time series processing stage. In the implementation of the SBAS algorithm by GMTSAR, whenever there is one NaN value in any of the unwrapped phases, it will skip the whole column before doing the inversion and hence the pixel gets masked out. However, LiCBAS only removes that NaN value from the observation vector and will do the inversion using the remaining valid unwrapped phases. It should also be noted that even in the case of gaps in the small-baseline network (mainly due to decorrelation associated with factors such as vegetation or snow cover), LiCSBAS is able to do the inversion by imposing a linear temporal constraint across pixels with network gap. This explains the larger coverage of the LiCSBAS estimated velocity compared to the one obtained using GMTSAR. Particularly, in areas with relatively lower coherence, LiCSBAS seems to provide better estimate of the land subsidence.

While LiCSBAS is able to accurately detect large-scale land deformation (~ km) of the plain, its accuracy may be considerably affected by the land cover and use type of the area under investigation. Land subsidence measurements of the Ardebil plain obtained using the LiCSBAS and GMTSAR packages were used to evaluate how land use type and cover affect land subsidence estimation using the two packages. For this evaluation, the plain was classified into five different classes based on their landcover types (Fig. 12): Class A is a cropland area that has the highest potential for land subsidence; Class B is the City of Ardebil, classified as urban area with the largest population and located at the central part of the plain; Class C is the airport and its runways, which are located very close to the area with the most subsidence; Classes D and F mainly consist of rangeland and dry farming respectively. The displacement rates of these class regions were then estimated using LiCSBAS and compared with those obtained using GMTSAR.

Figure 12. Land cover type and classification of the Ardebil plain.

According to Table 2, the displacement velocities of LiCSBAS at cropland area for A-101 and D-6 track show subsidence at rates of 52 mm/yr and 47 mm/yr, respectively, which are consistent with GMTSAR. The Std of differences between the derived results in cropland area for both tracks are 0.09 mm/yr and 2.46 mm/yr, respectively. To evaluate consistency between the LiCSBAS and GMTSAR packages, a time series sample of subsidence behavior in Class A region was used.

Table 2

Information generated for mean LOS velocity maps of the S1A images at the landcover.
Area\Frame (S1A	A-101						D-6
	LiCSBAS			GMTSAR			LiCSBAS			GMTSAR
	Pixel	Vel Std (mm/yr)		Pixel	Vel Std (mm/yr)		Pixel	Vel Std (mm/yr)		Pixel	Vel Std (mm/yr)
Plain	128736	-52.73	6.67	99003	-53.22	6.58	43708	-47.51	6.64	20942	-46.44	4.18
A (Cropland)	38010	-51.83	7.65	9346	-51.19	8.16	10107	-47.51	8.34	1195	-46.44	7.94
B (Urban)	6363	-8.43	2.21	6508	-7.91	2.29	1606	0.93	1.31	1440	-4.52	3.84
C (Airport)	2590	-23.43	2.73	2407	-18.73	3.55	NaN	NaN	NaN	NaN	NaN	NaN
D (Rangeland)	15218	-14.21	3.88	8979	-20.59	5.95	14868	-8.87	2.96	4347	-13.26	3.11
F (Dry farming)	8019	-24.09	4.31	3277	-17.32	4.68	1518	-6.49	2.52	397	-13.57	3.84

Figure 13 shows the annual displacement time series of an area in Class A region located in the southeast part of the plain measured using LiCSBAS and GMTSAR for Tracks A-101 and D-6 for the period of March 2017 to January 2019. As presented in this figure, for the areas with expected high rates of subsidence, a good agreement is observed between the results of GMTSAR and LiCSBAS. The time series obtained from both packages indicate land subsidence with a relatively high rate of ~ 50 mm/yr. Both analysis packages show a reduction in the subsidence rate in 2019, that can be due to either implementation of some remediation measures including reduction of pumping draft and artificial recharge of aquifers from the land surface or the heavy rainfall in the region in 2019.

Figure 13. LOS displacement time series from S1A images (March 2017 to January 2019) for (a) A-101 track, (b) D-6 track, and correlation between the results of GMTSAR and LiCSBAS for (c) A-101 track, (d) D-6 track.

Figure 14 shows the Tropical Rainfall Measuring Mission (TRMM) and Standardized Palmer Drought Index (SPDI) of the region from 2017 to 2019. The SPDI is a monthly rainfall drought index based on calculation of the monthly cumulative and standardized rainfall anomalies [61]. SPDI is calculated using reference evapotranspiration and precipitation from the 4-km daily GRIDMET (Gridded Surface Meteorological) dataset and a static soil water holding capacity layer (https:// earth-engine/datasets/catalog/GRIDMET_DROUGHT). Dataset presented in Fig. 14 clearly show severe drought and reduced rainfall in 2018 with TRMM and SPDI of 0.02 mm/hr and − 4/8, respectively. Precipitation rate increased over the following years of 2019 and 2020, which resulted in a land subsidence decrease in the region.

Figure 14. LOS displacement time series from S1A images and average annual precipitation, (severe drought from March 2017 to January 2019).

In Class B region, as also presented in Table 2, differences in the displacement velocity and Std for A-101 and D-6 tracks are still small. To assess the LiCSBAS and GMTSAR results in Class B region, permanent GPS measurements were used (Fig. 15). The National Cartographic Center (NCC) of Iran (https://ncc.gov.ir), as a part of the Iranian Permanent Network for Geodynamic (IPGN) program, has established a permanent GPS station (ARDB) near Ardabil. Figure 12 shows the location of the ARDB GPS station by a red triangle. The green arrow indicates the velocity vector with an accuracy of 0.06 mm. Established in 2014, this station has been continuously measuring land deformation since then. Figure 15(a) shows a behavioral diagram of GPSUp (in vertical direction) measurements that indicates a maximum displacement rate of about 15 mm/yr from 2017 to 2019. The InSAR time-series for the location of the ARDB GPS station along with continuous GPSLOS (in LOS direction) measurements are illustrated in Fig. 15(b). Results indicated a very similar subsidence behavior for the analyzed time series, with a root mean square error (RMSE) of around 6 mm/yr. Also, a very good agreement was observed between the InSAR and GPSLOS time-series fluctuations, in terms of both range and data as observed in Fig. 15(b).

Figure 15. Comparison of the GPS and InSAR time-series of LiCSBAS and GMTSAR results: (a) Displacement trends from 2017 to 2019 of GPS permanent station (https://ipgn.ncc.gov.ir/pggn/). (b) RMSE between the InSAR and GPSLOS time-series from 2017 to 2019.

As presented in Table 2, the LiCSBAS displacement of Class C region for the A-101 track shows subsidence at a rate of 23 mm/yr which is in an acceptable range when compared to the 19 mm/yr subsidence rate measured using GMTSAR. In Class D and F regions, which are the areas with the lowest subsidence rates, the difference between the displacement rates estimated using the two packages were more pronounced. In these regions, displacement rate measurements estimated using LiCSBAS and GMTSAR differed about 6 to 7 mm/yr for the A-101 and D-6 tracks, respectively.

In order to evaluate the subsidence maps obtained from LiCSBAS and their compatibility with real data, it is necessary to adapt it to ground observations. In this study, relationship between land subsidence and groundwater changes, and a geotechnical-based numerical model were used to validate the results of land subsidence obtained using Envisat and S1A images by LiCSBAS.

Relationship between land subsidence and groundwater level (GWL) changes

One of the most important parameters for hydrogeological assessment of the land subsidence is the study of changes in groundwater level. Figure 16 shows contour maps of the groundwater levels and their comparison with land subsidence measurements in the region in 2006 and changes in groundwater level in Khalilabad during the years of 1991 to 2020. According to Fig. 16(a), GWL drops from east to west of the plain, with its lowest elevation (about 50 m below the ground surface) in the southeast part of the plain. Also results indicate a 27 m drop in GWL in Khalilabad during the years of 1991 to 2017 (Fig. 16b). As a result of these changes in GWL along with low rates of precipitation, great rates of subsidence were recorded that particular area of the plain. According to the results of GWL change, Envisat data (2004 to 2010) and S1A data (2014 to 2020) InSAR measurements, every meter drop in GWL has resulted in about 13 mm of annual land subsidence during the years of 2003 to 2020. Similar observations were reported by [39].

Figure 16. (a) Map for GWL changes and land subsidence of InSAR measurements. (b) The trend of groundwater level changes in one of the observation wells in the Khalilabad in with the time series of subsidence.

Geotechnical numerical modeling and comparison with SBAS-InSAR results

The land subsidence measurements obtained using SBAS-InSAR results were compared with the results of a rigorous geotechnical numerical analysis of land subsidence performed using FLAC2D in Ardabil plain between the years of 2006 and 2020. FLAC2D (Fast Lagrangian Analysis of Continua in 2 Dimensions) is a numerical modeling software that implements an explicit finite volume formulation (analogous to finite-difference for 2D geometries) to capture the complex behavior of geotechnical systems and soil-structure interaction. For numerical modeling, soil layers were divided into saturated and unsaturated layers based on the location of the groundwater level and then appropriate constitutive models were assigned to each layer to properly capture the stress-deformation behavior of the layers. Terzaghi's (1943) one-dimensional consolidation theory was used for the saturated layer and the constitutive model proposed by [63] was used for the unsaturated level. More details about [63] constitutive model were provided in Appendix A [62]. The input parameters for the constitutive models and methods that were implemented for their measurements are presented in Table. 3.

Table 3

The required parameters for the numerical modeling
Parameters
λ	SI*	k	SD*	λs	ν0	ks	Cν	k1	Sr0	k2	γd	P'c	Suction Stress	Groundwater Level
λ & k		These two parameters are obtained from the consolidation curve in which λ is the normalized consolidation slope in the loading mode and κ is the normal consolidation slope is in the case of loading in saturation condition of ν-ln p'.
λs & ks		λs is the initial wet and dry route slope and ks is the elastic section slope in Wheeler characteristic curve.
SI&SD		related to the air inflow (AEV) of the modified SWRC curve of the Wheeler model in two wet and dry modes.
P'c		known as soil pre-consolidation stress, whose value is obtained for consolidation test for different layers for soil saturation, but it should be calculated, for soils that are in the unsaturated state since the hardness or pre-consolidation stress of the soil increases by increasing the suction (decreasing the percentage of soil moisture).
k1 & k2		These two parameters are used to model the interaction between mechanical and hydraulic behavior, which moves the LC curve when it comes to SI or SD curves, and vice versa.
ν0 & γd		These two parameters are measurable according to the soil Index characteristics and two parameters of porosity ratio and soil Gs.
Cν		When a time-dependent settlement estimation is required, the coefficient of consolidation should be known. The amount of this parameter can be obtained using the one-dimensional consolidation test (edometer).

As presented in Table. 3, the input parameters include a series of hydraulic and mechanical properties that were mainly determined by performing a series of simple lab experiments on soil samples taken from the plain. In this study, parameters were calibrated based on properties of soil taken from a site near Khalilabad. This area has the lowest GWL and is expected to experience the highest rate of subsidence over the years. The location of the site is presented in Fig. 17. A rotatory wash-boring method was adopted to drill a 105 meter borehole in the field. During boring, a SPT (standard penetration test) split sampler was used to take soil samples at 5 meter depth intervals. The samples were then packed in rigid boxes and carefully transferred to the advanced soil mechanics laboratory, where index and strength tests were performed on the collected soil samples to measure input parameters for the numerical model. Five Decagon dielectric sensors were also inserted into the drilled borehole to capture and monitor profiles of soil temperature and moisture during an extended time period. The sensor measurements were then used to obtain profiles of suction and calibrate predictions of the numerical model following a procedure described by [64]. Based on the groundwater level measurements presented in Fig. 17, a groundwater level drop of 1 m/yr was assumed for the numerical analysis. Based on field measurements and sensor readings, the groundwater table was assumed to be at the depth of 50 m and a 150 m depth was considered as the depth of the bedrock. Details of soil sampling and model calibration are provided in [15].

Figure 17. Drilling site of a borehole 105 meters deep in the field.

Figure 18 compares predictions of the model over the years of 2006 to 2020 with land subsidence measurements obtained using SBAS-InSAR results. Based on land subsidence measurements, three time periods can be considered for the study: The first period is between the years of 2006 to 2015, where, according to Fig. 18(a), severe reduction in groundwater level (nearly 1.5 m/yr), resulted in considerable land subsidence in the region. From 2015 to 2018, the rate of changes in GWL decreased. However, due to long seasons of drought and significant increase in the SPDI drought index (Fig. 14b), the plain experienced significant land subsidence. Based on the results obtained using LiCSBAS, During the years of 2018 to 2020, a higher precipitation rate and some controlling measures that were taken by federal agencies decreased the rate of groundwater level drawdown to a value lower than 30 cm/yr. These changes justify lower rates of subsidence measured during the years of 2018 to 2020. The results of geotechnical numerical modeling showed land subsidence with an average annual rate of 3.8 cm between 2006 and 2020 which was close to what measured using SBAS-InSAR results. Results presented in this figure clearly show a good correlation between model predictions and InSAR measurements with an R² value of 0.77 (Fig. 18b).

Figure 18. Consistency evaluation of SBAS-InSAR results: (a) The subsidence values of each soil layer in the Ardabil plain between 2006 to 2020. (b) Correlation between model predictions and InSAR measurements.

In this study, an open-source interferometry time series analysis package, LiCSBAS that integrates with the automated Sentinel-1 InSAR processor was used to measure land subsidence of Ardebil plain, located at northwest, Iran, from October 2014 to January 2021. The results of this study show that LiCSBAS well performed in generating the Sentinel-1 time series over a full six years in the order of nearly 137 processed epochs for the A-101 and D-6 tracks. Using the LiCSBAS analysis package, a mean LOS velocity of 45 mm/yr was estimated for the Ardabil plain. Results also indicated subsidence rates as high as 50 mm/yr in the southeastern part of the plain. From 2018, due to higher rates of precipitation and some prevention measures, lower rates of subsidence were recorded. For consistency assessment of the land subsidence measurements using LiCSBAS, the annual displacement time series of an area near the city of Araloyebozorg (located at the southeast part of the plain) were obtained from March 2017 to January 2019 and then compared with the results of GMTSAR. LiCSBAS provided acceptable results of land subsidence in areas experiencing high rates of displacement and relatively low coherence. Comparison of the results of subsidence estimated using InSAR measurements and those obtained from GPSLOS measurements and a geotechnical based numerical modeling technique showed that LiCSBAS is able to accurately measure land subsidence, especially in areas that experience high rates of subsidence.

DATA AVAILABILITY

All the analyzed groundwater level data, climate change and geotechnical data, and GPS observations collected during this study are included in this paper. The InSAR products that support the findings of this study include the following:

LiCSAR products based on a predefined LiCSAR frame )from COMET-LiCS portal: https://comet.nerc.ac.uk/COMET-LiCS-portal/) (~ 25 GB in total): ID 101A_05193_131313 include 453 interferograms from (20141129 to 20200910):https://gws-access.jasmin.ac.uk/public/nceo_geohazards/LiCSAR_products/101/101A_05193_131313/interferograms/ and ID 006D_05111_131313 include 365 interferograms from (20141006 to 20210114): https://gws-access.jasmin.ac.uk/public/nceo_geohazards/LiCSAR_products/6/006D_05310_131313/interferograms/
Sentinel-1A data download for InSAR time series analysis by GMTSAR (from Copernicus Open Access Hub data portal: https://scihub.copernicus.eu/dhus/#/home) (~ 240 GB in total): Ascending_101 (20170330 to 20190131) and Descending_6 (20170324 to 20190113).
DEM files for use with GMTSAR (SRTM DEM 30 m): https://topex.ucsd.edu/gmtsar/demgen/). The ESA website has the orbit files: (https://qc.sentinel1.eo.esa.int/aux_poeorb/).

The generated datasets along with any models or codes that support the findings of this study are available from the corresponding author on reasonable request.

Calderhead, A., Therrien, R., Rivera, A., Martel, R., & Garfias, J. Simulating pumping-induced regional land subsidence with the use of InSAR and field data in the Toluca Valley, Mexico. Advances in Water Resources, 34(1), 83-97 (2011).
Foroughnia, F., Nemati, S., Maghsoudi, Y. & Perissin, D. An iterative PS-InSAR method for the analysis of large spatio-temporal baseline data stacks for land subsidence estimation. Int. J. Appl. Earth Obs. Geoinf, 74, 248–258 (2019).
Khorrami, M., Alizadeh, B., Tousi, E.G., Shakerian, M., Maghsoudi, Y. & Rahgozar, P. How Groundwater Level Fluctuations and Geotechnical Properties Lead to Asymmetric Subsidence: A PSInSAR Analysis of Land Deformation over a Transit Corridor in the Los Angeles Metropolitan Area. Remote Sensing, 377 (2019).
Qu, F., Zhao, C., Zhang, Q., Lu, Z., Yang, C. & Zhang, J. Land subsidence and ground ﬁssures in Xi'an, China 2005–2012 revealed by multi-band InSAR time-series analysis. Remote Sens. Environ. 155, 366–376 (2014).
Bhattacharya, A. K., Basak, S., & Patra, M. N. Land subsidence in Calcutta under the effect of hydrogeological conditions and over-extraction of groundwater. EJGE, 9, 1-16 (2004).
Wood, D.B. Water use and associated effects on ground-water levels, Las Vegas Valley and vicinity, Clark County, Nevada, 1980–95, Nevada Division of Water Resources Information Report 35 (1999).

Zhou, C., Gong, H., Chen, B., Li, J., Gao, M., Zhu, F. & Liang, Y. InSAR time-series analysis of land subsidence under different land use types in the Eastern Beijing Plain, China. Remote Sensing, 9(4), 380 (2017).

Poland, J. F., ed., Guidebook to studies of land subsidence due to groundwater withdrawal: Paris, France, UNESCO Studies and Reports in Hydrology, 305 p (1984).
Amighpey M, Arabi S, Talebi A, Djamour, Y Elevation changes of the precise leveling tracks in the Iran leveling network. Scientific report published in National Cartographic Center (NCC) of Iran, Tehran (In Persian) (2006).
Motagh, M., Djamour, Y., Walter, T.R., Wetzel, H.U., Zschau, J. & Arabi, S. “Land subsidence in Mashhad Valley, northeast Iran: results from InSAR, levelling and GPS.” Geophysical Journal International, vol. 168, pp. 518-526 (2007).
Baniani, S. R., Chang, L., & Maghsoudi, Y. Mapping and analyzing land subsidence for Tehran using Sentinel-1 SAR and GPS and geological data (No. EGU21-295). Copernicus Meetings (2021).
Rezaei, A., Mousavi, Z., Khorrami, F. & Nankali, H. Inelastic and elastic storage properties and daily hydraulic head estimates from continuous global positioning system (GPS) measurements in northern Iran. Hydrogeology Journal, 28(2), 657-672 (2020).
Ng, C. W. W., Sadeghi, H., Hossen, S. B., Chiu, C. F., Alonso, E. E. & Baghbanrezvan, S. Water retention and volumetric characteristics of intact and re-compacted loess. Canadian Geotechnical Journal, 53(8), 1258-1269 (2016).
Sadeghi, H. & Ali Panahi, P. Saturated hydraulic conductivity of problematic soils measured by a newly developed low-compliance triaxial permeameter. Engineering Geology. 278, 105827 (2020).
Mojtahedi, F. F., Ali, K., Nazari, A., Rezvani, S., Khatami, A. & Ahmadi, N. Measurement of Moisture and Temperature Profiles in Different Layers of Soil. IFCEE. pp. 266-278 (2018).
Galloway, D.L. & Burbey, T.J. Review: Regional land subsidence accompanying groundwater extraction. Hydrogeol. J. 19, 1459–1486 (2011).
Zhu, L., Gong, H.L., Jing, L.L., Su, Y.M., Li, X.J. & Jiang, J. Research on evolution of land subsidence induced by nature and human activity by utilizing remote sensing technology. In Proceedings of the Urban Remote Sensing Event 2009, Shanghai, China, 20–22. pp. 1331–1335 (2009).
Shanker, P., Casu, F., Zebker, H. A. & Lanari, R. Comparison of persistent scatterers and small baseline time-series InSAR results: A case study of the San Francisco bay area. IEEE Geoscience and Remote Sensing Letters, 8(4), 592-596 (2011).
Ferretti A, Prati C, Rocca F. Permanent scatterers in SAR interferometry. IEEE Transactions on Geoscience & Remote Sensing. 39(1), 8–20 (2001).
Lanari, R., Lundgren, P., Manzo, M. & Casu, F. Satellite radar interferometry time series analysis of surface deformation for Los Angeles, California. Geophys. Res. Lett. 31 (23) (2004).
Berardino, P., Fornaro, G., Lanari, R., Sansosti, E. A new algorithm for surface deformation monitoring based on small baseline diﬀerential SAR interferograms. IEEE Trans. Geosci. Remote Sens. 40 (11), 2375–2383 (2002).
Chen, X., Tessari, G., Fabris, M., Achilli, V., & Floris, M. Comparison Between PS and SBAS InSAR Techniques in Monitoring Shallow Landslides. In Workshop on World Landslide Forum (pp. 155-161) (2020).
Zhou, L., Guo, J. M. & Li, X. Monitoring and analyzing surface subsidence based on SBAS-InSAR in Beijing region, China. In International Conference on Intelligent Earth Observing and Applications 2015 (Vol. 9808, p. 98081Y). International Society for Optics and Photonics (2015).
Raspini, F., Loupasakis, C., Rozos, D., Adam, N., & Moretti, S. Ground subsidence phenomena in the Delta municipality region (Northern Greece): Geotechnical modeling and validation with Persistent Scatterer Interferometry. International Journal of Applied Earth Observation and Geoinformation. 28, 78-89 (2014).
Bozzano, F., Esposito, C., Franchi, S., Mazzanti, P., Perissin, D., Rocca, A., & Romano, E. Analysis of a Subsidence Process by Integrating Geological and Hydrogeological Modelling with Satellite InSAR Data. In Engineering Geology for Society and Territory-Volume 5 (pp. 155-159) (2015).
Catalão, J., Nico, G., Lollino, P., Conde, V., Lorusso, G., & Silva, C. Integration of InSAR Analysis and Numerical Modeling for the Assessment of Ground Subsidence in the City of Lisbon, Portugal. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing. 9(4), 1663-1673 (2015).
Rezaei, A. & Mousavi, Z. Characterization of land deformation, hydraulic head, and aquifer properties of the Gorgan confined aquifer, Iran, from InSAR observations. Journal of Hydrology, 579, 124196 (2019).
Pepe, A., Bonano, M., Zhao, Q., Yang, T. & Wang, H. The use of C-/X-band time-gapped SAR data and geotechnical models for the study of Shanghai’s ocean-reclaimed lands through the SBAS-DInSAR technique. Remote Sensing, 8(11), 911 (2016).
Iannacone, J. P., Falorni, G. & Macdonald, B. The role of InSAR in detecting and evaluating geotechnical risk from ground deformation. In Conference: Risk and Resilience Mining Solutions (2016).
Ciampalini, A., Solari, L., Giannecchini, R., Galanti, Y., & Moretti, S. Evaluation of subsidence induced by long-lasting buildings load using InSAR technique and geotechnical data: The case study of a Freight Terminal (Tuscany, Italy). International Journal of Applied Earth Observation and Geoinformation. 82, 101925 (2019).
Maghsoudi, Y., Amani, R. & Ahmadi, H. A Study on Land Subsidence in West of Tehran Using Sentinel-1 Images and Persistent Scatterers Interferometry, 299-313 (2019).
Chatrsimab, Z., Alesheikh, A., Vosoghi, B., Behzadi, S. & Modiri, M. Land Subsidence Modelling Using Particle Swarm Optimization Algorithm and Differential Interferometry Synthetic Aperture Radar. ECOPERSIA. 8(2), 77-87 (2020).
Hao, Q. N., & Takewaka, S. Detection of land subsidence in Nam Dinh coast by DInSAR analyses. In International Conference on Asian and Pacific Coasts (pp. 1287-1294) (2019).
Sandwell, D., Mellors, R., Tong, X., Wei, M. & Wessel, P. GMTSAR: Open radar interferometry software for mapping surface deformation. Eos Trans. AGU. 92(28), doi:10.1029/2011EO280002 (2011).
Hooper, Andrew. A multi-temporal InSAR method incorporating both persistent scatterer and small baseline approaches. Geophysical Research Letters. 35(16), 96-106 (2008).
Morishita, Y., Lazecky, M., Wright, T. J., Weiss, J. R., Elliott, J. R. & Hooper, A. LiCSBAS: An Open-Source InSAR Time Series Analysis Package Integrated with the LiCSAR Automated Sentinel-1 InSAR Processor. Remote Sensing. 12(3), 424 (2020).
Loibl, D., Bookhagen, B., Valade, S. & Schneider, C. OSARIS, the open source SAR investigation system for automatized parallel InSAR processing of sentinel-1 time series data with special emphasis on cryosphere applications. Frontiers in Earth Science. 172 (2019).
Lazecký, M., Spaans, K., González, P. J., Maghsoudi, Y., Morishita, Y., Albino, F. et al. LiCSAR: An automatic InSAR tool for measuring and monitoring tectonic and volcanic activity. Remote Sensing. 12(15), 2430 (2020).
Ghorbani, Z., Joodaki, Gh. & Javadnia, E. Land subsidence analysis in the Ardabil plain using sbas InSAR technique. In Conference: 3rd TRIGGER Conference 2019, TRIGGERS_abstract:https://iasbs.ac.ir/~trigger/pages/TRIGGERS_abstract.pdf. 111, https://www.researchgate.net/publication/347976526_Land_subsidence_analysis_in_the_Ardabil_plain_using_SBAS_InSAR_technique (2019).
Ulma, T., Anjasmara, I. M. & Hayati, N. Atmospheric phase delay correction of PS-InSAR to Monitor Land Subsidence in Surabaya. In IOP Conference Series: Earth and Environmental Science (Vol. 936, No. 1, p. 012033). IOP Publishing (2021).
Senturk, S., Cakir, Z. & Ustundag, B. B. The potential of Sentinel-IA interferometric SAR data in monitoring of surface subsidence caused by over drafting groundwater in agricultural areas. In 2016 Fifth International Conference on Agro-Geoinformatics (Agro-Geoinformatics) (pp. 1-4). IEEE (2016).
Anjasmara, I. M., Yulyta, S. A., Cahyadi, M. N., Khomsin, Taufik, M., & Jaelani, L. M. Land subsidence analysis in Surabaya urban area using time series InSAR method. In AIP Conference Proceedings (Vol. 1987, No. 1, p. 020071). AIP Publishing LLC (2018).
Aslan, G., Foumelis, M., Raucoules, D., De Michele, M., Bernardie, S., & Cakir, Z. Landslide mapping and monitoring using persistent scatterer interferometry (PSI) technique in the French Alps. Remote Sensing, 12(8), 1305 (2020).
Isya, N. H., Niemeier, W. & Gerke, M. 3D estimation of slow ground motion using InSAR and the slope aspect assumption, a case study: the puncak pass landslide, indonesia. ISPRS Annals of Photogrammetry, Remote Sensing & Spatial Information Sciences. 4 (2019).
Suresh, D. & Yarrakula, K. InSAR based deformation mapping of earthquake using Sentinel 1A imagery. Geocarto International. 35(5), 559-568 (2020).
Sreejith, K. M., Sunil, P. S., Agrawal, R., Saji, A. P., Rajawat, A. S. & Ramesh, D. S. Audit of stored strain energy and extent of future earthquake rupture in central Himalaya. Scientific reports. 8(1), 1-9 (2018).
GAMMA Remote Sensing, GAMMA Software Information, https://gamma-rs.ch (2021).
Kord, M., Asghari-Moghaddam, A. & Nakhaeei, M. Investiga-tion of hydrogeological characteristics of Ardabil plain aquifer, Northwest of Iran. ISESCO Journal of Science and Technology 9(15), 63–69 (2013).
Ardabil Regional Water Authority (ARWA) Evaluation of groundwater in Ardabil plain. Ardabil Regional Water Authority, https://www.arrw.ir (2014).
Aghazadeh, N., M. Chitsazan, & Y. Golestan. Hydrochemistry and quality assessment of groundwater in the Ardabil area, Iran. Applied Water Science 7, no. 7, 3599-3616 (2017).
Chen, C.W. & Zebker, H.A. Phase unwrapping for large SAR interferograms: Statistical segmentation and generalized network models. IEEE Trans. Geosci. Remote Sens, 40, 1709–1719 (2002).
Doin, M.-P., Lodge, F., Guillaso, S., Jolivet, R., Lasserre, C., Ducret, G., et al. Presentation of the small baseline NSBAS processing chain on a case example: The Etna deformation monitoring from 2003 to 2010 using Envisat data. In Proceedings of the ESA Fringe 2011 workshop (pp. 19–23) (2011).
López-Quiroz, P., Doin, M.-P., Tupin, F., Briole, P. & Nicolas, J. M. Time series analysis of Mexico City subsidence constrained by radar interferometry. Journal of Applied Geophysics, 69, 1–15 (2009).
Morishita, Y. Nationwide urban ground deformation monitoring in Japan using Sentinel-1 LiCSAR products and LiCSBAS. Progress in Earth and Planetary Science, 8(1), 1-23 (2021).
Biggs, J., Wright, T., Lu, Z., & Parsons, B. Multi-interferogram method for measuring interseismic deformation: Denali Fault, Alaska. Geophysical Journal International. 170(3), 1165–1179 (2007).
De Zan, F., Zonno, M., & López-Dekker, P. Phase Inconsistencies and Multiple Scattering in SAR Interferometry. IEEE Transactions on Geoscience and Remote Sensing, 53(12), 6608–6616 (2015).
Hooper, A. & Zebker, H. Phase Unwrapping in Three Dimensions with Application to InSAR Time Series, J. Optical Soc. of America. 24, 2737-2747 (2007).
Hanssen, R. F. Radar Interferometry: Data Interpretation and Error Analysis. Vol. 2, (Springer Science & Business Media, Berlin (2001).
Aalipour, E. M., Malekmohammadi, B., & Jafari, H. R. Risk Zoning of Land Subsidence Due to Groundwater Level Declining using Fuzzy Analytical Hierarchy Process (Case Study: Ardabil Plain). Iranian journal of watershed management science and engineering . Volume 11 , Number 38. 25-34 (2017).
Sandwell, D., Mellors, R., Tong, X., Xu, X., Wei, M. & Wessel, P. GMTSAR: An InSAR Processing System Based on Generic Mapping Tools (2016).
McKee, T.B., Doesken, N.J. & Kleist, J., The relation-ship of drought frequency and duration to time scales,Proc. Eighth Conf. on Applied Climatology. 17–22, 179–184 (1993).
Terzaghi K. Theory of consolidation. Theor Soil Mech. 265–296 (1943).
Wheeler, S. J., Sharma, R. S., & Buisson, M. S. R. Coupling of hydraulic hysteresis and stress–strain behaviour in unsaturated soils. Géotechnique. 53(1), 41-54 (2003).
McCartney, J. S. & Khosravi, A. Field-monitoring system for suction and temperature profiles under pavements. Journal of performance of constructed facilities. 27(6), 818-825 (2013).

No competing interests reported.

Use of InSAR Data for the Measurement of Land Subsidence in Response to Groundwater Fluctuations and Climate Change (Case Study: Ardabil Plain, Iran)

Status:

Version 1

Abstract

Figures

Introduction

Description Of Study Area And Research Methodology

Description of Study area

Research methodology

InSAR data

LiCSBAS time-series analysis

Results

Interpretation of the S1A data by LiCSBAS velocities

Consistency assessment

Discussion

Relationship between land subsidence and groundwater level (GWL) changes

Geotechnical numerical modeling and comparison with SBAS-InSAR results

Conclusion

Declaration

References

Additional Declarations

Supplementary Files

Status:

Version 1