# On Time-Series InSAR by SA-SVR Algorithm: Prediction and Analysis of Mining Subsidence

Authors: Yun Shi; Qianwen Li; Xin Meng; Tongkang Zhang; Jingjian Shi
Journal: Journal of Sensors
 (2020)

---

## Abstract

Given the increasingly serious geological disasters caused by underground mining in the Hancheng mining area in China and the existing problems with mining subsidence prediction models, this article uses the small baseline subset interferometric synthetic aperture radar (SBAS-InSAR) technology to process 109 Sentinel-1A images of this mining area from December 2015 to February 2020. The results show that there are three subsidences: one in Donganshang, one in south of Zhuyuan village, and one in Shandizhaizi village. In the basin, the maximum annual average subsidence rate is 300 mm/a, and the maximum cumulative subsidence is 1000 mm. The SBAS-InSAR results are compared with Global Positioning System (GPS) observation results, and the correlation coefficient is 74%. Finally, a simulated annealing (SA) algorithm is used to estimate the optimal parameters of a support vector regression (SVR) prediction model, which is applied for mining subsidence prediction. The prediction results are compared with the results of SVR and the GM (1, 1). The minimum value of the coefficient of determination for prediction with SA-SVR model is 0.57, which is significantly better than that those of the other two prediction methods. The results indicate that the proposed prediction model offers high subsidence prediction accuracy and fully meets the requirements of engineering applications.


---

## Body

## 1. Introduction

Ground subsidence in a mining area is a type of vertical deformation of the ground, which is prone to slow regional changes due to the destruction of the structure of the rock mass caused by mining. Long-term underground mining results in movement and deformation of the overlying rock layer and the ground surface and can cause continuous or discontinuous sinking, tilting, curvature, stretching, compression, etc., of buildings and other structures located within the mining area. Deformation, cracks, collapses, collapse pits, landslides (movement), and similar types of damage degrade the land resources and ecological environment in a mining area, significantly impacting the productivity and daily lives of residents and also affect future engineering construction in the area [1]. In a word, the ground subsidence caused by mining is a very complicated process. To effectively control the ground subsidence in mining areas caused by mining and to reduce losses, it is necessary to accurately understand the causes of mining subsidence and its development process over time.Traditional levelling and Global Positioning System (GPS) measurements have the disadvantages of sparse measurement points, a small range, and a high cost for monitoring surface subsidence. As an alternative, differential interferometric synthetic aperture radar (D-InSAR) has been proven to be able to probe small surface deformations [2]; it is widely used in the monitoring of landslides [3], urban areas [4], mining area ground subsidence [5], and other geological disasters, but it is susceptible to limitations such as decoherence and atmospheric delay, and it cannot provide continuous time-series deformation information. As a sequential InSAR technology, small baseline subset InSAR (SBAS-InSAR) mitigates the shortcomings of D-InSAR technology by adopting a short time-space baseline [6–8]. Zhang et al. [9, 10], Hu et al. [11], and Jiang et al. [12] have used time-series InSAR technology to investigate land subsidence in Los Angeles, California, USA; Ningbo, China; and the Lost Hills Oil Field in California. Yin et al. [13], Yang et al. [14], Saygin et al. [15], and Mark et al. [16] have applied time-series InSAR technology for mining area monitoring in Lengshuijiang, China; Datong, China; Zonguldak Province, Turkey; and Springfield, Illinois, USA, and have shown the applicability of SBAS-InSAR technology in mining surface subsidence monitoring applications. He et al. [17] used ALOS-1 data to identify large-scale surface deformations near China’s Hancheng coal mine. Through a literature review, it can be found that most applications of InSAR technology in mining subsidence monitoring have focused on data monitoring and analysis, and the research time span has been relatively short. However, for mining subsidence monitoring, it is necessary to combine monitoring and prediction methods to establish a complete mining subsidence prediction model with long time-series capabilities.The current mining subsidence prediction methods can be divided into four categories: influence function methods, empirical methods, theoretical methods, and other methods [18, 19]. The most commonly used influence function methods are the probability integral method in China and the Budryk-Knothe method elsewhere. The formulas for these two methods are the same. By adjusting the prediction parameters, this method can be applied for prediction under different geological and mining conditions, but the expected values of the parameters must be known and are difficult to obtain. Empirical methods include method based on typical curves, profile functions, and other approaches. A large amount of measured data is used to determine the laws governing the observed deformation. The prediction formulas are simple, but this approach is suitable only for areas with similar geological mining conditions, so its applicability is poor. Theoretical methods rely on complex mechanical parameters and calculations. Previous studies have shown that such methods are very difficult to apply. With the rapid development of computer technology, a large number of complex calculations have become possible, but it is still difficult to obtain mechanical parameters on a large scale. Other methods include neural network backpropagation (BP) [20], grey system theory (GM) [21], and support vector regression (SVR) [22, 23]. Researchers use computer algorithms to establish prediction models to predict mining subsidence in mining areas. Various new ideas are emerging in the field of mining subsidence prediction, but they are still in the exploration stage. Although the BP algorithm considers the nonlinearity of ground subsidence, when few training samples are available, the prediction results obtained are unreliable. The GM algorithm transforms messy original data into regular time-series data in a prescribed way. It is suitable only for medium- and short-term forecasts, and some data fitting results are far from the original data. In this paper, SBAS-InSAR technology is combined with the SA-SVR algorithm to predict mining subsidence for the first time. The SA-SVR algorithm combines SA and SVR. The SVR algorithm can handle any nonlinear situation, has a strong generalization ability and good theoretical support, and has been widely used in research on the regression fitting of nonlinear functions [24].This article is organized into five parts. The first part summarizes the current research progress and shortcomings of mining subsidence monitoring and prediction. The second part briefly describes the basic idea of the combined SA-SVR algorithm. The third part presents the extraction of the 4-year average annual deformation rate and the time-series cumulative deformation of a mining area in Hancheng based on SBAS-InSAR technology, analyses the overall deformation behaviour in the mining area before and after mining, provides a decision basis for the subsequent mining management of the mining area, and verifies the SBAS-InSAR monitoring results against GPS data. The fourth part reports the prediction results achieved by using the SBAS-InSAR monitoring results as the training and test sets for the SA-SVR algorithm and finally compares the SA-SVR prediction results with those of SVR and GM. The fifth part summarizes the whole paper and discusses the advantages and disadvantages of the SA-SVR algorithm.


## 2. Basic Idea of the SA-SVR Algorithm

The SVR algorithm performs parameter regression estimation on the basis of a support vector machine (SVM). This algorithm is suitable for low-information, linear and nonlinear models and for the recognition of high-dimensional spatial patterns from small samples. Research has shown that the SVR algorithm enables stable time-series analysis and statistical forecasting [25]. The basic principle of the method is to find the optimal division plane in a high-dimensional space and to minimize the sum of the energy of the regression coefficients under the constraint that the error of the given sample data with respect to the division boundary is within a certain value. The original SVR algorithm can also be used as a regularization method to balance the error when fitting sample data and the energy of the regression coefficients.Time-series InSAR technology is used to obtain time-series data on mining area settlement, denoted byx1,y1,x2,y2⋯xn,yn, xi∈X=Rn,yi∈Y=R. To study the time-space relationship of settlement in a mining area and use it to make settlement predictions, we need to analyse xi,yii=1,2,⋯,n quantitatively. More specifically, we hope to find a suitable mapping Φ to describe the relationship between xi and yi and use SVR for data learning and prediction. Under the assumption that the mapping Φ between x,y follows an approximately linear relationship, we use the following formula:(1)y=ω,x+bwhere·,· represents the vector dot product, x is the value of the input variable, and y is the corresponding output value. According to Smola and Scholkopf [26], the SVR algorithm can be applied to solve the following optimization problem:(2)min:12ω2+C∑i=1nξi+ξi∗s.t.yi−ω,xi−b≤εi+ξiω,xi+b−yi≤εi+ξi∗ξi,ξi∗≥0where C is a penalty or equilibrium factor and ξi,ξi∗ are the upper and lower bounds of the regression line, respectively, where an excess violation of εi is allowed.The parameters of the SVR model include the accuracyεi describing the regression model and the data column, the penalty factor C, and the mapping Φ. The kernel function Kx,x′=Φx,Φx′ can be obtained from this model. These parameters and the kernel function together determine the overall performance of the prediction model. From the constraint conditions of the optimization model given in (2), it can be seen that when the difference between the theoretical value and the actual observed settlement value is less than εi, SVR will not penalize such a violation; that is, the corresponding values of ξi,ξi∗ are zero. Therefore, SVR is also called ε-insensitive regression parameter estimation. Generally, if the observed settlement values are the same in terms of accuracy, all εi have the same value ε. If the accuracy value εi is too small, that is, the accuracy of the data column is too high, then under-fitting or over-fitting can easily occur. Smola and Scholkopf [26] defined the allowable kernel function and presented a detailed discussion of the structure of the allowable kernel function problem. For details, please see Smola and Scholkopf [26].In the literature on SVR, to solve for the optimal parameters of the optimization problem defined in (2), the original problem is generally converted into the corresponding dual problem using the Lagrange multiplier method, namely,(3)min:12α−α∗TXTXα−α∗+∑i=1nεiαi+αi∗−∑i=1nyiαi−αi∗s.t.α−α∗,ε=0αi,αi∗∈0,CIn SVMs, sequential algorithms are usually used to solve for the optimal coefficientsα^ and α∧∗ of such a dual problem (see [26]). α^ and α∧∗ essentially correspond to the Lagrangian coefficients of the first and second constraints, respectively, of the optimization problem given in (2). It can be proven (Smola and Scholkopf [26]) that the estimate ω^ of the unknown parameter in the original regression model defined in (1) is uniquely determined by the Lagrangian multiplication coefficients α^ and α∧∗:(4)ω^=Xα^−α∧∗The unknown parameterb is then given by the Karush-Kuhn-Tucker condition [27, 28]. Thus, the optimal prediction model is(5)y^=∑i=1nα^ι−α^ι∗xi,x+b^In fact, by comparing the original optimization problem (2) and its corresponding dual problem (3), it can be seen that the original optimization problem (2) has more parameters (ω and b) than the dual problem (3). Nevertheless, if the number of unknown parameters in the original regression model (1) is low, then the number of calculations will not be very different between the original optimization problem (2) and its dual problem (3).In this paper, an algorithm that combines SA with SVR is adopted to achieve the best prediction effect for mining subsidence. The SA algorithm was proposed by Metropolis et al. in 1953 [27] and started to be widely used after the publication of an article in Science by Kirkpatrick et al. in 1983. In the SA algorithm, random factors are introduced into the search process. In the iterative update process, a solution that is worse than the current solution may be accepted with a certain small probability, allowing the algorithm to jump out of locally optimal solutions to reach a better optimal solution. In this paper, an SA optimization algorithm is used in combination with SVR to estimate the parameters of the time-series forecasting model.The basic principle of the SA algorithm is to imitate the ideal crystal generation process in physics. It is assumed that at a high temperature, solid matter will be in a liquid state, in which the order of the molecules is free and random. Through a careful cooling process, the molecules in this random state will be arranged in an orderly manner to generate crystals in the ideal state. To use an optimized mathematical language to describe the process of forming an ideal crystal state through physical cooling, we need to formulate the ideal crystal state and the physical cooling process in terms of an optimized objective function and an optimization process, respectively. Specifically, the ideal crystal state is the minimum-energy state, that is, the state with the global optimal value of the objective function, and the physical cooling process is equivalent to the iterative process of optimizing the objective function, thus completing the mathematical simulation ideal crystal generation. In particular, Metropolis et al.’s original text from 1953 proposes random sampling from a uniform probability distribution as an iterative process for optimization. If the objective function is fully understood, random sampling from the corresponding probability distribution can instead be performed for optimization (see Peiliang Xu’s private message from February 2020) [29].The basic step of the optimization process in the SA algorithm is to compare the objective function value at a sampled point to the current optimal value. If the former is better than the latter, then the objective function value of the sampled point replaces the current optimal value, and the algorithm enters the next iteration. If the objective function value at the sampled point is worse than the current optimal value, to prevent the SA algorithm from converging to a locally optimal value too quickly, Metropolis and others have suggested that this worse solution should nevertheless be accepted with a certain probability during the optimization process. For an optimization functionfx, if the current solution is denoted by x0 and the search point for the next round is denoted byx′, Metropolis and others suggest using the following probability:(6)Px0⇒x′=1,fx′<fx0exp−fx′−fx0T,fx′≥fx0Accordingly, we decide whether to acceptx′ as the new solution. In formula (6), T represents the temperature variable, and ⇒ represents a state transition.The SA algorithm obtains the optimal solution through repeated iterations of the above process. However, ifT changes too quickly, then the SA algorithm will quickly converge to a locally optimal solution. To ensure that the algorithm will converge in a limited time, it is necessary to adjust the temperature parameter T. The steps of the adjustment process are as follows:(i)
The initial temperatureT (0) is set high enough that all states are acceptable(ii)
The annealing rate generally decreases exponentially(7)Tn=λTn,n=1,2,3⋯Here,λ takes a value in the range [0.8, 0.99] so that each temperature will have a certain probability to be tried.(iii)
A termination temperature is specified. Once the temperature reaches this preset threshold during the iterative process, annealing is completedThe flow chart of the SA-SVR algorithm is shown in Figure1.Figure 1
SA-SVR algorithm flow chart.


## 3. Time-Series InSAR Monitoring and Analysis

### 3.1. Research Background

The study area is a mine in the southwest of the Hancheng mining area, which is located on the southeastern edge of the Weibei Uplift in the Ordos Basin. The formation is generally a monoclinic structure inclined to the northwest. The geographical location is 110°17′–110°36′ east in longitude and 35°21′–35°51′ north in latitude. The mine is adjacent to the Dongze village structural belt and the Long-ting structural belt in the west and to the Huaxian fault zone in Hancheng and the Yellow River in the east. The geographical location of the mine is shown in Figure 2 [30, 31]. The main coal-bearing strata in the mining area belong to the Shanxi Formation and the Taiyuan Formation; coal seam No. 3 belongs to the Shanxi Formation, and coal seam No. 5 belongs to the Taiyuan Formation. The long-arm comprehensive mechanized coal mining method or the caving mining method has been adopted at all working faces, and the caving method is used to manage the roof. The mining thickness of each working face is 1.5–2.2 m, and the coal seam inclination is 1–3°, being nearly horizontal. This article mainly studies data taken from the 21308 working face of the No. 3 coal seam from December 2015 to February 2020. The strike length is 930 m, and the incline length is 210 m.Figure 2
Geographical location of the study area.An investigation performed in May 2017 showed that only the Donganshang group exhibited multiple cracks in houses and 3 cracks in the road pavement, with a crack length of approximately 3 m, a crack width of approximately 0.05–0.1 m, and an angle between the crack and the road surface of approximately 90°, which caused the ground to form a bulge of 0.1–0.2 m. There were 3 large-scale cracks in the ground, which were approximately 500 m long and 0.2–0.6 m wide, and 3 small-scale cracks, which were approximately 40~80 m long and 0.1–0.3 m wide. In the village, there were 7 landslides, with lengths of 8–15 m, widths of 5–12 m, and drops of 2–5 m.As research data, 109 scenes of Sentinel-1A radar satellite data in the study area were selected, which cover a time span from December 24, 2015, to February 25, 2020. The Sentinel-1A data are in the C-band, with a wavelength of 5.66 mm; the orbit revisit period is 12 days, the polarization mode is VV, the incidence angle is approximately 33.7°, the resolution of these data is5×20m, and the width of a single image is 250 km. The data used in this experiment are all products of the interferometric wide-width mode (IW) and are single-view complex data with slant distance information [32].The SRTM1 digital elevation model (DEM) was chosen as a source of external DEM data; its ground resolution is30×30m. For the precise orbital parameters, the precise orbit determination ephemeris data provided by the European Space Agency (ESA) (https://qc.sentinel1.eo.esa.int/) were adopted, which can be downloaded 21 days after acquisition of the Sentinel-1A data and have a positioning accuracy of within 5 mm [33]. External DEM and orbital data can be used to remove the levelling effect and orbital errors.


### 3.2. Data Processing

The SBAS-InSAR data process was realized based on the GAMMA Software. The main steps were as follows: extracting public bursts from the 109 Sentinel-1A scenes to improve the subsequent processing speed and efficiency, performing precise registration with an accuracy of 0.001 pixels, and setting a short time baseline (three consecutive scene images). The research area was cropped out based on the corresponding row and column numbers, and multiview interference processing was performed at a visual number ratio of 5 : 1 (range: azimuth). Differential interference, filtering (Goldstein), and unwinding (minimum cost flow method) were performed based on the initial baseline and the precise baseline. The high coherence coefficient method was used to extract points with high coherence, and elevation and atmospheric errors were removed through regression analysis of the high-coherence points. Finally, complete and accurate deformation rate and time-series cumulative shape variables were obtained [34–36].


### 3.3. SBAS-InSAR Results and Analysis

Following the above data processing method, the annual average subsidence rate (Figure3) and the time-series cumulative subsidence map (Figure 4) along the line-of-sight (LOS) direction were obtained for the mining area of Hancheng from December 2015 to February 2020. Figure 4 shows the timing accumulation of 11 scenes selected from among all 109 scenes at equal intervals. Figure 3 shows that the entire study area contains three relatively large subsidence areas: the Donganshang subsidence area, Zhuyuan village south subsidence area, and the Shandizhaizi subsidence area. A field investigation showed that there is good consistency between the subsidence area on the Donganshang and the distribution of the working face in the south No. 1 mining area in Hancheng. Both the subsidence area to the south of Zhuyuan village and the Shandizhaizi subsidence area have undergone coal seam mining, and mountain slippage, road cracks, and house cracks have occurred.Figure 3
Average subsidence rates from December 2015 to February 2020.Figure 4
Time-series cumulative subsidence from December 2015 to February 2020.The results presented in Figures3 and 4 show that the maximum annual average settlement rate in the study area is 300 mm/a and the maximum cumulative settlement is 1000 mm from December 2015 to February 2020. The area on the Donganshang presents the most serious subsidence, with a subsidence area of approximately 3.66 km2 and a cumulative subsidence of approximately 150–1000 mm. All of the working faces in this mining area show different degrees of subsidence. The 21506, 21307, and 21308 working faces stopped being mined in February 2017, and their subsidence rates are approximately 40–295 mm/a. The 21309 working face corresponds to the area with the most severe deformation during the study period. As of February 2020, the settlement rate of the 21309 mining face is the highest. The settlement area to the south of Zhuyuan village is approximately 2.56 km2, the settlement rate is approximately 40–257 mm/a, and the cumulative settlement is approximately 100~1000 mm. The Shandizhaizi settlement area is approximately 0.68 km2, the settlement rate is approximately 40–129 mm/a, and the cumulative settlement is approximately 100–570 mm.It can be seen from Figure4 that each mining face in the mining area shows a different degree of nonlinear settlement. Let the first image (December 24, 2015) be considered the initial image, in which the settlement amount is 0. The three settlement areas had formed by February 16, 2017. Over time and with the mining of the working faces, the subsidence area and the cumulative subsidence continued to increase. As of November 2, 2018, the subsidence area had gradually expanded westward, forming three obvious subsidence basins. As of February 25, 2020, the cumulative settlement had reached its maximum; the impact of the settlement had continued to expand from east to west, and the settlement area and cumulative settlement had gradually increased, indicating that the settlement area had not reached a stable state by that time.To further analyse the settlement information of the study area, the 21308 working face in the settlement area on the Donganshang is taken as the object of investigation, and profile lines A1A2 and B1B2 are established along the direction and inclination (the position of the section line is shown in Figure4). The working face trend is advancing along the direction from A2 to A1, and the working face tendency is advancing along the direction from B2 to B1. For image readability, the data from 10 scene images corresponding to some of the time-series cumulative settlement figures were selected to draw the corresponding strike and longitudinal time-series profile lines using the Origin software (Figure 5). Figure 5 shows that the main settlement areas are between points 11 and 41, while the maximum sinking position is always at point 26. The maximum accumulated settlement amounts to 748 mm and 880 mm. In general, over time, the settlement of each point along the direction and inclination of the working surface is gradually increasing, and an obvious settlement funnel is forming on the spatial scale, which is consistent with the surface mining subsidence in the mining area.Figure 5
Subsidence curves of the working face profiles: (a) A1A2 sinking curve; (b) B1B2 sinking curve.
(a)(b)Because this certain mining area in Hancheng has been repeatedly mined at multiple coal seams, to better analyse the pattern of the mining subsidence, four characteristic points (the positions of these characteristic points are shown in Figure3) on the Donganshang were selected for time-series cumulative settlement analysis, as shown in Figure 6. P1 and P2 are located at the edge of the deformation area on the Donganshang. They are affected by the mining of the 21504 and 21506 working faces, respectively. The cumulative settlement at points P1 and P2 from December 24, 2015, to February 25, 2020, amounted to 255 mm and 392 mm, respectively. P3 and P4 are located in the centre of the settlement area and are mainly affected by mining at the 21308 working face. The cumulative settlement at these points from December 24, 2015, to February 25, 2020, was relatively large, reaching 888 mm at P3 and 838 mm at P4. In general, the cumulative settlement at these four characteristic points has increased over time. As of February 25, 2020, the settlement had not stabilized.Figure 6
Time-series cumulative settlement at characteristic points.


### 3.4. Result Verification and Analysis

To quantitatively verify the ground subsidence monitoring results for the mining area as obtained via SBAS-InSAR technology, GPS observation points established in the subsidence area on the Donganshang were used to obtain ground displacement monitoring data, and the monitoring results were compared and analysed. GPS observation data collected from December 24, 2015, to January 12, 2017, at 8 observation points were selected. Because the surface deformation data obtained via SBAS-InSAR were obtained along the radar LOS direction, the three-dimensional deformation data obtained from the GPS observations were obtained in the vertical, east-west, and north-south directions. To compare the results of these two methods, three steps were required. First, it was necessary to project the vertical, east-west, and north-south deformations observed from the GPS data to the LOS direction. The projection formula is(8)LOS=cosθ−sinθcosαsinθsinαdvdEdNwhere θ is the local incidence angle of the satellite; ɑ is the azimuth of the satellite’s base; and dv, dE, and dN are the vertical, east-west, and north-south deformations, respectively, obtained from the GPS observations.The above formula was used to convert the GPS-observed deformations into deformations in the LOS direction. Then, the coordinates of the GPS observation points were converted into InSAR data coordinates through Gaussian projection to make the two sets of data consistent in their spatial dimensions. Finally, the extracted values were used via the Point tool of the ArcGIS software to extract the InSAR data deformation points corresponding to the GPS observation points in accordance with the neighbouring pixel values using bilinear interpolation to make the GPS results correspond to the InSAR results. Then, the two sets of results were compared and analysed. As shown in Figure7, the overall trends of the deformation values obtained from the SBAS-InSAR and GPS observations in the LOS direction are relatively consistent. Only the GPS results and SBAS-InSAR results of 8 observation points at point 4 are abnormal. Because these results were obtained through a combination of multiple SBAS-InSAR processing experiments, this abnormality may have been caused by human error transmission in the GPS measurements. The smallest deviation of 0.1 cm is observed at point 3. The average absolute error of the 8 points is 1.11 cm, which is within the acceptable range. A correlation analysis of the two sets of data shows that the correlation coefficient is 86%, indicating that they are highly correlated. On the whole, the two sets of results are basically consistent, indicating that SBAS-InSAR technology is reliable for the monitoring of ground subsidence in mining areas and can accurately reflect the corresponding information.Figure 7
Comparison of GPS and SBAS-InSAR results.


## 3.1. Research Background

The study area is a mine in the southwest of the Hancheng mining area, which is located on the southeastern edge of the Weibei Uplift in the Ordos Basin. The formation is generally a monoclinic structure inclined to the northwest. The geographical location is 110°17′–110°36′ east in longitude and 35°21′–35°51′ north in latitude. The mine is adjacent to the Dongze village structural belt and the Long-ting structural belt in the west and to the Huaxian fault zone in Hancheng and the Yellow River in the east. The geographical location of the mine is shown in Figure 2 [30, 31]. The main coal-bearing strata in the mining area belong to the Shanxi Formation and the Taiyuan Formation; coal seam No. 3 belongs to the Shanxi Formation, and coal seam No. 5 belongs to the Taiyuan Formation. The long-arm comprehensive mechanized coal mining method or the caving mining method has been adopted at all working faces, and the caving method is used to manage the roof. The mining thickness of each working face is 1.5–2.2 m, and the coal seam inclination is 1–3°, being nearly horizontal. This article mainly studies data taken from the 21308 working face of the No. 3 coal seam from December 2015 to February 2020. The strike length is 930 m, and the incline length is 210 m.Figure 2
Geographical location of the study area.An investigation performed in May 2017 showed that only the Donganshang group exhibited multiple cracks in houses and 3 cracks in the road pavement, with a crack length of approximately 3 m, a crack width of approximately 0.05–0.1 m, and an angle between the crack and the road surface of approximately 90°, which caused the ground to form a bulge of 0.1–0.2 m. There were 3 large-scale cracks in the ground, which were approximately 500 m long and 0.2–0.6 m wide, and 3 small-scale cracks, which were approximately 40~80 m long and 0.1–0.3 m wide. In the village, there were 7 landslides, with lengths of 8–15 m, widths of 5–12 m, and drops of 2–5 m.As research data, 109 scenes of Sentinel-1A radar satellite data in the study area were selected, which cover a time span from December 24, 2015, to February 25, 2020. The Sentinel-1A data are in the C-band, with a wavelength of 5.66 mm; the orbit revisit period is 12 days, the polarization mode is VV, the incidence angle is approximately 33.7°, the resolution of these data is5×20m, and the width of a single image is 250 km. The data used in this experiment are all products of the interferometric wide-width mode (IW) and are single-view complex data with slant distance information [32].The SRTM1 digital elevation model (DEM) was chosen as a source of external DEM data; its ground resolution is30×30m. For the precise orbital parameters, the precise orbit determination ephemeris data provided by the European Space Agency (ESA) (https://qc.sentinel1.eo.esa.int/) were adopted, which can be downloaded 21 days after acquisition of the Sentinel-1A data and have a positioning accuracy of within 5 mm [33]. External DEM and orbital data can be used to remove the levelling effect and orbital errors.


## 3.2. Data Processing

The SBAS-InSAR data process was realized based on the GAMMA Software. The main steps were as follows: extracting public bursts from the 109 Sentinel-1A scenes to improve the subsequent processing speed and efficiency, performing precise registration with an accuracy of 0.001 pixels, and setting a short time baseline (three consecutive scene images). The research area was cropped out based on the corresponding row and column numbers, and multiview interference processing was performed at a visual number ratio of 5 : 1 (range: azimuth). Differential interference, filtering (Goldstein), and unwinding (minimum cost flow method) were performed based on the initial baseline and the precise baseline. The high coherence coefficient method was used to extract points with high coherence, and elevation and atmospheric errors were removed through regression analysis of the high-coherence points. Finally, complete and accurate deformation rate and time-series cumulative shape variables were obtained [34–36].


## 3.3. SBAS-InSAR Results and Analysis

Following the above data processing method, the annual average subsidence rate (Figure3) and the time-series cumulative subsidence map (Figure 4) along the line-of-sight (LOS) direction were obtained for the mining area of Hancheng from December 2015 to February 2020. Figure 4 shows the timing accumulation of 11 scenes selected from among all 109 scenes at equal intervals. Figure 3 shows that the entire study area contains three relatively large subsidence areas: the Donganshang subsidence area, Zhuyuan village south subsidence area, and the Shandizhaizi subsidence area. A field investigation showed that there is good consistency between the subsidence area on the Donganshang and the distribution of the working face in the south No. 1 mining area in Hancheng. Both the subsidence area to the south of Zhuyuan village and the Shandizhaizi subsidence area have undergone coal seam mining, and mountain slippage, road cracks, and house cracks have occurred.Figure 3
Average subsidence rates from December 2015 to February 2020.Figure 4
Time-series cumulative subsidence from December 2015 to February 2020.The results presented in Figures3 and 4 show that the maximum annual average settlement rate in the study area is 300 mm/a and the maximum cumulative settlement is 1000 mm from December 2015 to February 2020. The area on the Donganshang presents the most serious subsidence, with a subsidence area of approximately 3.66 km2 and a cumulative subsidence of approximately 150–1000 mm. All of the working faces in this mining area show different degrees of subsidence. The 21506, 21307, and 21308 working faces stopped being mined in February 2017, and their subsidence rates are approximately 40–295 mm/a. The 21309 working face corresponds to the area with the most severe deformation during the study period. As of February 2020, the settlement rate of the 21309 mining face is the highest. The settlement area to the south of Zhuyuan village is approximately 2.56 km2, the settlement rate is approximately 40–257 mm/a, and the cumulative settlement is approximately 100~1000 mm. The Shandizhaizi settlement area is approximately 0.68 km2, the settlement rate is approximately 40–129 mm/a, and the cumulative settlement is approximately 100–570 mm.It can be seen from Figure4 that each mining face in the mining area shows a different degree of nonlinear settlement. Let the first image (December 24, 2015) be considered the initial image, in which the settlement amount is 0. The three settlement areas had formed by February 16, 2017. Over time and with the mining of the working faces, the subsidence area and the cumulative subsidence continued to increase. As of November 2, 2018, the subsidence area had gradually expanded westward, forming three obvious subsidence basins. As of February 25, 2020, the cumulative settlement had reached its maximum; the impact of the settlement had continued to expand from east to west, and the settlement area and cumulative settlement had gradually increased, indicating that the settlement area had not reached a stable state by that time.To further analyse the settlement information of the study area, the 21308 working face in the settlement area on the Donganshang is taken as the object of investigation, and profile lines A1A2 and B1B2 are established along the direction and inclination (the position of the section line is shown in Figure4). The working face trend is advancing along the direction from A2 to A1, and the working face tendency is advancing along the direction from B2 to B1. For image readability, the data from 10 scene images corresponding to some of the time-series cumulative settlement figures were selected to draw the corresponding strike and longitudinal time-series profile lines using the Origin software (Figure 5). Figure 5 shows that the main settlement areas are between points 11 and 41, while the maximum sinking position is always at point 26. The maximum accumulated settlement amounts to 748 mm and 880 mm. In general, over time, the settlement of each point along the direction and inclination of the working surface is gradually increasing, and an obvious settlement funnel is forming on the spatial scale, which is consistent with the surface mining subsidence in the mining area.Figure 5
Subsidence curves of the working face profiles: (a) A1A2 sinking curve; (b) B1B2 sinking curve.
(a)(b)Because this certain mining area in Hancheng has been repeatedly mined at multiple coal seams, to better analyse the pattern of the mining subsidence, four characteristic points (the positions of these characteristic points are shown in Figure3) on the Donganshang were selected for time-series cumulative settlement analysis, as shown in Figure 6. P1 and P2 are located at the edge of the deformation area on the Donganshang. They are affected by the mining of the 21504 and 21506 working faces, respectively. The cumulative settlement at points P1 and P2 from December 24, 2015, to February 25, 2020, amounted to 255 mm and 392 mm, respectively. P3 and P4 are located in the centre of the settlement area and are mainly affected by mining at the 21308 working face. The cumulative settlement at these points from December 24, 2015, to February 25, 2020, was relatively large, reaching 888 mm at P3 and 838 mm at P4. In general, the cumulative settlement at these four characteristic points has increased over time. As of February 25, 2020, the settlement had not stabilized.Figure 6
Time-series cumulative settlement at characteristic points.


## 3.4. Result Verification and Analysis

To quantitatively verify the ground subsidence monitoring results for the mining area as obtained via SBAS-InSAR technology, GPS observation points established in the subsidence area on the Donganshang were used to obtain ground displacement monitoring data, and the monitoring results were compared and analysed. GPS observation data collected from December 24, 2015, to January 12, 2017, at 8 observation points were selected. Because the surface deformation data obtained via SBAS-InSAR were obtained along the radar LOS direction, the three-dimensional deformation data obtained from the GPS observations were obtained in the vertical, east-west, and north-south directions. To compare the results of these two methods, three steps were required. First, it was necessary to project the vertical, east-west, and north-south deformations observed from the GPS data to the LOS direction. The projection formula is(8)LOS=cosθ−sinθcosαsinθsinαdvdEdNwhere θ is the local incidence angle of the satellite; ɑ is the azimuth of the satellite’s base; and dv, dE, and dN are the vertical, east-west, and north-south deformations, respectively, obtained from the GPS observations.The above formula was used to convert the GPS-observed deformations into deformations in the LOS direction. Then, the coordinates of the GPS observation points were converted into InSAR data coordinates through Gaussian projection to make the two sets of data consistent in their spatial dimensions. Finally, the extracted values were used via the Point tool of the ArcGIS software to extract the InSAR data deformation points corresponding to the GPS observation points in accordance with the neighbouring pixel values using bilinear interpolation to make the GPS results correspond to the InSAR results. Then, the two sets of results were compared and analysed. As shown in Figure7, the overall trends of the deformation values obtained from the SBAS-InSAR and GPS observations in the LOS direction are relatively consistent. Only the GPS results and SBAS-InSAR results of 8 observation points at point 4 are abnormal. Because these results were obtained through a combination of multiple SBAS-InSAR processing experiments, this abnormality may have been caused by human error transmission in the GPS measurements. The smallest deviation of 0.1 cm is observed at point 3. The average absolute error of the 8 points is 1.11 cm, which is within the acceptable range. A correlation analysis of the two sets of data shows that the correlation coefficient is 86%, indicating that they are highly correlated. On the whole, the two sets of results are basically consistent, indicating that SBAS-InSAR technology is reliable for the monitoring of ground subsidence in mining areas and can accurately reflect the corresponding information.Figure 7
Comparison of GPS and SBAS-InSAR results.


## 4. Mining Subsidence Prediction

### 4.1. Model Building

For this experiment, the time-series settlement data based on SBAS-InSAR technology were selected as the training and test samples, and the combined SA-SVR algorithm introduced in the second part of this article was applied to construct a mining area subsidence prediction model and predict the time-series settlement in the mining area. Finally, the SA-SVR prediction results were compared with the prediction results of SVR and GM. To better predict the mining subsidence, a total of 109 sets of data from the Donganshang subsidence area collected from December 24, 2015, to February 25, 2020, were used as the research object in this experiment, and characteristic points were selected along the 21308 working face: points A17 and A25 were selected for the strike section line, and points B30 and B45 were selected for the inclined profile line. Thus, a total of four feature points were used to establish the mining subsidence prediction model.First, the parameters of the SVR model were optimally estimated based on the SA algorithm, and then, the optimal estimated parameters were substituted into the SVR model for settlement prediction. The settlement values at each point selected for the experiment constitute a nonlinear time seriesxi,yii=1,2,⋯,104. The first 104 sets of monitoring data (from December 24, 2015, to December 15, 2019) were used as SVR training samples, and the last five sets of monitoring data (from December 27, 2019, to February 25, 2020) were used as the prediction verification sample. To distinguish between these datasets, they are denoted by(9)xI,yI=x1,y1x2,y2⋮x104,y104,xII,yII=x105,y105x106,y106x107,y107x108,y108x109,y109The first set of dataxI,yI was used to determine the unknown parameters in equation (8), namely, ω^ and b^. The parameter estimates were then reused in the SVR model for settlement prediction using formula (8):(10)y^II=ω^,xII+b^where y^II is the predicted value of yII. The prediction results are listed in Table 1, along with the correlation coefficients. Meanwhile, the displacement map of training, test, and prediction results is shown in Figure 8.Table 1
Prediction results and correlation coefficientsR2 (unit: mm).
Point20191227 SBAS predict deviation20200108 SBAS predict deviation20200201 SBAS predict deviation20200213 SBAS predict deviation20200225 SBAS predict deviationRMSEMAER2A17-623.7-632.6-634.4-637.3-640.1-623.3-627.1-634.9-638.8642.72.82.10.75-0.5-5.50.51.52.5A25-768.1-772.0-772.7-781.8-783.1-763.0-767.8-777.5-782.3-787.24.13.80.51-5.2-4.24.80.64.1B30-740.6-744.7-748.0-756.0-759.2-735.0-740.6-751.7-757.3-762.83.93.70.68-5.6-4.23.71.33.6B45-521.2-523.9-526.2-528.0-528.1-522.8-523.8-525.8-526.8-527.80.90.70.881.6-0.1-0.4-1.2-0.3Note: RMSE (root mean square error) denotes the square root of the ratio of the sum of the squares of the deviations of the observations from their true values to the number of observations; MAE (mean absolute error) denotes the average value of the absolute error and is used to measure the deviation between the observed and true values.Figure 8
The displacement map of training, test, and prediction at 4 points.


### 4.2. Model Verification and Analysis

It can be seen from Table1 and Figure 8 that the prediction model based on the settlement information from the settlement area on the Donganshang is consistent with the SBAS-InSAR monitoring results. According to Table 1, the maximum absolute error between the predicted and measured values is 5.6 mm, and the maximum root mean square error is 4.1 mm. The maximum average absolute error is 3.8 mm, indicating that the prediction accuracy meets the engineering requirements. For the predictions generated by the mining subsidence prediction model established using the SA-SVR algorithm, the minimum value of the coefficient of determination reaches 0.51. The left side of the dotted line in Figure 8 is the original values and training results of the first 104 sets of training samples (December 24, 2015, to December 15, 2019), and the right side of the dotted line is the original values and prediction results of the last five sets of test samples (2019 from December 7, 2020, to February 25, 2020). It can be seen that the model trend is consistent with the training sample trend as a whole, the prediction error of the test sample is small, and the prediction result is good. It indicates that the subsidence prediction model constructed in this paper has high accuracy, fully meets the requirements of engineering applications, and provides a reliable theoretical basis for mining subsidence monitoring and prediction.To verify the superiority of the SA-SVR algorithm for the prediction of mining subsidence in the Hancheng mining area in the absence of mining parameter information for this area, characteristic points in time and space were selected to be used in a statistical analysis method for comparing the accuracy of the deformation prediction results of different algorithms. Figure9 presents a diagram comparing the prediction results for each feature point in the spatial dimension. It can be seen from Figure 9 that the values predicted with the SA-SVR algorithm based on the SBAS-InSAR results show the best fit; the prediction effect at all four feature points is obviously better than that of SVR and GM (1,1), mainly because of the parameter settings of the algorithm. In addition, as shown in Figure 10, the prediction accuracy at each selected feature point was analysed in the time dimension. It can be seen that although the SA-SVR predictions for feature points A17 and B45, with smaller settlement values, are very good, the prediction effect is obviously degraded for feature points A25 and B30, where the values are larger. These findings indicate that the prediction accuracy is affected by the settlement value, with the prediction error being larger when the settlement value is larger.Figure 9
Comparison chart of the results of each point in the spatial dimensions: (a) 20191227. (b) 20200108. (c) 20200201. (d) 20200213. (e) 20200225.
(a)(b)(c)(d)(e)Figure 10
Comparison of the trend at each point in the time dimension: (a) A17. (b) A25. (c) B30. (d) B45.
(a)(b)(c)(d)


## 4.1. Model Building

For this experiment, the time-series settlement data based on SBAS-InSAR technology were selected as the training and test samples, and the combined SA-SVR algorithm introduced in the second part of this article was applied to construct a mining area subsidence prediction model and predict the time-series settlement in the mining area. Finally, the SA-SVR prediction results were compared with the prediction results of SVR and GM. To better predict the mining subsidence, a total of 109 sets of data from the Donganshang subsidence area collected from December 24, 2015, to February 25, 2020, were used as the research object in this experiment, and characteristic points were selected along the 21308 working face: points A17 and A25 were selected for the strike section line, and points B30 and B45 were selected for the inclined profile line. Thus, a total of four feature points were used to establish the mining subsidence prediction model.First, the parameters of the SVR model were optimally estimated based on the SA algorithm, and then, the optimal estimated parameters were substituted into the SVR model for settlement prediction. The settlement values at each point selected for the experiment constitute a nonlinear time seriesxi,yii=1,2,⋯,104. The first 104 sets of monitoring data (from December 24, 2015, to December 15, 2019) were used as SVR training samples, and the last five sets of monitoring data (from December 27, 2019, to February 25, 2020) were used as the prediction verification sample. To distinguish between these datasets, they are denoted by(9)xI,yI=x1,y1x2,y2⋮x104,y104,xII,yII=x105,y105x106,y106x107,y107x108,y108x109,y109The first set of dataxI,yI was used to determine the unknown parameters in equation (8), namely, ω^ and b^. The parameter estimates were then reused in the SVR model for settlement prediction using formula (8):(10)y^II=ω^,xII+b^where y^II is the predicted value of yII. The prediction results are listed in Table 1, along with the correlation coefficients. Meanwhile, the displacement map of training, test, and prediction results is shown in Figure 8.Table 1
Prediction results and correlation coefficientsR2 (unit: mm).
Point20191227 SBAS predict deviation20200108 SBAS predict deviation20200201 SBAS predict deviation20200213 SBAS predict deviation20200225 SBAS predict deviationRMSEMAER2A17-623.7-632.6-634.4-637.3-640.1-623.3-627.1-634.9-638.8642.72.82.10.75-0.5-5.50.51.52.5A25-768.1-772.0-772.7-781.8-783.1-763.0-767.8-777.5-782.3-787.24.13.80.51-5.2-4.24.80.64.1B30-740.6-744.7-748.0-756.0-759.2-735.0-740.6-751.7-757.3-762.83.93.70.68-5.6-4.23.71.33.6B45-521.2-523.9-526.2-528.0-528.1-522.8-523.8-525.8-526.8-527.80.90.70.881.6-0.1-0.4-1.2-0.3Note: RMSE (root mean square error) denotes the square root of the ratio of the sum of the squares of the deviations of the observations from their true values to the number of observations; MAE (mean absolute error) denotes the average value of the absolute error and is used to measure the deviation between the observed and true values.Figure 8
The displacement map of training, test, and prediction at 4 points.


## 4.2. Model Verification and Analysis

It can be seen from Table1 and Figure 8 that the prediction model based on the settlement information from the settlement area on the Donganshang is consistent with the SBAS-InSAR monitoring results. According to Table 1, the maximum absolute error between the predicted and measured values is 5.6 mm, and the maximum root mean square error is 4.1 mm. The maximum average absolute error is 3.8 mm, indicating that the prediction accuracy meets the engineering requirements. For the predictions generated by the mining subsidence prediction model established using the SA-SVR algorithm, the minimum value of the coefficient of determination reaches 0.51. The left side of the dotted line in Figure 8 is the original values and training results of the first 104 sets of training samples (December 24, 2015, to December 15, 2019), and the right side of the dotted line is the original values and prediction results of the last five sets of test samples (2019 from December 7, 2020, to February 25, 2020). It can be seen that the model trend is consistent with the training sample trend as a whole, the prediction error of the test sample is small, and the prediction result is good. It indicates that the subsidence prediction model constructed in this paper has high accuracy, fully meets the requirements of engineering applications, and provides a reliable theoretical basis for mining subsidence monitoring and prediction.To verify the superiority of the SA-SVR algorithm for the prediction of mining subsidence in the Hancheng mining area in the absence of mining parameter information for this area, characteristic points in time and space were selected to be used in a statistical analysis method for comparing the accuracy of the deformation prediction results of different algorithms. Figure9 presents a diagram comparing the prediction results for each feature point in the spatial dimension. It can be seen from Figure 9 that the values predicted with the SA-SVR algorithm based on the SBAS-InSAR results show the best fit; the prediction effect at all four feature points is obviously better than that of SVR and GM (1,1), mainly because of the parameter settings of the algorithm. In addition, as shown in Figure 10, the prediction accuracy at each selected feature point was analysed in the time dimension. It can be seen that although the SA-SVR predictions for feature points A17 and B45, with smaller settlement values, are very good, the prediction effect is obviously degraded for feature points A25 and B30, where the values are larger. These findings indicate that the prediction accuracy is affected by the settlement value, with the prediction error being larger when the settlement value is larger.Figure 9
Comparison chart of the results of each point in the spatial dimensions: (a) 20191227. (b) 20200108. (c) 20200201. (d) 20200213. (e) 20200225.
(a)(b)(c)(d)(e)Figure 10
Comparison of the trend at each point in the time dimension: (a) A17. (b) A25. (c) B30. (d) B45.
(a)(b)(c)(d)


## 5. Conclusion

This paper has combined SBAS-InSAR technology with the SA-SVR algorithm for the first time for long-term mining subsidence monitoring and prediction for a mining area in Hancheng, China. The surface deformation information of the mining area was extracted based on SBAS-InSAR technology, GPS data were used to verify the monitoring results, and the SBAS-InSAR monitoring results were used as training samples to develop a SA-SVR model for predicting mining subsidence. The SVR and GM algorithms were also used to generate predictions for use as a baseline to compare and analyse the SA-SVR results. The conclusions of this study are as follows:(i)
The SBAS-InSAR monitoring results are found to be basically consistent with the GPS observations. The minimum deviation is 0.1 cm. The correlation coefficient between the two is 86%, indicating high correlation. This finding illustrates the reliability of the long-term sequence of SBAS-InSAR monitoring data for the surface subsidence in the mining area and reduces cross-validation error(ii)
The SVR algorithm represents the application of the SVM concept for function approximation and regression estimation. The SVM model is an optimization model proposed by Vaprik [37] in 1995 for regression estimation. To solve the SVM problem, a variety of algorithms have been developed. In 2011, Chang and Lin [38] systematically introduced various mathematical optimization models for solving the SVM problem. The SVM algorithm essentially transforms the original model to be solved into a quadratic mathematical programming optimization model with equality and inequality constraints. Essentially, the original model is expanded; the number of variables changes from the number of unknowns in the regression model to the total number of observations. If there are more observations than unknown parameters of the model, the amount of calculation will increase accordingly. In this paper, the SA algorithm has been directly applied to the original SVM model. By optimizing the parameters of the original model, it becomes more convenient to directly optimize problems with fewer parameters than the number of observations, which may have advantages in terms of convenient programming and ease of calculation. However, in essence, the SA-SVR algorithm obtains an improved solution that is better than the local optimal solution, also called an approximate global optimal solution, at the cost of a large number of calculations(iii)
The method adopted in this paper improves the accuracy of mining subsidence prediction, indicating that the combination of SBAS-InSAR and SA-SVR can be effectively used for mining subsidence monitoring and prediction



---

*Source: 8860225-2020-11-04.xml*