Removing InSAR Topography-Dependent Atmospheric Effect Based on Deep Learning

Abstract

Atmospheric effects are among the primary error sources affecting the accuracy of interferometric synthetic aperture radar (InSAR). The topography-dependent atmospheric effect is particularly noteworthy in reservoir areas for landslide monitoring utilizing InSAR, which must be effectively corrected to complete the InSAR high-accuracy measurement. This paper proposed a topography-dependent atmospheric correction method based on the Multi-Layer Perceptron (MLP) neural network model combined with topography and spatial data information. We used this proposed approach for the atmospheric correction of the interferometric pairs of Sentinel-1 images in the Baihetan dam. We contrasted the outcomes with those obtained using the generic atmospheric correction online service for InSAR (GACOS) correction and the traditional linear model correction. The results indicated that the MLP neural network model correction reduced the phase standard deviation of the Sentinel-1 interferogram by an average of 64% and nearly eliminated the phase-elevation correlation. Both comparisons outperformed the GACOS correction and the linear model correction. Through two real-world examples, we demonstrated how slopes with displacements, which were previously obscured by a significant topography-dependent atmospheric delay, could be successfully and clearly identified in the interferograms following the correction by the MLP neural network. The topography-dependent atmosphere can be better corrected using the MLP neural network model suggested in this paper. Unlike the previous model, this proposed approach could be adjusted to fit each interferogram, regardless of how much of the topography-dependent atmosphere was present. In order to improve the effectiveness of DInSAR and time-series InSAR solutions, it can be applied immediately to the interferogram to retrieve the effective displacement information that cannot be identified before the correction.

1. Introduction

InSAR technology is a geodetic measurement technology [1] based on spaceborne Earth observation, which has been developed rapidly in the past three decades. Satellite radar remote sensing has prominent strengths, such as wide-range coverage, high accuracy, and all-day and all-weather applicability [2]. At present, although InSAR technology has been widely used in landslide identification and monitoring [3,4,5,6,7,8], seismic research [9,10,11,12], ground subsidence monitoring [13,14,15], volcano monitoring [16,17,18], glacier movement [19,20] and polar research [21], atmospheric delay error is still one of the major errors limiting the high-accuracy application of InSAR technology [22,23]. When the atmospheric relative humidity changes by 20% in time and space, the atmospheric effect can contribute to a deformation error of 10–14 cm, which is not allowed for the high-accuracy (millimetre level) deformation detection of InSAR [24].

Many scholars have proposed various methods to reduce the influences and errors caused by atmospheric effects. The correction methods are divided into two categories: atmospheric correction with external data or by the model without external data. Atmospheric corrections without external data include atmospheric delay field models [25,26], linear superposition methods [27], topography-dependent turbulence models (referred to as GTTM) [22], and mesoscale meteorological models [28] and linear models [29,30]. The second is the correction method based on external data, which mainly includes Global Navigation Satellite System (GNSS) water vapour data and ground-based meteorological observations correction method [25,31,32], Moderate-resolution Imaging Spectroradiometer (MODIS) and Medium Resolution Imaging Spectrometer (MERIS) water vapour product data correction method [33,34]. Recently, the global medium-resolution numerical atmospheric reanalysis products from the European Centre for Medium-Range Weather Forecasts (ECMWF) have been widely used for the tropospheric delay estimation and correction of individual interferograms and time-series InSAR [35,36,37]. Based on an iterative tropospheric decomposition model that combines continuous GNSS tropospheric estimates with high-resolution ECMWF data, Yu et al. [38] created a generic atmospheric correction online service (GACOS) product with global coverage, all-day, and all-weather data.

The above-generalized method is more appropriate for a broad area, while the study area presents a topography-dependent atmospheric effects with vertical distribution. Among the topography-dependent atmospheric effect correction research, Beauducel et al. (2000) [39] obtained better correction results for topography-dependent atmospheric effects using atmospheric vertical stratification coefficients. Li et al. (2004) [40] proposed a model that integrated continuous GPS (CGPS) and ground meteorological observation data to simulate and correct topography-dependent atmospheric effects in InSAR measurements. Li et al. (2006) [22] established an atmospheric correction model (abbreviated as GTTM) using GPS data, considering the vertical distribution of tropospheric water vapour in terms of topography, which reduces water vapor’s influence on interferograms. External data are often restricted by spatial resolution. Thus, the correction of the topography-dependent atmosphere at a small scale without the help of external data still warrants further study.

This paper proposes a topography-dependent atmospheric correction method based on the MLP neural network model integrated with topographic and spatial data (described in Section 3). The method is applied to the Baihetan dam, comparing the correction results with the traditional linear model and GACOS correction model to evaluate the correction effects of the three models (Section 4). Meanwhile, the method was practically applied to two pairs of interferograms for atmospheric correction in the Baihetan reservoir area, where the topography-dependent atmospheric effects are severe. Potential landslides that could not be detected before the atmospheric correction was extracted using the method (Section 5), which is of practical significance to improve the accuracy of monitoring disaster sites through InSAR technology.

2. Study Area and Data Source

2.1. Study Area

The study area in this research was Baihetan dam, which is located in Ningnan County, Sichuan Province (27.04°N, 102.91°E) and Qiaojia County, Yunnan Province (26.78°N, 103.12°E). It is located at the southeast edge of the Qinghai-Tibet Plateau, with elevations between 1000 and 3000 m, and belongs to the high mountain and plateau geomorphological unit and the Transverse Mountain System in southwest Sichuan and northeast Yunnan. Due to the significant elevation difference, the study area has distinct vertical stereoclimatic features covering all climate types from southern subtropical to cold temperate zones. The area is mainly composed of river erosion landforms, tectonic landforms and ice erosion landforms, with deep river valleys and severe weathering and denudation [41]. Figure 1b shows that the reservoir area is located in the lower reaches of the Jinsha River, and the terrain is high in the east and low in the west. Moreover, the slopes along both sides of the river are steep, the mountains are high, the valleys are deep, and the river valley landforms mostly comprise “V”-shaped development, which is a typical high mountain deep canyon area. Due to impoundment, the water elevation rises dramatically, leading to landslides, mudflows, collapses and other geological hazards in the study area, so a more accurate displacement monitoring of the study area is especially crucial [42]. The real topography of the study area is shown in Figure 1c,d.

2.2. Data Source

A total of 17 Sentinel-1 images, including 10 images of the descending orbit and 7 images of the ascending orbit, were used in this study. The Sentinel-1 satellite incorporates four modes of operation: stripmap mode (SM), interferometric wide swath mode (IW), extra-wide swath mode (EW), and wave mode (WV). In this paper, we selected single-look complex (SLC) images of VV polarization in IW mode, spanning from 21 April 2021 to 4 July 2021, all after the impoundment of the Baihetan Reservoir. The Sentinel-1 single satellite uses topographic observation in progressive scans (TOPS) mode to acquire images on a global scale every 12 days [43]; that is, the satellite has a revisit period of 12 days. Other parameters, such as the incidence angle and resolution of the Sentinel-1 satellite, are shown in

Table 1. Sentinel-1 satellite SAR sensor major parameters.

SAR System Parameters	Values	SAR System Parameters	Values
Date of launch	April 2014	The angle of incidence	29.1–46.0°
Operation band	C	Resolution	5 m × 20 m
Revisit period	12d	Width	250 km
Proposed shooting mode	IW	Polarization mode	HH + HV/VV + VH/HH/VV

.

3. Methodology

3.1. DInSAR Processing

Differential Interferometry Synthetic Aperture Radar (DInSAR) monitor surface displacement was performed using the phase information of 2 SAR images acquired from the SAR sensor at different times in the same area [44]. As shown in Figure 2, SAR images covering the same area were acquired and co-registered with the precise orbit data to form interferograms. Shuttle radar topography mission digital elevation model (SRTM DEM) data, with a 30 m resolution, were used to remove the topographic and flat effect from the interferogram. The Goldstein adaptive filtering method [45] was used for filtering to remove the noise. Finally, a displacement map was acquired after phase unwrapping using the minimum cost flow (MCF) method [46,47].

3.2. MLP-Based Topography-Dependent Atmospheric Correction Method

In the deep learning community, we have witnessed the rapid development of neural networks. The architecture of neural networks developed from the MLP (Multi-Layer Perceptron), to the CNN (Convolutional Neural Network), to the Transformer, and back to the MLP [48]. In this study, we proposed an MLP with ten layers, referred to as FatMLP, for removing the atmospheric phase in InSAR. FatMLP was designed to learn high-order non-linear interactions of the digital elevation model (DEM), longitude and latitude.

The architecture of FatMLP is shown in Figure 2, which consists of ten layers. In FatMLP, an MLP layer was used to learn the linear correlations of the input variables. We used a wider MLP to capture more linear correlations. We designed deeper layers for higher-order nonlinear interactions of the input variables. The input layer included a fully connected layer with 3 neurons, which represent the DEM, longitude and latitude, respectively. Each pixel of the DEM was treated as a sample. First, the input data were lifted into a high-dimensional space, in which the representation of the original data was limited [49]. To this end, the first two hidden layers were fully connected layers with 4096 neurons. Such a high-dimensional space is useful for learning complex interactions of the input variables. Second, the intermediate sparse data were gradually compressed. The third and fourth hidden layers were two fully connected layers with 2048 neurons. The fifth and sixth layers were two fully connected layers with 1024 neurons. The seventh and eighth layers were two fully connected layers with 512 neurons. Finally, the output layer was fully connected with 1 neuron, representing the atmospheric phase. For each hidden layer, the activation function was the Rectified Linear Unit (ReLU), which was used to learn non-linear interactions of the input variables. The depth of FatMLP was greater than conventional MLPs. Deeper neural networks are useful for learning high-order interactions of the input variables. Equations (1)–(5) show the implementation details of the proposed neural network, where FC represents a fully connected layer. The following hyperparameters and configurations are displayed: 9192 is the batch size; the Adaptive Moment Estimation (Adam) is the optimizer, which is one of the most widely used optimization algorithms in deep learning, with the advantages of fewer hyperparameters and good performance; the rate of learning is 0.001; and the proposed neural network was trained for 50 epochs.

$$X_{4096}=Relu\left(FC\left(Relu\left(FC\left(x_{dem},x_l,x_l\right)\right)\right)\right)$$

(1)

$$X_{2048}=Relu\left(FC\left(Relu\left(FC\left(X_{4096}\right)\right)\right)\right)$$

(2)

$$X_{1024}=Relu\left(FC\left(Relu\left(FC\left(X_{2048}\right)\right)\right)\right)$$

(3)

$$X_{512}=Relu\left(FC\left(Relu\left(FC\left(X_{1024}\right)\right)\right)\right)$$

(4)

$$Y_p=FC\left(Relu\left(FC\left(X_{512}\right)\right)\right)$$

(5)

3.3. Linear Model Correction

The distribution of atmospheric water vapour is the primary cause of tropospheric delays, with impacts of atmospheric vertical stratification introducing topography-dependent delayed phases and turbulence introducing local delayed phases. In contrast to delays caused by turbulent mixing, vertical stratification effects are regarded as static stratification in a particular location and time scale. In the case of changes in surface elevation, the vertical stratification tropospheric delay increases as the water vapour concentration and pressure to temperature ratio rise [50]. The concentration of water vapour declines exponentially with altitude, and the theoretical delay curve is an exponential function of altitude [29]. Suppose the interferogram is majorly affected by the static tropospheric delay. In this case, the delay signal is the difference between the corresponding delay curves of the primary and secondary SAR images, and the resulting exponential function is a Taylor-expanded and second-order function. Higher-order terms are neglected, and a linear relationship can be obtained as follows [51]:

$$\Delta Ø=b+Kh$$

(6)

where b is a deviation term, h is the topography, $\Delta Ø$ is the phase, and K is the transfer function between the topography and the phase. The correction of tropospheric delay using a linear model is applicable to most weather conditions, except for special cases, such as reverse or non-monotonic troposphere. Since the elevation–phase relationship may have a local trade-off with the elevation displacement relationship, the transfer function K is best determined through global linear regression [52].

3.4. GACOS Correction Method

GACOS, an InSAR atmospheric correction model developed by Zhenhong Li’ (2006) [22] at Newcastle University, United Kingdom, is a numerical atmospheric product based on a high-resolution European Centre for Medium-Range Weather Forecasts grid of up to 0.1° × 0.1°. The model utilises an iterative tropospheric decomposition (ITD) model to separate the topography-dependent vertically stratified delay component and the turbulent mixing delay component from the total tropospheric delay in order to generate a high spatially resolved (default 90 m) zenith total delay map for InSAR atmospheric delay correction [53]. InSAR users can obtain ZTD products for the corresponding time period by submitting the location of the study area and the time of SAR image acquisition through (http://www.gacos.net) (accessed on 5 January 2022).

4. Results and Discussion

4.1. Analysis of Correction Effect

As shown by the analysis in Figure 3, the interferograms a–f suffered from more serious atmospheric delay effects associated with elevation, and the atmospheric delay phase was related to the elevation change. The atmospheric simulated interferogram (i.e., MLP-modelled atmospheric correction that needs to be removed from the original interferograms) is shown in Figure 3g,h. The topography-dependent atmospheric delay effects in Figure 3m–r were reduced to different degrees after correction by the MLP neural network method, and the phases of the interferograms were smoother. Overall, the effect of atmospheric delay phase removal was more obvious after correction using the MLP neural network method. Meanwhile, the method is not required to build models, does not introduce any external data, has certain applicability for topography-dependent atmospheric corrections, and the atmospheric correction effect is relatively significant.

Table 2. The compared standard deviation of the unwrapping phase before and after correction of MLP neural network.

Interferometric Pair	Phase Standard Deviation (stdDev)
	Original Phase	Corrected Phase	Rate of Change
20210421–20210503	2.5847	1.0807	58%
20210421–20210515	2.5763	0.6896	73%
20210503–20210515	2.4871	0.6379	74%
20210517–20210604	2.5069	0.7856	69%
20210523–20210604	3.1698	0.6192	80%
20210604–20210610	4.6820	0.7508	84%

shows the statistics of the change in the phase standard deviation of the original deconvolution of each interferometric pair. The analysis revealed that the phase standard deviation of all interferometric pairs after correction by the MLP neural network model was significantly reduced compared to the phase standard deviation before correction, with the 20210604–20210610 interferometric pair having the best correction effect and the highest change rate of up to 84%.

4.2. Linear Model and GACOS Correction Method

Figure 4 shows the corrective effect of the Sentinel-1 (21 April 2021 to 4 July 2021) interferogram with the three correction methods. As seen in Figure 4, there was a strong association between the interferogram phase spatial distribution and the topographic relief of the study area, indicating that the delayed vertical stratification of the atmosphere predominates in terms of the overall trend. After correction using the conventional linear model, GACOS, and MLP neural network model, the combination of the second and third rows of Figure 4 demonstrated that the topography-dependent atmospheric delay in the linear model and GACOS correction findings were attenuated to varying degrees. The analysis of the fourth row of Figure 4 revealed, however, that after the correction of the MLP neural network model, the simulated phase of the atmosphere was visibly lowered in comparison to the phase value before correction, and the phase value in the majority of places was close to 0. Thus, the MLP neural network correction had a greater impact than the linear model and GACOS.

To optimize the use of different interferometric SAR atmospheric correction methods and to avoid potential uncertainty, it is necessary to establish appropriate statistical metrics and evaluate the correction performance [54]. The phase standard deviation was applied to measure the effect of atmospheric delay estimation and correction. In particular, the reduction in the interferometric phase was used to demonstrate the reduction in the corrected atmospheric delay phase [36,37]. Using the three aforementioned techniques, atmospheric effects were removed from 27 interferograms, and the standard deviation of the interferometric phase before and after correction was calculated to evaluate the overall performance of the correction results (Figure 5).

Figure 5 shows the magnitude of the reduction in phase standard deviation of the 27 interferometric pairs in the Baihetan reservoir region after atmospheric correction using linear model correction, external data GACOS and MLP neural network, respectively. The analysis revealed that after linear correction and GACOS correction, the phase standard deviation of some interference pairs decreased, which was not significant, while the phase standard deviation of some interference pairs increased instead. Nevertheless, with MLP neural network correction, the phase standard deviation of all 27 interference pairs was reduced, and the reduction was more significant.

Table 3. The number of, percentage of, and average reduction in interference pairs with reduced phase standard deviation and average increment of interference pairs with increased phase standard deviation of the interferogram after the correction using the three correction methods.

The Interferometric Pairs with a Reduced Standard Deviation	Linear Model	GACOS	MLP Neural Network
Number	16	20	27
Percentage	59%	88%	100%
Average reduction	11.1%	17.4%	64%
The average increment of interferometric pairs with increasing standard deviation.	22.1%	10.5%	0

shows that after correction using the linear model, the interferometric pairs with reduced phase standard deviation accounted for 59% of the 27 interferometric pairs, which is the lowest among the three methods. The GACOS atmospheric correction had 88% of the interferometric pairs with reduced phase standard deviation of the unwinding. The percentage of interferograms with reduced standard deviation was larger than that of the linear model. However, the percentage of interferometric pairs with reduced phase standard deviation after correction using the MLP neural network was 100%, which was significantly higher than that of the linear model and GACOS. Hence, for the topography-dependent atmospheric delay in the Baihetan reservoir area, the atmospheric correction results of the MLP neural network were both superior to those of the linear model and the external data GACOS, and the GACOS was superior to the linear model. Compared with the atmospheric delay correction accuracy, the phase standard deviation of all interferometric pairs was reduced by 11.1% on average after using the linear model for atmospheric delay correction, which was the lowest among the three methods. However, after atmospheric delay correction using GACOS, the average reduction in the standard deviation of the unwrapping phase was 17.4%, which was higher than the linear model and located in the second rank. Finally, after the correction of atmospheric delay using the MLP neural network, the average reduction in phase standard deviation was approximately 64%, which was much higher than the first two correction methods, and the correction effect was more obvious for the topography-dependent atmospheric delay in the Baihetan reservoir area.

4.3. Comparative Analysis of Topographic-Phase Correlations

Jolivet et al. (2014) [36] estimated the correlation coefficient between the interferometric phase and elevation as an alternative to the ratio between turbulent and vertical stratified delays. Here, we used the correlation coefficient to describe the elevation dependence in the interferogram. The smaller difference from zero indicates the lower proportion of signals related to the terrain in the interferometric phases. If the absolute value of the correlation decreases after atmospheric correction, it represents the mitigation of the vertical stratified components [36].

As demonstrated by the example of the analysis of the correlation atmospheric delay correction for the elevation of the Baihetan reservoir area in Figure 6, the correlation between phase and elevation for the original interferograms 20210604–20210610 (Figure 6a) was −0.092. After the linear model and GACOS correction, the correlation between phase and elevation was reduced to −0.082 and −0.002, respectively, as shown in Figure 6b,c. However, there was still some correlation between phase and elevation. Figure 6d illustrates that after rectification using the MLP neural network model, the correlation was reduced to −0.000036, which is quite close to 0. The correlation between phase and elevation was rather weak at this stage.

5. Displacement Identification after MLP Correction

When InSAR technology is applied to early displacement monitoring, some minor displacement information is prone to being obscured due to the atmospheric delay effect in the study area, which poses an invisible threat to human life and property. Therefore, it is essential to effectively identify the deformers in some areas through atmospheric correction. In this paper, the interferometric pair 20210517–20210604 with severe atmospheric delay in the Baihetan reservoir area was atmospherically corrected using the above three methods. The results of the atmospheric correction were processed for phase turn displacement, and it was found that there were signs of displacement in the local area. The specific displacement results are shown in Figure 7. The corrected displacement results of the original displacement results, linear model, GACOS, and MLP neural network model are represented in Figure 7a–d, respectively.

As seen in Figure 7, the original displacement results map can be seen in the slope area without signs of displacement, surrounded by atmospheric phase error. It was seriously impacted, and was unable to determine the exact displacement location and displacement magnitude accurately. Combined with Figure 7b,c, the linear model and GACOS correction were ineffective, and the displacement region was not easily identified. However, from Figure 7d, it is clear that the MLP neural network model had a remarkably effective correction effect. In addition to the obvious displacement in the slope ranges, the displacement volume in the non-displacement area was approximately 0mm, and it was able to analyse the displacement area more accurately.

The displacement area, located near Sanjia Village, Yogu Country, Huidong County, Liangshan Yi Autonomous Prefecture, Sichuan Province, with a central longitude of 103°1′19″ and a central latitude of 26°36′59″, was further analysed. This slope elevation was approximately 772–1028 m, with a differential height of 256 m and a relatively steep slope. Combined with Figure 8a,b, the strong displacement area of this slope was more obvious (as shown by the blue dotted line), which was located in the median–lower part of the slope, with a maximum line-of-sight displacement of up to approximately 28 mm.

Meanwhile, the interferometric pairs 20210604–20210610 with severe atmospheric delay in the Baihetan reservoir area were corrected for atmospheric correction using the three methods mentioned above, and the results of the atmospheric correction were processed for phase turn displacement. From the specific displacement results of Figure 9, we found that the local area showed signs of displacement. The original displacement results, linear model, GACOS, and MLP neural network model corrected displacement results are represented in Figure 9a–d, respectively.

Figure 9a showed that the atmospheric delay effect seriously influenced the original displacement results, and the displacement location and displacement magnitude could not be accurately identified. The linear model correction was not effective, and the atmospheric phase remained serious after the correction (Figure 9b). From Figure 9c, it can be seen that the GACOS correction had a certain effect but failed to accurately identify the location of the landslide. Figure 9d reveals that the MLP neural network model had the optimal correction effect. Moreover, with the exception of the obvious displacement in the slope ranges, the displacement volume in the non-displacement area was approximately 0 mm, allowing for a more accurate analysis of the displacement area.

The further analysis of the displacement area showed that the area is located near Dapingzi village, Nagu Town, Huize County, Qujing City, Yunnan Province, with a central longitude and latitude of 103°4′24″ and 26°31′3″, respectively. The elevation of the slope is approximately 1160–1492 m, with a differential height of 332 m, and the slope is relatively steep. Based on the displacement monitoring results shown in Figure 10a, it is evident that the area of strong displacement in this region was more apparent (as indicated by the blue dotted line), with the maximum line-of-sight displacement reaching approximately 40 mm/year. Figure 10b shows the remote sensing interpretation of the displacement area with clear displacement signs. This active slope with displacement was also validated by the SBAS-InSAR result and our fieldwork.

Among the advantages of GACOS are its versatility and global availability. GACOS has advantages in dealing with long-wave atmospheric delay errors and may have some limitations in local applications with large undulating mountains. Therefore, atmospheric errors can be further and more effectively removed in local applications with the help of more refined models and methods (such as the deep learning methods proposed in this paper). It should be mentioned that in order to maintain a high coherence of interferometry, the temporal baseline range was within one month in this study. The effect of the temporal atmospheric water content variations from wet/dry seasons should be evaluated if interferograms with long temporal baselines were used.

6. Conclusions

This study proposed a topography-dependent atmospheric correction method based on the MLP neural network model integrated with topographic and spatial data. The method was used to correct the atmosphere of an interferometric pair composed of Sentinel-1 ascending and descending orbit SAR in 17 views of the Baihetan dam. Moreover, the correction results were compared with the traditional linear model and GACOS to evaluate the effect of correction and the rate of change of the phase standard deviation of the three models.

The results demonstrated that the MLP neural network method could effectively estimate and correct the topography-dependent atmospheric delay, and the corrected Sentinel-1 wide-field InSAR interferogram phase standard deviation mean reduction showed a 64% improvement compared to that of the traditional linear model and GACOS. Furthermore, through fitting the phase-elevation linear relationship, it was clear that the linear coefficient of the phase and elevation correlation after the MLP neural network model correction was almost 0, which means that the phase and elevation are relatively uncorrelated. Eventually, two practical applications confirmed that the MLP neural network model could successfully identify two landslides in Sanjia and Dapingzi villages in the Baihetan reservoir area. The MLP neural network model proposed in this paper is more credible for topography-dependent atmospheric correction results, which can achieve the high-accuracy application of InSAR technology in reservoir area displacement monitoring.