Pre- and Postcollapse Ground Deformation Revealed by SAR Interferometry: A Case Study of Foshan (China) Ground Collapse

On the evening of 7 February 2018, a deadly collapse of a metro tunnel under construction in the Southern China city of Foshan caused 11 deaths, 8 injuries, and 1 missing person. For disaster prevention and mitigation, the spatiotemporal ground deformations before and after the collapse event were derived from 55 Sentinel-1A synthetic aperture radar (SAR) images spanning from March 2017 to January 2019. The results showed that prominent ground subsidence in the shape of a funnel with a maximum rate of 42mm/year was observed in the vicinity of the collapse area before the accident. After the accident, the area and magnitude of subsidence decreased compared with precollapse subsidence. This decrease is related to the progress of tunnel excavation and groundwater changes. In the temporal domain, continuous subsidence was observed over a year before and after the accident, and accelerated subsidence appeared one month before the collapse accident. Soft soil consolidation and tunnel-induced soil losses were the main reasons for the subsidence over the study area. The leakage of groundwater accounted for the collapse event. The leaked groundwater eroded the soil, resulting in the formation of an arched hole. The connection between the arched hole and the tunnel reduced the bearing capacity of the soil layer above the arched hole, triggering the collapse event. The findings provide scientific evidence for future collapse monitoring and early warning due to tunnel excavation.


Introduction
As an advanced space observation technique, synthetic aperture radar (SAR) interferometry has demonstrated advantages for monitoring the ground deformation for collapse sinkholes compared with some ground-based instruments and techniques, such as high-precision leveling [1], robotic total stations [2], global positioning systems [3], and terrestrial laser scanner [4]. The sinkholes along the Dead Sea coast in Israel and Jordan have been successfully identified and monitored by SAR interferometry since 2002, providing rich references and experiences for studying the deformation of collapse sinkholes [5][6][7][8][9][10][11]. The continuous subsidence of sinkholes in Wink, Texas, USA, was observed from timeseries SAR interferometry, suggesting that continuous monitoring of subsidence in the vicinity of sinkholes could help prevent and mitigate catastrophic outcomes [12][13][14]. A dramatic increase in vertical deformation in the last few years preceding the sinkhole formation at a shopping mall in Heer-len, the Netherlands, was recorded in SAR data between 1992 and 2011, showing the feasibility of using satellite radar interferometry to detect a migrating cavity [15]. Similar techniques were applied to the collapse sinkhole in the Ebro Valley, Spain [16]; South Africa [17]; Poland [18]; and Arizona, USA [19].
On the evening of 7 February 2018, a deadly ground collapse of a metro tunnel under construction occurred in Foshan, Guangdong Province, China. The accident caused 11 deaths, 8 injuries, 1 missing person, and a direct economic loss of more than 53 million yuan (USA $7.8 million) [20]. The metro construction is an east-west line including 12 underground stations and 5 elevated stations. The project has 5 shield tunnel intervals, and some of the intervals are constructed in strata that feature rich underground water and complex geological conditions. The collapse area is part of the shield tunnel interval between Huchong Station and Lvdaohu Station, which is 1804 m long and 17.2-34.4 m deep [20]. After the collapse event, a sinkhole 81 m long, 65 m wide, and 8 m deep was formed, as shown in Figure 1. Since the collapse area was under metro construction, the event was attributed to the underground tunnel excavation. The results of a field investigation indicated that the accident was due to the groundwater leakage, which let the water pour into the tunnel and triggered the collapse [21]. Based on this investigation, the leakage of groundwater was the crucial factor causing the collapse event. Accordingly, several researchers analyzed the reasons for the groundwater leakage, such as complex geological strata [21], surface loading [22], and the leakage at the tail of the shield machine [23]. For the ground deformation, existing research focused on the whole city and there was no information available for deformation monitoring and analysis regarding to the collapse area [22,24,25].
To better understand the Foshan collapse event, multitemporal SAR interferometry was used in this study to obtain the spatiotemporally evolution of the ground deformation in the vicinity of the collapse sinkhole. For this, 55 C-band Sentinel-1A images were acquired in an ascending orbit. The annual deformation velocity was estimated to characterize the spatial pattern of deformation over the collapse area. Subsequently, the time-series deformation was interpreted in detail. Finally, the cause of the accident was discussed based on the observed deformation.

Study Area and Data Collection
2.1. Geological Setting of the Study Area. The study area is located at Foshan, which lies in the northern part of the Pearl River, in Southern China. Crossed by the Tanzhou, Dongping, and Jili Rivers, the collapse area is a large alluvial delta, where it is mostly flat with few hills and has an altitude from 1.3 to 4.6 m above sea level [21]. Belonging to the subtropical maritime monsoon climate, the study area is characterized by short spring and autumn, long summer, warm winter, abundant sunshine, and plentiful rainfall [26].
Quaternary strata, with unconsolidated sediments from 0 to 60 m thick, are widespread in the study area. Figure 2(a) shows that the thickness of quaternary strata decreases from west to east. Figure 2(b), displaying the vertical stratigraphic structure around the constructed tunnel, shows that the stratum from top to bottom is artificial filling soil, silty clay, mucky soil, powdered fine sand, medium sand, gravel, sandstone, silty mudstone, and siltstone [21]. The artificial filling soil with a thickness of 0-5 m is distributed along the road and mainly from the weathered residual soil, weathered rock, and gravel. Below the artificial filling soil, it is distributed by the silty clay characterized by high water content and high compressibility. The experiment indicates that the water content is about 10%-48% and compression modulus is about 3-13 MPa for silty clay over the study area [26]. The mucky soil, distributed below the silty clay, is mainly from the sea alluvial silt and mucky soil, which has properties of high water content, high compressibility, high void ratio, and low strength. For the study area, the mucky soil has the water content of 24.3%-68.8%, compression modulus of 1.6-2.4 MPa, and void ratio of 0.8-2.0 [26]. Powdered fine sand and medium sand belong to the sandy soil, which is characterized by high permeability and distributed below the mucky soil. Gravel, sandstone, silty mudstone, and siltstone are classified as sandstone soil, which is distributed at the bottom of quaternary and has the properties of high permeability, macropores, and high shear strength.
Groundwater of the study area is mainly sourced from quaternary pore water and bedrock fissure water, where the pore water is further divided into unconfined and confined water. The unconfined water mainly lies in the soft soil layer, while the confined water is mainly found in the sandy and sandstone soil layers [21]. The changes of groundwater are affected by the precipitation and anthropogenic activities [23]. The annual average precipitation is about 1,690 mm, of which 50% concentrates in summer, 40% concentrates in spring and autumn, and 10% concentrates in winter [21]. Therefore, changes of groundwater due to the precipitation present the long-term periodical variations: high groundwater level in summer, middle groundwater level in spring and autumn, and low groundwater level in winter. Changes of groundwater due to the anthropogenic activities show the short-term fluctuations with domestic and industrial    3 Journal of Sensors water consumptions, as well as underground construction. The domestic and industrial water consumptions are not dominated over the study area due to the rich surface water [26]. The operations of pumping and discharging groundwater during the underground construction, such as excavations of a metro tunnel, may cause the significant groundwater changes.

Data Collection.
A total of 55 C-band Sentinel-1A images were collected to derive the ground deformation over the Foshan collapse area. Table 1 and Figure 3 provide the parameters and coverage of these SAR images, respectively. The Sentinel-1A images in the ascending orbit were acquired between March 2017 and January 2019, which allowed us to obtain the ground deformation over one year before and after the formation of the collapse sinkhole. A TerraSAR-X addon for Digital Elevation Measurements (TanDEM-X) with a spatial resolution of 90 m was acquired as an external digital elevation model (DEM) to remove the topographic phase from the differential interferograms.

Methodology
Using the collected 55 C-band Sentinel-1A images, the annual deformation velocity and time-series deformation were calculated using the multitemporal SAR interferometry technique. During this procedure, we followed the methods of Chen et al. [27]. Interferometric pairs were generated by setting small temporal and spatial baselines. On the basis of the experiment, a spatial baseline below 100 m and a temporal baseline less than 30 days were designed to generate the interferograms. A total of 96 interferograms were produced, among which 42 interferograms were used for precollapse deformation estimation and 54 interferograms were used for postcollapse deformation estimation, as shown in Figure 4. The topographic and orbital phases were simulated   Journal of Sensors and subsequently removed from the differential interferometric phase using the collected TanDEM-X DEM and SAR orbital data [27]. Then, adaptive spectral filtering with a window size of 32 was executed to suppress the interferometric noise [28]. After that, phase unwrapping using the minimum cost flow (MCF) method was conducted to retrieve the absolute phase [29]. Finally, the atmospheric delay errors were reduced through the cascade of a lowpass filtering implemented in the two-dimensional spatial domain, followed by a temporal high-pass filtering [30]. After these processes, the unwrapped interferograms were used to estimate the annual deformation velocity and timeseries deformation before and after the formation of the collapse sinkhole. The annual deformation velocity, which allowed the characterization of the spatial pattern of ground deformation over the study area, was estimated by the method in [31]: where Ph rate is the annual deformation velocity, N is the number of interferograms involved in the estimation, Δt j is the time interval between the master and slave pair, and φ j is the unwrapped interferograms. The standard deviation of Ph rate was estimated by the method in [32]: where λ is the radar wavelength and VarðPh rate Þ is the standard deviation of Ph rate . VarðPh rate Þ was used to validate our observed annual deformation velocity due to lack of additional observations, such as ground-based leveling measurements.
The time-series deformation demonstrates the temporal variations in ground deformation and can be estimated by the method in [33]: where matrix A consists of time intervals between consecutive SAR acquisition, θ is the incidence angle, V u represents the unknown vertical velocities that are to be determined, and φ represents the unwrapped interferometric phase. The unknown parameters V u for each pixel are solved by applying the singular value decomposition (SVD) [34]. Once V u was estimated, the time-series deformation d = was reconstructed from the computed deformation rates by numerical integration, where M is the number of SAR images involved in calculation.

Precollapse Ground Deformation.
The annual deformation velocity was estimated by Equation (1) using 26 Sentinel-1A images acquired between 12 March 2017 and 30 January 2018, which revealed the spatial pattern of ground deformation over the year prior to the collapse accident, as shown in Figure 5(a). All measurement values are relative to a common reference point that was considered to be stable over the entire time period, and positive and negative values represent uplift and subsidence, respectively. The prominently subsiding region was detected in the vicinity of the collapse sinkholes, as indicated by the black outline polygon in Figure 5(a). This subsiding region was further divided into three subregions. The first subregion contained the collapse Pre-collapse interferograms P -o a n r e o r m Post-collapse interferograms -o a s t r g am Figure 4: The distribution of the temporal and spatial baselines before and after the collapse event.
5 Journal of Sensors region and showed the maximum subsiding rate of 42 mm/year. The second subregion, with a maximum subsiding rate of 38 mm/year, was located at about 600 m north of the collapse sinkhole. The third subregion was located at about 900 m northwest of the collapse sinkhole, where the maximum subsiding rate was up to 34 mm/year. To evaluate the observed result, the error distribution of deformation velocity was estimated using Equation (2), as shown in Figure 5(b). It was observed that the errors were below 3 mm for most regions and up to 5 mm for the severe subsiding area, suggesting the reliability of our observed subsidence.
Based on this observation, it was summarized that prominent ground subsidence with a maximum rate of 42 mm/year appeared in the vicinity of the sinkhole before the collapse.
To investigate the temporal variation of precollapse deformation, the time-series deformation was estimated using multitemporal SAR interferometry processing. Figure 6 shows the time-series deformation from 12 March 2017 to 30 January 2018, which was derived from 26 Sentinel-1A images acquired before the formation of the collapse sinkhole. Continuous subsidence was observed in the vicinity of the collapse, where the cumulative maximum  Point P2, which is about 900 m away from the center of the collapse sinkhole, showed similar variations in time-series deformation as point P1. The cumulative subsidence was about 28 mm for point 2, which is less than at point P1. As for the deformation in the month prior to the collapse event on 7 February 2018, accelerated subsidence was observed at both points: the deformation decreased by 7 mm and 8 mm between 6 January and 30 January 2018 for P1 and P2, respectively. In summary, continuous subsidence occurred in the vicinity of the collapse before the accident, and accelerated subsidence was observed when metro construction started and one month before the formation of the sinkhole.

Postcollapse
Ground Deformation. Figure 8 shows the annual deformation velocity and its error distribution after the formation of the collapse sinkhole, which was derived from 29 Sentinel-1A images acquired between 11 February 2018 and 25 January 2019. Similar to the precollapse deformation velocity map in Figure 5(a), the subsiding region occurred in the vicinity of the collapse, as shown by the black outline in Figure 8(a). However, the area and magnitude of subsidence are smaller compared with Figure 5(a). The statistics showed that the maximum subsiding rates were about 30, 31, and 32 mm/year for the first, second, and third subsiding subregions, respectively. Meanwhile, it was observed that the subsidence started to extend to the east, where subsidence area and magnitude increased. The lower subsiding rates were due to the completion of the metro construction over the study area, and the construction continuing to the east extended the subsidence. This is similar to the precollapse deformation where the deformation velocity errors were below 3 mm for most regions and up to 5 mm for the severely subsiding area, as shown in Figure 8(b). Based on this observation, we found decreases in the area and magnitude of subsidence after the collapse, which may be related with the groundwater level changes. Figure 9 shows the time-series deformation from 11 February 2018 to 25 January 2019, which was derived from 29 Sentinel-1A images after the collapse of the sinkhole. Like the precollapse time-series deformation in Figure 6, continuous subsidence was observed in the vicinity of the collapse sinkholes, where the cumulative maximum subsidence reached 30 mm. However, the area and magnitude of subsidence are smaller compared with Figure 6. Figure 10 displays the time-series deformation at P1 (blue rectangle) and P2 (red circle), both of which are marked in the last image of Figure 9. Point P1 experienced a cumulative subsidence of 20 mm between 11 February 2018 and 25 January 2019, which is a decrease of 15 mm compared with the precollapse subsidence. Although the wave-shaped deformation variation was retained, the amplitude of variation was lower when compared with Figure 7(a). Point 2 experienced a cumulative subsidence of 28 mm, which was the same as the precollapse subsidence. In summary, continuous but decelerated subsidence was observed after the incident, and accelerated subsidence was observed over 48 days after the formation of the sinkhole.

Discussion
5.1. Subsidence due to Soft Soil Consolidation. Located on a large alluvial delta, unconsolidated sediments are widespread in the study area (Figure 2(a)). The stratum from top to bottom is artificial filling soil, silty clay, mucky soil, powdered fine sand, medium sand, gravel, sandstone, silty mudstone, and siltstone (Figure 2(b)). Considering the presented characteristics of high water content and high compressibility, the upper silty clay and mucky soil are classified as soft soil, which has many properties that are unfavorable for their use in projects [26]. The lower soil layers, including powdered fine sand, medium sand, gravel, sandstone, silty mudstone, and siltstone, are classified as sandy and sandstone   Journal of Sensors soils, which are characterized by high permeability, low compressibility, and high shear strength. In this context, pore water in soft soils is easily discharged along the sandy and sandstone soils in the presence of natural and human loadings, resulting in the decrease of pore fluid volume, compression of soft soils, and subsequent ground subsidence. The relevant materials show that the magnitude of subsidence is related to the thickness of soft soils [27,34]. Thus, the relationship between ground deformation and soft soil thickness was investigated and analyzed over the study area. Figure 11(a) displays the thickness of soft soils over the study area, where the thickness between Huchong Station and Lvdaohu Station is from 20 to 30 m. Figure 11(b) shows the relationship between ground deformation and soft soil thickness. The green circle in Figure 11(b) represents the ground deformation between Huchong Station and Lvdaohu Station, which is extracted from Figure 5(a). The red and blue lines are, respectively, the fitting deformation and soft soil thickness between Huchong Station and Lvdaohu Station. The comparison results show the approximate consistency between the ground deformation and the soft soil thickness in the spatial distribution; the greater the subsidence, the   Journal of Sensors thicker the soft soils, particularly in the vicinity of the collapse sinkhole, where the maximum subsidence and soft soil thickness are observed. The literature indicates that the subsidence is due to the soft soil consolidation mainly from the surface loadings [22] and groundwater exploitation [25]. With respect to the collapse area, surface loadings come from the vehicle loadings and some buildings and infrastructure. Groundwater exploitation to meet the needs of domestic and industrial water consumption was rare over the study area due to the rich surface water [26]. However, this process occurs during the tunnel excavation, where the operations of pumping and discharging groundwater are needed to reduce the external water pressure [35].

Subsidence due to Soil Losses.
The studies indicate that the soil losses during the tunnel excavation may cause the subsidence, which is located at the center of the tunnel [36,37]. This subsidence occurs when a difference exists in the volumes of the excavated soil and tunnel. In this case, soil losses contributed to the subsidence according to the adopted shield tunnel excavating method and the local geological conditions [37]. Figure 12 shows the annual deformation velocity along the planned metro line. It is observed that prominent subsidence appears along the completed metro line, which is caused by the tunnel-induced soft soil consolidation and soil losses. However, the deformation is not significant for the region where the tunnel excavation has not started. This observation suggests the effect of tunnel excavation on ground subsidence.
In order to further show the ground subsidence induced by soil loss, Peck model, an empirical model to calculate the deformation, was used in this study [36]. A schematic diagram of the Peck model is shown in Figure 13, where point O is the center of the tunnel, R is the tunnel radius, H is the depth of the tunnel centerline, S max is the maximum subsidence, SðxÞ is the deformation due to soil loss, and i is the width coefficient of the subsiding trough, which can be empirically estimated as a function of the depth of the tunnel centerline H, and the horizontal direction of the tunnel cross section is the x axis. The model is expressed as [36] where V s and V l are the soil loss and soil loss rate, respectively, both of which are difficult to accurately estimate since they are related to many factors, such as stratigraphic condition, hydrogeology, geotechnical parameters, rock/soil porosity, water content, granulometry, construction technique, and excavating method. In this study, R and H in Equations (4), (5), and (6) were directly derived from [21] and parameters V l = 3% and i = 0:76 * H were derived from [37]. Based on these parameters, the ground subsidence due to soil loss was estimated, as shown in Figure 14. It was found that the subsidence due to soil losses was in the shape of a funnel, with the subsiding center was located at the center of the tunnel. Compared with the deformation map in Figure 6, the three detected subsiding subregions were approximately shaped as funnels, suggesting consistency in the deformation shape. For the subsiding center, the first subsiding subregion showed a center along the metro line. However, the other two subregions did not overlap the metro line. This phenomenon was due to the reduction in the subsiding area after the tunnel excavation, which was identified in the time-series deformation map. In terms of magnitude of subsidence, the modeled maximum subsidence due to soil losses in this case was about 25 mm. The observed subsidence due to soil losses at point P1 in Figure 7 was about 28 mm when considering 80% of total subsidence from soil loss [37]. The modeled and observed subsidence at point P1 is comparable in magnitude. Based on the discussion, it was summarized that the observed subsidence in this case was due to soft soil consolidation and soil losses. Subsidence due to soft soil consolidation occurred before and after tunnel excavation, whereas subsidence due to soil losses occurred during tunnel excavation.

Collapse due to Groundwater
Leakage. The result from an investigation report showed that the Foshan collapse occurred due to leakage of groundwater above the tunnel [21]. The groundwater poured into the tunnel, causing the tunnel structure damage and triggering the collapse. In view of filed investigation, the damage of underground pipelines is considered to be related with the leakage of groundwater.
Here, three aspects are shown to support this point. The first is the aging of underground pipelines over the study area. Like most other cities, there are more than 30-year-old underground pipelines in Foshan region [38]. This aging of underground pipelines is easily damaged in the case of ever-increasing urbanization and the continuous construction, development, and expansion of urban areas, such as 11 Journal of Sensors the abrupt collapse due to a 44-year-old sewer pipeline in Fraser, USA [39]. The second is the observed inhomogeneous ground subsidence over the study area. The filed investigation suggests that some infrastructures have been affected by ground subsidence in the vicinity of the collapse area, such as ground crack, wall crack, and pipeline rupture, as shown in Figures 15(a)-15(d). Although there is no direct evidence to show this effect on the underground pipelines, plenty of research has proved that the inhomogeneous ground subsidence can cause the underground pipeline leakage [40][41][42].    Journal of Sensors The leakage in underground pipelines mainly refers to the leakages, bursts, or blockages in sewer, drain, and/or water pipelines [43,44]. The third is the similar collapse event in the vicinity of the collapse area. On November 5, 2017, a small ground collapse event was observed near the station of Huchong, which caused the formation of a collapse sinkhole with an area of 60 m 2 [2]. It is confirmed that the collapse is caused by the underground pipeline leakage, as shown in Figures 15(e) and 15(f). Based on these three aspects, it is the damage of underground pipelines causing the leakage of groundwater. Therefore, the detection of underground pipeline leakage should be conducted to prevent the collapse events, particularly for subsiding areas. Related studies investigated the mechanisms and necessary conditions for collapse sinkhole formation due to the leakage of groundwater [6,16,[45][46][47]. Based on these studies, the widely accepted process of collapse sinkhole formation is shown in Figure 16. Firstly, the leakage of groundwater erodes the soil, forming an arched hole. This arched hole is generally related to the size and location of the later formed collapse sinkhole and is therefore prerequisite to the formation of a collapse sinkhole [48,49]. Anthropogenic influences and the presence of urban infrastructure can significantly expedite the process of arched hole formation [50,51]. With the excavation of the tunnel, the arched hole and the tunnel are connected. This connection reduces the bearing capacity of the soil layer above the arched hole; hence, a ground collapse sinkhole forms. Thus, the collapse event occurred in this region due to groundwater leakage and arched hole formation.

Conclusions
On the evening of 7 February 2018, a deadly collapse of a metro tunnel under construction in the Southern China city of Foshan caused 11 deaths, 8 injuries, and 1 missing person. To better understand this event, the ground deformations before and after the collapse was derived using a multitemporal InSAR technique. The driving factors for the ground deformation and the collapse were discussed based on the investigation report, local geological conditions, and observed deformation. The findings provide scientific evidence for future collapse monitoring and early warning due to tunnel excavation. Based on this study, the following conclusions were made: (1) Precursory subsidence was observed in the vicinity of collapse sinkhole. Using 55 C-band Sentinel-1A images, the annual deformation velocity and timeseries deformation maps before and after the event were retrieved through multitemporal SAR interferometry processing. The results showed that prominent ground subsidence with a maximum rate of 42 mm/year occurred in the vicinity of the collapse before the formation of the sinkhole. Compared with precollapse deformation, the area and magnitude of subsidence decreased during postcollapse deformation. In the temporal domain, accelerated subsidence was observed during metro construction and one month before the formation of the collapsed sinkhole. These observations indicated the presence of precursory subsidence for collapse sinkhole, which allowed us to forecast the collapses and thus reduce the risk for people and property (2) Soft soil consolidation and soil losses were the main factors affecting the subsidence over the study area. The quaternary stratum from top to bottom was soft, sandy, and sandstone soils. In this context, the ground subsidence easily occurred due to soft soil consolidation and tunnel excavation. The timeseries deformation maps showed the continuous subsidence before and after the metro construction. This subsidence was mainly from the soft soil consolidation due to surface loadings and groundwater exploitation. The modeled and observed subsidence is comparable in magnitude, suggesting the subsidence due to soil losses during tunnel excavation (3) The leakage of groundwater accounted for the collapse event. The result from an investigation report showed that the collapse occurred due to the leakage of groundwater above the tunnel. Considering the  14 Journal of Sensors mechanisms and necessary conditions for collapse sinkhole formation, the process of collapse was retrieved. The leaked groundwater eroded the soil, forming an arched hole. The connection between the arched hole and the tunnel reduced the bearing capacity of the soil layer above the arched hole, thereby triggering the collapse event These findings contribute to the understanding of the spatial-temporal evolution of ground deformation before and after the collapse event. However, clearly describing the reasons for collapse is challenging due to a lack of related data, such as continuous ground-based deformation monitoring data. Therefore, sufficient in situ measurements will be collected in the future to analyze the reasons for the collapse.

Data Availability
Sentinel-1A/B images were from the European Space Agency (ESA) project, and these data were downloaded from the Sentinel-1 Scientific Data Hub. ALOS-2 images were provided by the Japan Aerospace Exploration Agency (JAXA).

Conflicts of Interest
The authors declare no conflict of interest.