InSAR-derived land subsidence in Wuhan between 2015 and 2020

ABSTRACT Land subsidence has become a challenging problem with urbanization and underground exploitation. The land subsidence in Wuhan, China, was inferred from 47 Sentinel-1A images acquired between June 2015 to October 2020 based on Synthetic Aperture Radar (SAR) Interferometric Point Target Analysis (IPTA) technology. Results show that significant subsidence mainly occurred in Houhu, Qinshan, Baishazhou, and Zhuankou areas, and the most significant subsidence occurred in Baishazhou, reaching an annual subsidence rate of 41.0 mm/a with accumulated subsidence of 219.4 mm., The subsidence was more significant in the carbonate solution zone with the overlying Holocene strata compared to the rest carbonate karst belts. Analysis shows thatland building construction is a major reason for land subsidence in Wuhan. The construction of subway underground stations is another factor causing land subsidence. The effect of water level fluctuations of the Yangtze River on long-term trends is negligible. Our study provides further insight into the temporal-spatial pattern of Wuhan’s displacements and suggests that the land subsidence should be continuously monitored for alleviation of potential destructive hazards, especially during the period of construction of a land building, such as residential, public facilities, and subway.


Introduction
Land subsidence is a gradual or sudden sinking of the Earth's surface due to theremoval or displacement of subsurface earth materials, such as underground mining, groundwater depletion, etc. (Galloway et al., 2000). A recent study showed that land subsidence had become a serious global environmental geological problem, which was threatening urban infrastructure and normal human activities (Herrera-García et al., 2021). Herrera-García et al., (2021) also pointed out that subsidence caused by groundwater depletions had occurred at 200 sites in 34 countries (Iran, Mexico, India, etc.) over the past century. Haghshenas Haghighi & Motagh, (2019) obtained the surface displacement time series of Tehran, Iran from 2003 to 2017 by using InSAR technology, and the results showed that Tehran had an obvious subsidence trend. Chaussard et al., (2021)    also obtained land subsidence in Mexico City based on the 2014-2020 Sentinel-1 data, and found that more than 457,000 properties and 1.5 million residents in the ZMVM were at high risk of surface rupture due to the land subsidence. Based on 60 Sentinel-1 data Karanam et al., (2021) investigated the subsidence in the Jariya Coalfield India at a rate of 10 cm/a to 20 cm/a and discussed the impact of coal fires on the subsidence, resulting in the subsidence velocity had a strong positive correlation with the temperature anomalies caused by coal fires. Garg et al., (2022) used sentinel data from 2014 to 2020 to obtain the land subsidence rate in Delhi, India and found that the subsidence in the kapashela and Faridabad districts was increasing. In addition, China also faces huge risks of land subsidence. In 2015, the area affected by land subsidence in China (more than 10 mm/a) exceeded 120,000 km 2 , including more than 57 prefecture-level (or above) cities, and the area affected by serious subsidence (more than 50 mm/a) was 5860 km 2 (L. Zhang, 2017).
In the monitoring techniques of land subsidence, although the traditional geodetic methods (e.g. leveling survey, GPS) have high accuracy (Abidin et al., 2008;Serrano-Juan et al., 2017), they are still timeconsuming and labor-consuming methods, and are not efficient in large-scale continuous monitoring . As an active remote sensing technology, InSAR is a reliable tool for large-scale surface deformation monitoring because of its largescale, all-weather, remote sensing capability (Samsonov et al., 2013;Simons, 2002;Wright, 2004). Differential Interferometric Synthetic Aperture Radar (D-InSAR) was first used to obtain the land subsidence of an irrigated area in California, USA (Bai, 2017;Gabriel et al., 1989). However, the D-InSAR is highly susceptible to decoherence phenomena due to spatiotemporal varying features of observation targets, and atmospheric delay. In order to overcome the spatiotemporal decoherence effect, Permanent Scatterers InSAR technology (PS-InSAR), a multi-temporal technique, was firstly proposed in the late 1990s based on a set of SAR acquisitions which allowed to identify high coherence point targets, mitigate atmospheric effects and retrieve time series displacements (Ferretti et al., 2001). PS-InSAR technology can largely reduce the influence of temporal and spatial decoherence on signals by only using the high coherence of PS points. The method improves the InSAR-derived displacement by estimating and mitigating the phase contribution of atmosphere change based on the differential spatio-temporal features of deformation and atmospheric phases. The PS-InSAR technology has a potential to detect millimeter surface change(D.R. Li et al., 2004;Liao & Lin, 2002;Liao & Wang, 2014), which has been widely used to detect changes due to surface deformation (e.g. urban subsidence), artificial buildings(e.g. dams and railways), and natural objects (e.g. bare rocks; Gu et al., 2014;Li, 2021;J. Zhang et al., 2021). To identify possible future land subsidence areas, susceptibility mapping is a widely used method to hazard mitigation using mathematical or statistical models (Arabameri et al., 2022;Dong, 2006;Ghorbanzadeh et al., 2020Ghorbanzadeh et al., , 2018, such as artificial neural network (ANN; Lee et al., 2012), decision tree (Lee & Park, 2013), k-nearest-neighbor analysis (K-NN; Pradhan & Jebur, 2017), and adaptive neuro-fuzzy inference systems (ANFIS; Arabameri et al., 2022;Basser et al., 2014;Chen et al., 2017;Polykretis et al., 2019;Sezer et al., 2011;Tien Bui et al., 2012). Recent studies have shown more accurate susceptibility prediction using artificial intelligence (AI) models (e.g. ANFIS) in combination with highquality training and validation data from InSAR (Chatrsimab et al., 2020;Ghorbanzadeh et al., 2020).
As a mega city in central China, Wuhan City has experienced rapid infrastructure development in recent years, such as the subway railway, which has undoubtedly led to large-scale land subsidence (Y. Zhang, 2019). By 2016, land subsidence has become one of the main geological disaster types in Wuhan (Liao et al., 2019). There are several studies on land subsidence monitoring in Wuhan using multi-temporal InSAR technology. Bai et al., (2019) obtained the monitoring results of land subsidence in Wuhan from 2013 to 2015 by using StaMPS based on TerraSAR-X RADAR data. The results showed that the most serious land subsidence occurred in the Houhu area of Hankou, and the subsidence rate exceeded 78 mm/a. The main causes of land subsidence were urban construction, groundwater level, and land loading on the soft base surface. Sun et al., (2019) used the PS-InSAR method to obtain land subsidence rate inferred from COSMO-SkyMed data covering the main urban area of Wuhan from June 2012 to June 2017. They found that the most serious land subsidence occurred in the Houhu area of Hankou, reaching a rate of 60.9 mm/a. (Tien Bui et al., 2012) used ALOS-PALSAR data from 2007 to 2011 and Radarsat-2 data from 2015 to 2019 to extract the corresponding land subsidence rate and time series in the main urban area of Wuhan during these two periods, and quantitatively analyzed the spatio-temporal pattern of land subsidence. The occurrence and evolution of land subsidence were also investigated in combination with precipitation, soft soil, and urban construction. This study still showed that the land subsidence in Wuhan between 2007 and 2011 was mainly concentrated in the Hankou area, the coastal area of Shahu Lake, the Baishazhou area, and the coastal area of Yangtze River to the north of Shahu Lake. This study concluded that the main reasons for the land subsidence were the widespread distribution of quaternary and soft soil, underground space development as well as engineering construction. From 2015 to 2019, the widespread distribution of the quaternary system and construction in Hankou, northern Shahu Lake, and the Baishazhou area were the main driving factors of land subsidence in Wuhan City. Han et al., (2020) analyzed the causes of land subsidence in Wuhan with ALOS-1, Envisat, and Sentinel-1 images obtained from 2007-2010, 2008-2010, and 2015-2019 respectively. The result showed that human activities such as intensive urban construction and groundwater extraction were the main causes of subsidence and the subsidence center would change significantly from the change in temporal and spatial characteristics. This study also suggested that the seasonal variation of water levels may affect the land subsidence along the Yangtze River in Wuhan. Shi et al., (2021) used the SBAS-InSAR (Small Baseline Subset InSAR) method to obtain the land subsidence rate in Wuhan based on Sentinel-1A data from 2015 to 2019, and considered that the current rapid urban construction was the main cause of land subsidence. Jiang et al., (2021) found that Houhu, Xinrong, and Guanggu were the most serious area of land subsidence in Wuhan by using long-time COSMO-SkyMed data. The study also pointed out that the subsidence of Houhu was caused by the soft soil consolidation, while the subsidence in Xinrong and Guanggu wascaused by human factors, such as subway construction, and urban construction, by analyzing the influence of natural factors and human factors on land subsidence and their correlation.
However, the current studies focus only on the land subsidence occurred before 2019, and other factors (e.g. subway construction of elevated station and underground stations) which have a potential to lead to land subsidence are not comprehensively investigated. Therefore, we use inbuilt Interferometric Point Target Analysis (IPTA) method of GAMMA software to obtain the time series displacements and estimate average annual rate of land subsidence in Wuhan using 47 scenes Sentinel-1A data acquired from 2015 to 2020, and compared with the subsidence results of benchmark along the second phase of Line 8. Several main subsidence centres have been comprehensively analysed, and corresponding impact factors of land subsidence are also investigated, especially the water level of the Yangtze River and subway construction. In the study of the water level of the Yangtze River, in order to reduce the influence of time delay (Han et al., 2020),the accumulative deformation was used to analyse.

IPTA principle
IPTA is one of the concrete realisation methods of PS-InSAR technology (Werner et al., 2003). This method inherits the idea of PS-InSAR technology, and in IPTA the interferograms are only interpreted for the selected coherent point targets. Compared with traditional interferometry, IPTA uses vector data structure rather than raster data structure for differential interferometry and regression analysis, which allows a drastic reduction of required disk space, resulting in improved computing efficiency and disk space savings (H.B. Zhang et al., 2016;Liu et al., 2019;Wang & Xu, 2017). The important aspect of IPTA is point targets in long-baseline interferometric phase of pairs can be also analysed, even above the critical baseline. A flow chart of typical IPTA processing is shown in Figure 1.
The core algorithms consist of the generation of a point list, differential interferogram, interferometric point analysis, and model refinement, including atmospheric phase estimation and removal. The main steps include the following: (1) Point list generation. Based on the co-registered Single Look Complex (SLC), a candidate list of point targets is determined according to the stable and pointlike scattering characteristics. The initial selection of point target candidates is based on criteria high backscattering and low spectral phase diversity (GAMMA, 2003; 2) Differential interferogram. In this stage, SLC values are extracted and written into the point list data stack. Initial interferometric baselines are estimated from the orbital state vectors, and differential interferograms are computed. The interferometric phase ɸ can be decomposed into: where ɸ flat is the flat phase; ɸ def is the deformation phase; ɸ noise is the noise phase; ɸ atm is the atmospheric phase. Then, the initial interference baselines and Digital Elevation Model (DEM) data are used to simulate the terrain phase, and subtract it from the interferometric phase to get the initial differential phase.
(3) Initial phase regression analysis. Establish a regression model with the differential interferometric phase as the dependent variable, the linear deformation phase and the elevation correction as the independent variables, and the residual phase as a constant (G.Y. Li et al., 2018). The least square method is used to solve the correction values elevation and deformation rate, and the phase model is updated.
(4) Regression analysis iteration. Due to the influence of terrain error, elevation error, orbit error, etc., the residual phase is usually included in the initial differential interferograms. Therefore, iterative regression analysis is needed for the updated phase model to improve the accuracy of solving deformation parameters. Firstly, the minimum cost flow (MCF) method is used for multi-batch unwrapping in the space domain, and two-dimensional regression analysis is used to solve the correction of elevation error and deformation rate iteratively in the time domain. A least square method is used to correct baseline errors.
(5) Residual phase separation. The residual phase is the sum of the nonlinear deformation phase, atmospheric phase, and noise phase, i.e. ɸ res ¼ɸ non def þɸ atm þɸ noise (2) Where ɸ res represents the residual phase; ɸ non def represents the nonlinear deformation phase; ɸ atm represents the atmospheric phase; ɸ noise represents the noise phase. The atmospheric phase can be divided depending on its different spatio-temporal characteristics.

Validation of the result with leveing
In the IPTA analyses, all displacements were in the direction of Line-Of-Sight (LOS), and the displacement using leveling was in the vertical direction. In order to compare the results of IPTA and leveling, the movements measured in the IPTA should be converted to vertical. The transformation formulate (Fryksten & Nilfouroushan, 2019) was shown in Equation 3 in which no horizontal displacements are assumed, where the m v was the movements in the vertical direction, the m LOS was the movements in the LOS direction, and the θ was the incidence angle. The incidence angle is 39.3, and the ratio between in a length in the LOS direction and vertical direction was 1.29.

Study area
Wuhan is situated in the east of Jianghan Plain and is the capital city of Hubei Province in China. The Yangtze River and Han River divide Wuhan into three parts: Wuchang, Hanyang, and Hankou, which are generally known as Wuhan's Three Towns. The area of carbonate karst is distributed widely, accounting for about 13% of the total urban area (Luo, 2014), resulting in frequent karst ground collapses (Fan, 2006). The karst geology in Wuhan is dominated by covering and buried types (Guan et al., 2008). The geographical location of the study area is 113.8396° to 114.7059° E and 30.2567° to 30.7941° N, covering the main urban areas of Wuhan (Hanyang District, Qiaokou District, Jianghan District, Jiangan District, Qingshan District, Wuchang District, and Hongshan District), part of the suburban and surrounding areas. Figure 2 shows the geographical location and karst distribution of the study area. There are six carbonate rock belts with five types of karst geology in the main urban areas of Wuhan (Luo, 2014(Luo, , 2013Luo et al., 2018;X Wang et al., 2020). The Baishazhou karst belt is the area with the most frequent karst collapse accidents, and 27 karst collapses were recorded in 2018 in the region (Y. .

Data sources and pre-processing
We used 47 Sentinel-1A TOPS (Terrain Observation with Progressive Scans) SAR images in interferometric wide (IW) swath mode acquired from June 2015 to October 2020. The data is Ascending SLC images with VV polarisation. The registration images are multi-looked in both range and azimuth directions (5 × 1 pixels) to improve the SAR image quality by reducing the speckle and to obtain a square pixel on the output image. The image acquired on 18 June 2018 is selected as the master image and the rest 46 images are slave images. Finally, 46 interferometric pairs are generated for the subsequent time series analysis. Figure 3 shows their spatio-temporal baseline distribution. Most perpendicular baselines are less than 50 m, which suggests a higher coherence. In addition, the DEM with a spatial resolution of 30 m directed from the Shuttle Radar Topography Mission (SRTM) is used to remove the influence of terrain effects in the interferometric phase.

IPTA analysis and displacementresults
There are 793,644 points targets to be selected in the point list generation stage by the thresholds of the amplitude deviation which is measured by the Mean to Sigma Ratio (MSR) (MSR = 1.7) and the coherence coefficient (γ = 0.5). The point target at the new archway of Wuhan University is selected as the reference   point, and the displacement of IPTA is LOS direction. Figure 4 shows the time series displacements of the main urban area of Wuhan from June 2015 to October 2020. The warm color (positive) represents the land subsidence relative to the master image, while the cold color (negative) represents the land uplift relative to the master image. The land subsidence mainly occurred at the intersection area of the Han River and the Yangtze River, and on the coast of the Yangtze River on the right side. The displacement at the right bank of the junction of the Han River and the Yangtze River decreased in time series, the right bank of the Yangtze River also shows significantly continuous land subsidence in the time period.
To further analyze the spatio-temporal distribution characteristics of the study area, we obtained the average rate of land subsidence in the main urban area of Wuhan between 2015 and 2020 ( Figure 5). Similar to previous studies, significant subsidence occurred in areas where subsidence has been found Y. Zhang, 2019), and the subsidence rates range from −41.0 mm/a-12.0 mm/a. The measurements of the land subsidence rate exceeding 5 mm/a are less than 3%, which shows that the main urban area of Wuhan is relatively stable as a whole during the time period. Compared with previous studies Sun et al., 2019), Houhu in Hankou, once a severe subsidence area, now shows a slowing trend of subsidence rate, and the severe subsidence center has moved northward to the area of active construction. Other subsidence areas are mainly situated in Baishazhou, Zhuankou, Qingshan District, and other places. The severest subsidence occurred in Qinglinghe Road in Baishazhou, with an average annual subsidence rate of −41.0 mm/a and accumulative subsidence of −219.4 mm. During the time period, Baishazhou is in a period of large-scale construction, for example, the construction of the Yangsigang Yangtze River Bridge and the construction of commercial residential areas (e.g. Zhongjian Platinum residence), etc.
In addition, land subsidence still occurred in some suburbs. Some uplift areas in the previous study have also become subsidence (Y. Zhang, 2019), and were shown marked by the red box in Figure 5, which is due to the large-scale urban construction from August 2019 to October 2020. We further analyze the land subsidence results over several main subsidence areas (marked by blue boxes in Figure 5) and investigate the potential impact factors, including the water level of the Yangtze River, and the construction of Subway Line 16.

The validation of InSAR-dervied displcaments
In order to assess the accuracy of displacement measurements, IPTA results obtain from 22 June 2015 to 14 October 2016 were compared to land subsidence which was measured by leveling along Subway Line 8 observed from 30 June 2015 to 13 September 2016 (Y. Zhang, 2019). These leveling measurements are surveyed twice in the second period of Subway Line 8. Due to the spatial discontinuity of IPTA, Kriging interpolate was used in ArcGIS to obtain the IPTA results at the location of the leveling. Figure 6 shows the IPTA results and their corresponding leveling results. The differences between IPTA and leveling were less than 10 mm, and most of them were between ±5 mm. For quantitation of the validation, Spearman's rank correlation coefficient was used for correlation analysis, since the leveling results do not conform to a normal distribution. Unlike the Pearson correlation coefficient, this method can work for data that do not conform to a normal distribution. The statistical information on differences between IPTA, PS-InSAR, and leveling is shown in Table 1. The IPTA results were weakly correlated with leveling results, but there is a correlation of 0.457(P-value<0.001) with PS-InSAR results in the paper by (Y. Zhang, 2019). The mean and standard deviation were better than the results in Zhang.
The significant difference between IPTA results and leveling may be attributed to four factors. (1) The IPTA values are computed by use of the interpolation method, which shows the average change over a certain area, while the leveling measurements are obtained on a certain site. (2) The penetration of SAR makes the IPTA measurements represent integrated change over a certain depth from the surface, while the leveling only measures the surface change.
(3) The time asynchrony between the InSAR observation and leveling survey also contributes to the difference, the time differences are generally more than 30 days. (4) The error of the two techniques may be another factor. In summary, the displacements of this experiment are considered reliable (Shi et al., 2021).

Relationship between subsidence and karst ground
The risk of karst collapse is different in carbonate karst belts, due to the different overlying rock and soil mass. In order to understand the influence of different carbonate karst belts on displacement, we analyzed the annual average displacement rate of the PS points according to the karst distribution at the location of the PS points.
It can be seen that the areas along the river in the Baishazhou karst belt and Hannan karst belts (Figure 2) show significant land subsidence (hatched areas in Figure 7). In Figure 8, the median and mean values of PS points on the carbonate karst belts were close to 0, and the values between the first and third quarter were found between −0.4 and 1.9 mm, except for the Baishazhou karst belt and Hannan karst belts which were between −1.3 ~ 0.3  and −0.9 ~ 1.0, respectively. It was worth noting that, compared with the overlying geology of other belts, the banks of the Baishazhou karst belt and Hannan karst belts were Holocene strata coverage areas, and the silt sand is directly above the carbonate strata. The sandy soil can be lost directly through the channels such as dissolution gaps and holes, which will cause the ground to collapse and cause karst geological disasters.

Houhu
The Houhu was covered mainly by lakes in history, and the soil is mostly mucky clay or silt. In recent years, the reclamation of Houhu has covered its original soil type. Large-scale infrastructure and high-rise building construction increased the land loading in this area, resulting in groundwater decline, which further exacerbate the land subsidence in the Houhu area(Bai et al., 2019;   Shi et al., 2021;Zhen et al., 2003). According to the previous monitoring results of movement in Wuhan Shi et al., 2021;Zhen et al., 2003), the subsidence rate of Houhu was relatively large, reaching 78.1 mm/a, during the construction period of the Subway Line 3 between 2012 to 2015. In addition, from the end of 2013 to the beginning of 2015, Houhu also showed a tendency for land subsidence rapidly (Jiang et al., 2021). Figure 9 shows the annual subsidence rate distribution and movement time series of three sites (denoted by HH1, HH2, and HH3) in the Houhu area. After the subway construction, the land subsidence of Houhu tends to be moderated during the monitoring period from 2015 to 2020. The subsidence center has moved towards the areas with more active urban construction. Figures 10(a) and Figure 10(b) are the satellite images acquired in July 2015 and February 2020 respectively, corresponding to the area marked by the blue frame in Figure 9. The area with building change delineated by the white lines is often relative to more serious land subsidence.

Qingshan
The soil in the Qingshan area is composed of loose clay, sandy clay, sandy soil, and a small part of weathered limestone. Figure 11 shows the distribution of annual average land subsidence and displacement time series of two sites (denoted by QS1 and QS2) in the Qingshan Area. The displacement results reveal a significant trend of subsidence. The distribution  patterns of subsidence centers are generally consistent with that of a previous study from 2015 to 2019 (Shi et al., 2021), and the subsidence occurred in the areas with silt or soft soil types. The area at the QS1 point   from May 2018 to January 2020 shows accelerated subsidence. Figure 12a,b are satellite images acquired in December 2017 and February 2020 respectively, and the white frames delineate a new area of high-rise building construction. The accelerated subsidence is well in agreement with the period of the construction, which suggests a significant effect of building construction.

Baishazhou
Baishazhou is an area where ground collapse frequently occurs owing to geological conditions and large-scale construction. Figure 13 shows the annual average land subsidence distribution of Baishazhou and the Displacement time series of four sites (denoted by BSZ1, BSZ2, BSZ3, and BSZ4). The site of BSZ2 is located in the Baishazhou carbonate band (Fan, 2006), it can be seen from Figure 7, that its land subsidence is more likely related to karst geology. The area at BSZ1, BSZ3, and BSZ4 occurred larger land subsidence and the most serious subsidence occurred in the areas marked by BSZ1 and BSZ4, which is possibly contributed by large-scale construction during the monitoring period. Figures 14(a) and Figure 14(b) are satellite images acquired in August 2015 and December 2019. There are many new high-rise buildings along the Qingling River to the south during 2015-2019, which is in good agreement with the distribution of land subsidence along the Qingling River. It should be noted that excessive land subsidence can cause collapse. Due to the poor and unstable geological conditions of the miscellaneous fill, and illegal construction, a sudden collapse event occurred in the area south of Baisha Fifth Road (close to BSZ4) on 13 June 2018(Investigation Report,).

Zhuankou
The Zhuankou area is located in the lacustrine engineering geological subregion. The surface of the land is covered by loose artificial fill, and the bottom consists of soft-flowing mucky silt and soft plastic clay, which has poor engineering properties and poor selfstability (C.A. Li et al., 2019). Figure 15 shows the annual average land subsidence and the displacement time series at the sites of ZK1 and ZK2 in the Zhuankou area. The area is generally relatively stable, or with slight uplift, while the area around ZK1 shows an obvious trend of subsidence compared with the surrounding area. Figure 16 shows satellite images acquired in August 2015 and October 2020 in coverage of the vicinity of ZK1. The building area shown by the white line was located in the subsidence area and there was no apparent construction activity during the monitoring period. We found that there are two logistics centers and Wuhan Jiantou Railway Transportation Co., LTD in the area with subsidence. The logistics companies and rail transport co., LTD. are mainly engaged in warehousing, cargo handling, and transportation during the period, resulting in frequent loading on the ground. Since the land type in this area is lacustrine facies, which suggests that its land loading capacity is weak. Therefore, we consider land subsidence in this area is related to its karst geology and frequent loading.

Influence of subway construction on land subsidence
Subway Line 16 has a length of 33.0 kilometers with 12 stations. The underground stations are situated between Guobo Center South station and Tanjun Road station, and elevated stations are located between Qingjiang station and Zhoujia River station. In order to investigate the impact of the construction of Subway Line 16, the two displacement datasets of the period before construction (From 9 August 2015 to 3 August 2016) and the period of construction (from 24 July 2018 to 19 July 2019) were identified respectively. Figure 17 shows the accumulative land subsidence during the non-construction period and the construction period of Subway Line 16 respectively. The maximum accumulative land subsidence during the non-construction period is about 32.8 mm, and the land subsidence occurred in the area between Guobo Center South station and Zhuankou Station. The Guobo Center South station and Laoguancun station are the transfer stations for the first phase of Subway line 6, and the construction period of line 6 was from 7 August 2013 to 28 December 2016. Therefore, in the non-construction period, the monitoring sensitivity of the Guobo Center South station to Yuejiazui Station is relatively strong, while land subsidence monitoring sensitivity is low where the station is not within the monitoring range in the southeast.
Similar to the situation in the non-construction stage phase, land subsidence between Guobo center South station to Zhuankou station still prevails during the construction period of subway line 16, but the area of subsidence is expanded significantly. Compared with the non-construction phase, increased subsidence is found in the construction phase of the underground station, and the maximum land subsidence between underground stations reaches37.8 mm. The land subsidence from Mayinghe station to Xiezihe station increased. However, no obvious land subsidence was found at the Zhijiaoyuan station, the station is an elevated station. Figure 17c   area changed significantly during this period. The land subsidence in this area was also considered to be related to the construction of high-rise buildings, while no obvious land subsidence is found during the construction period of the elevated station.
Phase 1 of Subway Line 3 and Phase 1 of Subway Line 6 were completed in December 2015 and December 2016 respectively (Figure 2). To further investigate the impact of underground subway station construction on land subsidence, we analyzed the land subsidence in subway line 3 (Wang Jiawan station -Longyangcun station), subway Line 6 (Jiangang Station -Qianjincun Station), and subway line 16 (Laoguancun Station-Laoguancun Depot Station) with similar the karst geologic types. The movement in the 100-meter buffer zone between the above stations was averaged and its time series of land subsidence is shown in Figure 18. The surface between Wangjiawan station and Langyangcun station (Line 3) showed a slow upward trend, while no obvious surface changes were observed between Jiangang station and Qianjincun station (Line 6). But in the area between Laoguancun station and Laoguancun Deport station, (Line 16) significant subsidence was found during the construction period, reaching up to 4 mm. It can be considered that the land subsidence between Laoguancun station and Laoguancun Deport station is attributed to the construction of the underground station of subway Line 16.

Influence of Yangtze River water level on surface movement
In order to assess the impact of water level change on beaches of the Yangtze River, we obtain the data of the water level on the Yangtze River from June 2015 to October 2020 at the Hankou hydrological station from the official website of China MSA (China MSA, 2022). PS points selected 400 meters away from Hankou hydrological station are used to investigate the relationship between land subsidence and water level changes. Figure 19 shows the average subsidence of within 200 m from the PS point as shown in Figure 5 and the water level of Yangtze River. In order to avoid the impact of construction on land subsidence, there is no obvious construction in the vicinity of PS point during the monitoring period, and the continuously accumulative subsidence of 8 mm was observed in the area during the observation period. The water level of the Yangtze River changes periodically, and the water level of the Yangtze River reaches its peak between July and September  every year. The ground shows a slight uplift of approximately 2 mm when the water level of the Yangtze River reaches its peak.

Discussion
The cause of land subsidence is the result of the combination of human activities and natural factors. This paper analyzed the effects of human factors such as urban construction and subway construction, as well as natural factors such as the distribution of carbonate karst geology and changes in the water level of the Yangtze River on land subsidence. Compared with previous studies, the subsidence center seems to have transferred to Baishazhou of Wuchang district from Houhu of Hankou, and the magnitude of land subsidence at Hankou became more moderately, which is consistent with the transfer of urban construction center in the time period.
We analyzed the relationship between subsidence and karst ground, and we found that highrisk areas for karst geological disasters occur in the areas of the carbonate strata directly covered by the Holocene strata. Unlike the Baishazhou karst belt, the carbonate strata of the Tianxingzhou karst belt and Zhuankou karst belt are directly covered with old clay or red bed, and there are old clay or red bed under the Holocene strata in some areas, so the risk of disasters is lower, which can be confirmed in the statistics of PS point displacement rate in Figure 8. The Tianxingzhou karst belt and Zhuankou karst belt had subsidence in a small area where was located at the Houhu and the point BSZ1 of Baishazhou in Figure 7. The areas of land subsidence were often accompanied by the progress of construction from the analysis of results of the four main subsidence centers also shows us the influence of construction on land subsidence.
The impact of changes in the water level of the Yangtze River can lead to a slight variation of displacements along its coastal areas rather than previous thoughts . The water level can only cause a transient uplift of approximately 2 mm after the water level of the Yangtze River reaches its peak, which suggests the change in the water level of the Yangtze River is not a key factor of land subsidence in the time period.
In addition, accurate validation of the IPTA result is a challenging issue in the study, because of the unavailability of synchronous leveling measurements, especially for continuous leveling observation. Furthermore, the leveling and IPTA measurements may mean different displacement features, because IPTA measurements present the change over a grid at the resolution of a pixel while leveling measurements are observed based on one site.

Conclusions
Land subsidence is caused by multiple factors. We obtained the time series and annual average subsidence rate of Wuhan City from June 2015 to October 2020 using 47 Sentinel-1A SAR satellite image data in Interferometric Wide Swath mode based on the inbuilt IPTA method of GAMMA software. The monitoring results show that the average annual subsidence rate in the main urban area of Wuhan can reach 41.0 mm/a, and accumulated subsidence of 219.4 mm during the study period. In the area where carbonate is directly covered with sandy soil was be found large areas of land subsidence. The largest subsidence center may have transferred to Baishazhou in the Wuchang district from Hankou compared with the previous study. It can be concluded that large-scale construction is a major factor leading to land subsidence in the urban area. In addition, long-term and frequent ground loading is also an alternative contribution to land subsidence in areas with weak ground bearing capacity. The construction of elevated stations has an insignificant contribution to land subsidence, but the construction of underground stations can lead to the subsidence of the surrounding ground. The change in the water level of the Yangtze River will cause a small periodic variation of the surrounding surface, while it does not change the overall trend of the surface displacements. Continuous monitoring of surface displacements inferred from InSAR technology is a reliable method for alleviation of potential destructive hazards, especially for the period of construction of land building and subway. preparation, W.L.; writing-review and editing, Q.S., H.W., C. K.S., L.J., B.Y.; supervision, Q.S., All authors have read and agreed to the published version of the manuscript.

Data availability statement
The Sentinel 1A data can be obtained freely from ESA. For detail please to visit Open Access Hub (https://scihub.coper nicus.eu/dhus/).