Landslide Mapping in Vegetated Areas Using Change Detection Based on Optical and Polarimetric SAR Data

Abstract

Mapping of landslides, quickly providing information about the extent of the affected area and type and grade of damage, is crucial to enable fast crisis response, i.e., to support rescue and humanitarian operations. Most synthetic aperture radar (SAR) data-based landslide detection approaches reported in the literature use change detection techniques, requiring very high resolution (VHR) SAR imagery acquired shortly before the landslide event, which is commonly not available. Modern VHR SAR missions, e.g., Radarsat-2, TerraSAR-X, or COSMO-SkyMed, do not systematically cover the entire world, due to limitations in onboard disk space and downlink transmission rates. Here, we present a fast and transferable procedure for mapping of landslides, based on change detection between pre-event optical imagery and the polarimetric entropy derived from post-event VHR polarimetric SAR data. Pre-event information is derived from high resolution optical imagery of Landsat-8 or Sentinel-2, which are freely available and systematically acquired over the entire Earth’s landmass. The landslide mapping is refined by slope information from a digital elevation model generated from bi-static TanDEM-X imagery. The methodology was successfully applied to two landslide events of different characteristics: A rotational slide near Charleston, West Virginia, USA and a mining waste earthflow near Bolshaya Talda, Russia.

1. Introduction

Large landslides are a global phenomenon, causing damage and casualties [1]. Landslides arouse emergency situations when urban areas or man-made constructions, such as buildings, bridges, railroads, and roads, are affected. Rapid mapping of landslides is crucial to detect the extent of the affected area, including grade and type of damage. Rapid mapping is a key element of fast crisis response, e.g., to support rescue, humanitarian, and reconstruction operations in the crisis area [2]. Therefore, Earth Observation (EO) based on satellite remote sensing plays a key role due to its fast response, wide field of view, and relatively low cost [3]. (Semi)-automatic landslide mapping based on satellite EO data provides an important information source to support field surveys. Furthermore, during rapid to very rapid events, i.e., deformation rates in the order of several meters per hour [4], access to the landslide area may be too difficult, making field surveys too dangerous [5].

The most frequently used EO data for rapid mapping of landslides is very high-resolution (VHR) optical satellite imagery [6,7,8]. A common way for landslide detection is the mapping of rapid changes of the vegetation layer derived from vegetation indices calculated for pre- and post-event optical EO imagery, e.g., [9,10,11,12,13,14,15].

However, as heavy rain events are one of the most frequent triggers for landslides, i.e., there is a high probability of cloud coverage right after the event, optical EO data is not always useful for rapid mapping applications [16,17]. The advantages of SAR compared to optical EO sensors are (I) the day and night availability of this active sensor, and (II) its almost complete weather independency due to its longer wavelength. In most cases, SAR EO data of a given crisis area is available earlier than cloud-free optical imagery. Therefore, faster disaster response is enabled by SAR-based rapid mapping procedures. However, for both optical and SAR sensors, the revisit time of the satellite has to be taken into account [18].

In addition to satellite-based methods for landslide mapping, there are also other approaches using for instance airborne laser altimetry (LiDAR), e.g., [19], or advanced field mapping techniques such as the combination of a laser rangefinder binocular combined with a GPS receiver [20], are described in the literature. Despite of their promising results, and contrary to satellite data-based landslide mapping methodologies, those approaches are not suited for a world-wide use due to their limited availability.

According to Czuchlewski et al. [21], one post-event single-polarized SAR image is insufficient for distinguishing and mapping landslides. The investigation of the temporal development of the interferometric coherence by analyzing a time-series of SAR imagery, including pre- and post-event imagery in the data stack, enables landslide detection, e.g., [22,23]. While multi-temporal SAR interferometry enables long-term monitoring of extremely slow and very slow movements, e.g., [24,25,26], speckle tracking approaches are able to measure higher deformation rates (up to tens of meters) [27].

Polarimetric SAR (PolSAR) data of at least two polarizations (dual-pol) provides more information on the ground, which enables a better land cover classification and landslide mapping. Quad-pol data, containing the full polarimetric backscattering (i.e., all four combinations of horizontal (H) and vertical (V) polarized waves) allows the most accurate land cover mapping [28] using SAR data.

Based on airborne L- and P-band quad-pol imagery, Rodriguez et al. [29] analyzed changes in the pedestal height and the Radar Vegetation Index (RVI) over time to detect landslides in Taiwan, which were triggered by the 1999 Chi-Chi earthquake.

Cui et al. [30] investigated landslides in earthen levees by means of a multi-classifier decision framework for textural features (grey level co-occurrence matrix) derived from multi-polarized SAR imagery.

Plank et al. [31] compared object-based landslide detection methods based on PolSAR (dual-pol TerraSAR-X) and VHR optical imagery for a case study in Taiwan. The PolSAR procedure is based on a textural analysis with focus on the Normalized Difference Standard Deviation (NDSD) of the calibrated intensities of both polarimetric channels.

Decomposition procedures based on quad-pol SAR imagery, such as the Freeman-Durden decomposition [32] and the further enhancement of it, the Yamaguchi decomposition [33,34], allow the derivation of surface (e.g., bare surfaces), volume (e.g., vegetation) and double-bounce (e.g., man-made objects and at tree trunks) scattering components. Watanabe et al. [35], Yamaguchi [36], Shibayama and Yamaguchi [16,37] report landslide detection procedures by detecting changes of the polarimetric scattering components. Landslides in vegetated areas cause a decrease of volume scattering, i.e., a loss of vegetation, and an increase of surface scattering (bare surfaces). Shibayama et al. [38] found that the local incidence angle has high influence on landslide detection based on polarimetric scattering analysis.

In addition to the aforementioned Freeman-Durden decomposition, Yonezawa et al. [39] investigated also the change of the entropy/anisotropy/alpha (H/A/α) decomposition [40] in pre- and post-event ALOS/PALSAR imagery. H showed lower values for landslide areas than for forested areas. However, farmlands showed similar low values of H as landslides, making the differentiation of these classes very difficult.

Except for [30,31], all aforementioned landslide detection procedures are based on change detection approaches of pre- and post-event VHR SAR imagery, requiring identical imaging geometries of both acquisitions. However, VHR archive SAR imagery acquired shortly before a landslide event, especially at the same image acquisition geometry as the next possible post-event imagery, is in most cases not available. Modern VHR SAR missions, such as TerraSAR-X, COSMO-SkyMed, or RADARSAT-2, do not systematically cover the entire world. Each acquisition has to be programmed manually. Furthermore, due to limited disk space on board the satellites and especially due to limited downlink transmission rates, these sensors are not able to provide worldwide coverage within a short time period—i.e., commonly, no archive image recorded shortly before the event is available. Here, we present a fast and transferable landslide detection methodology based only on post-event VHR PolSAR imagery supported by freely available and systematically-acquired pre-event high-resolution (HR) optical data. The post-event VHR PolSAR acquisition can be programmed before the next overpass of the satellite after the landslide event, independent of any geometrical restrictions by a pre-event SAR imagery. The proposed landslide mapping procedure is a semi-automatic change detection approach based on pre-event HR optical imagery of Landsat-8 or Sentinel-2 and post-event VHR PolSAR data (e.g., TerraSAR-X) acquired shortly after the event.

The methodology was successfully applied to two case studies of different characteristics: first, a rotational slide, which occurred on 12 March 2015 at the Yeager Airport near Charleston, West Virginia, USA was investigated. Second, the methodology was tested at a mining waste landslide event, which occurred on 1 April 2015 near Bolshaya Talda, Kemerovo Oblast, Russia (cf. Section 2).

Section 3 describes the developed landslide mapping methodology. Section 4 describes and discusses the results of both test sites. Finally, a conclusion and outlook is given in Section 5.

2. Study Sites and Data

Two landslide events are studied. The first one occurred at the Yeager Airport landslide, Charleston, West Virginia, USA. The second one is a mining waste landslide near Bolshaya Talda, Kemerovo Oblast, Russia.

On 12 March 2015, a large-scale landslide occurred at an artificial slope at the Yeager Airport. The airport, completed in 1947, was constructed atop seven semi-connected hilltops. In 2005, due to new Federal Aviation Administration (FAA) safety regulations, the construction of an Engineered Material Arrestor System (EMAS) was necessary. Therefore, the Yeager Airport had to be extended, leading to the construction of a large artificial slope, being the tallest geosynthetic reinforced slope in North America (horizontal/vertical ratio of 1:1). The construction of the slope was finished in 2007 [41]. The construction of this large artificial slope is described in detail in [42]. The functionality of the EMAS was successfully put to the test on 19 January 2010: a US Airways flight bound for Charlotte aborted takeoff. The CRJ 200 aircraft could not stop before the end of the runway. Fortunately, the jet was stopped approximately 45 m from the edge of the slope by the EMAS. All 34 passengers and crew survived the incident with only minor injuries reported.

First movements at the slope were noticed in June 2013. In the following time, the deformation increased. On 12 March 2015, the slope failed. A secondary failure of the slope occurred on 13 April 2015 (Figure 1). The landslide can be described as rotational debris slide [4]. Further details on this landslide event, as well as a video of the first slope failure, and a drone flight video recorded after the first slope failure are available at the AGU landslide blog of Dave Petley [43]. These videos are also attached to this article to guarantee long-term accessibility.

The available SAR data are two post-event TerraSAR-X HighResolution SpotLight (HS) dual-pol (HH/VV) imagery, acquired on 25 March 2015 (after the first slope failure) and on 16 April 2015 (after the second slope failure). In addition, one pre-event Landsat-8 imagery acquired on 15 January 2015 was used (cf.

Table 1. Timeline of the landslide events and used satellite imagery.

Study Site	Date of Event	Acquisition Date	Satellite 1	Relative Orbit 2	Polarization
Yeager Airport	-	2 April 2014	TanDEM-X	226/Asc.	HH
	-	15 January 2015	Landsat-8	18/33	-
	1st failure: 12 March 2015	25 March 2015	TerraSAR-X	44/Asc.	HH/VV
	2nd failure 13 April 2015	16 April 2015	TerraSAR-X	44/Asc.	HH/VV
Bolshaya Talda	-	14 September 2014	Landsat-8	146/22	-
	1 April 2015	26 April 2015	TerraSAR-X	14/Desc.	HH/VV
		14 August 2015	TerraSAR-X	14/Desc.	HH/VV

¹ TanDEM-X data was acquired in bi-static mode used for DEM generation; TerraSAR-X imagery was acquired in HighResolution SpotLight (HS) mode: 1.2 m (range) × 2.2 m (azimuth) spatial resolution; Landsat-8 provides imagery in 15 m (pan) and 30 m (multi-spectral) spatial resolution; ² relative orbit: For Landsat-8: WSR-Path/WSR-Row; for TerraSAR-X: Relative orbit/Path direction with ascending (Asc.) and descending (Desc.).

).

The second landslide, to which the developed landslide mapping procedure is applied on, is a large mining waste landslide, which occurred near a road between Novokuznetsk and Bolshaya Talda, Kemerovo Oblast, Russia. This landslide occurred on 1 April 2015 at 1 p.m. local time [47]. It can be described as a very rapid earthflow [4]. Unfortunately, no further details about this landslide event are available. However, an interesting video of the landslide is available and provided with this article (Source: [47]). Two post-event TerraSAR-X HS dual-pol (HH/VV) SAR imagery acquired on 26 April and 14 August 2015, as well as one pre-event optical Landsat-8 image acquired on 14 September 2014, were available for the analysis (cf. Table 1).

3. Method

The basic principle of the methodology proposed in this article is to detect landslides via change detection of freely available, systematically-acquired HR optical pre-event and VHR PolSAR post-event imagery. Figure 2 shows the flowchart of the procedure. Assuming land cover changes due to a landslide event, i.e., destruction and removal of the vegetation cover the first step (I) of the object-oriented procedure is the pre-selection of formerly-vegetated areas based on the Normalized Difference Vegetation Index (NDVI) of the multispectral pre-event imagery (e.g., Landsat-8 or Sentinel-2). (II) Next, after polarimetric speckle filtering using the edge-preserving refined Lee filter, (III) the H/α decomposition is applied to the post-event polarimetric SAR image to detect, within the pre-selected areas, regions characterized by low entropy (H) values, i.e., evidence of bare soil or rock (landslide material). (IV) Then, assuming a minimum slope value δ ≥ 20° (cf. Section 3.4) as a necessary requirement for a landslide event, the landslide detection map is refined accordingly. (V) Finally, to reduce the number of false classifications, all detected landslides smaller than a minimum mapping unit (MMU) are excluded (cf. Section 3.4).

3.1. Pre-Event Imagery: Selection of Vegetated Areas

Cloud-free optical (multispectral, MS) pre-event imagery is used to derive vegetated areas prior to the landslide event. Landsat-8 or the recently launched Sentinel-2 sensor are the preferred sources as their imagery is freely available and systematically acquired with a high repetition rate of five days (Sentinel-2 constellation) to 16 days (Landsat-8). In the ideal case, the optical imagery is acquired shortly before the landslide event. However, due to too high cloud coverage no useful imagery might be available, and optical imagery acquired one year before the event could be used. To minimize seasonal effects on the change detection procedure described below, it is important to use optical data acquired in the same season as the PolSAR imagery.

By using data of MS sensors working in the visible and near infrared (NIR) region of the electromagnetic spectrum, one is able to calculate the NDVI, being a proxy for the site’s vegetation density and greenness [48,49,50]. The NDVI uses the difference of the vegetation signature between the RED (0.6–0.7 μm) and NIR (0.7–1.1 μm) channel (Equation (1)) [51]:

$$NDVI=\frac{NIR-RED}{NIR+RED}$$

(1)

The rationing concept makes the NDVI independent of the illumination, atmospheric effects, topography, etc. Consequently, NDVI images acquired at different dates can be compared. The NDVI ranges from −1 to +1. As water has commonly no reflection in infrared, its NDVI is −1. The NDVI value of bare areas (rock, sand, and snow) is less than +0.1. The NDVI increases with vegetation density (NDVI range +0.1 to +0.7) [50]. Vegetated areas in the MS pre-event imagery are selected by setting NDVI > +0.1 as threshold.

3.2. Post-Event Imagery: Selection of Bare Areas

The selection of bare areas (i.e., possible landslide areas) by means of the post-event VHR PolSAR imagery requires several pre-processing steps, which are described in the following.

3.2.1. Polarimetric Speckle Filtering

The speckle effect, which is caused by the interference of the coherent reflected SAR waves of many individual scatterers within a resolution cell, complicates visual interpretation and classification of SAR images. The natural environment, characterized by distributed targets, is mainly affected by the speckle effect. To reduce this effect, polarimetric speckle filtering using the refined Lee filter is applied [52,53]. This filter aims to preserve the structure of the image, i.e., the edges, while filtering homogenous areas. The correlation between the different polarizations is conserved. The refined Lee filter searches for edges in eight directions: in the vertical, horizontal, and two diagonal directions. A kernel window of 7 × 7 pixels is used. Then, the covariance matrix is filtered.

3.2.2. Polarimetric Decomposition

Objects with different geometric and structural properties show different SAR backscatter. Based on physical assumptions, polarimetric decomposition procedures aim to separate these different backscatter types [54,55]. As dual-pol SAR imagery was available for this study, the H/α decomposition proposed by Cloude and Pottier [40] was applied.

This decomposition is based on the eigenvalues λ and eigenvectors of the covariance matrix C, which is shown in Equation (2) for the current dual-pol case (HH/VV):

$$\langle C_2\rangle =\langle \left[\begin{array}{cc}{\left|S_{HH}\right|}^2& S_{HH}{S}_{VV}^{*}\\ S_{VV}{S}_{HH}^{*}& {\left|S_{VV}\right|}^2\end{array}\right]\rangle$$

(2)

with the Sinclair-Matrix S_xy representing the two combinations of transmitted (index y) and received polarization (index x). The superscript * denotes the complex conjugate.

α describes the type of backscattering. α values close to zero indicate domination of surface scattering (single bounce scattering). α values around 45° show domination of volume scattering, caused by multiple scattering inside a volume, such as the crown of a tree or dense vegetation [56]. High α values (up to 90°) represent domination of double-bounce scattering (e.g., in urban area). α_m represents the mean of α₁ and α₂. The former describes the backscattering type of the dominant scatterer and the latter the backscattering type of the second dominant one (Equation (3)):

$${\mathsf{\alpha }}_m=\frac{1}{{\mathsf{\lambda }}_1+{\mathsf{\lambda }}_2}\left[\begin{array}{cc}{\mathsf{\lambda }}_1& {\mathsf{\lambda }}_2\end{array}\right]\left[\begin{array}{c}{\mathsf{\alpha }}_1\\ {\mathsf{\alpha }}_2\end{array}\right]$$

(3)

The heterogeneity of the scattering is represented by the entropy H, which ranges from 0 to 1 (Equation (4)). H = 0 indicates a dominant scatterer such as a corner reflector. H << 1 indicates natural areas free of vegetation, i.e., bare soil/rocks and landslide material. High H values with H close to 1 represent a random mixture of scattering mechanisms, e.g., multiple scattering inside the crown of a tree. Therefore, high H values are an indicator of vegetated areas such as forests.

$$H=\frac{-1}{{\mathsf{\lambda }}_1+{\mathsf{\lambda }}_2}\left[\begin{array}{cc}{\mathsf{\lambda }}_1& {\mathsf{\lambda }}_2\end{array}\right]{\text{log}}_2\left(\frac{1}{{\mathsf{\lambda }}_1+{\mathsf{\lambda }}_2}\left[\begin{array}{c}{\mathsf{\lambda }}_1\\ {\mathsf{\lambda }}_2\end{array}\right]\right)$$

(4)

Rapid landslides [4] remove the vegetation cover. Therefore, the entropy H can be used to detect possible landslides, which are characterized by low H values.

The use of α for landslide detection is more critical, especially in the dual-pol case, which is investigated in this study. As shown for the Yeager Airport landslide, the landslide is indistinguishable in the α image (Figure 3). Landslide material is very heterogeneous. Depending on the geological and environmental setting, landslide material may consist of rocks and debris of different size, as well as trunks and branches of fallen trees. Consequently, within a landslide body, α could show high variable values: α close to zero for rocks/debris and α close to 90° for tree trunks. Therefore, in the following only H is used to differ between landslides and areas not affected by landslides.

The change detection described in Section 3.3 requires all imagery to be in the same coordinate system. Therefore, the H image is orthorectified and map projected, i.e., transformed from the typical SAR geometry (range/azimuth) into a projected coordinate system.

3.3. Change Detection: Mapping of the Landslides

The landslides are detected by change detection of the aforementioned pre-event NDVI (cf. Section 3.1) and the post-event H derived from the VHR PolSAR imagery. The basic concept of the proposed methodology is to detect areas free of vegetation based on low H values at time t_post, which were previously covered by vegetation, i.e., NDVI > +0.1 at time t_pre. This change detection is executed in an object-based image analysis (OBIA) environment using the Cognition Network Language (CNL). First, the H image is segmented using the multiresolution approach based on the Fractal Net Evolution Approach (FNEA) [57,58]. The developed procedure described above uses a scale parameter of 10. The scale parameter is an abstract value to determine the maximum possible change of heterogeneity with no direct correlation to the object size measured in pixels [59]. A compactness (ranging from 0 to 1) value of 0.5 is chosen. The features of interest are natural ones. Therefore, the shape parameter, ranging from 0 to 1 was set to 0.1.

All thresholds mentioned in the following are mean values for the generated objects, with $\overline{X}$ representing the mean value of all pixels within a certain object for variable X. There are three possible cases (Equation (5)). Based on empirical tests, $\overline{H}$ ≤ 0.8 turned out to be best suited to detect landslide areas.

$$\{\begin{array}{cc}{\overline{NDVI}}_{t_{pre}}\le +0.1& \text{areas free of vegetation at time }{t}_{pre}\\ {\overline{NDVI}}_{t_{pre}}>+0.1\mathsf{\Lambda }{\overline{H}}_{t_{post}}\le 0.8& \text{landslide candidate}\\ {\overline{NDVI}}_{t_{pre}}>0.1\mathsf{\Lambda }{\overline{H}}_{t_{post}}>0.8& \text{no landslide candidate}\end{array}$$

(5)

Areas of $\overline{NDVI}$ ≤ +0.1 are free of vegetation at time t_pre and are not considered in the following. Based on a threshold of $\overline{NDVI}$ > +0.1 vegetated areas at the time t_pre before the landslide event are selected (cf. Section 3.1). Then, the polarimetric entropy $\overline{H}$ (cf. Section 3.2.2) is investigated for the areas, which were vegetated at time t_pre. Low H values indicate areas free of vegetation, such as landslide material. Consequently, objects are selected as landslide candidates if they are covered by vegetation at t_pre but are free of vegetation at t_post, i.e., if ${\overline{NDVI}}_{t_{pre}}>+0.1\mathsf{\Lambda }{\overline{H}}_{t_{post}}\le 0.8$ is true (Equation (5)).

3.4. Refinement of Classification by Topographic Information

Assuming a minimum slope δ value as a necessary requirement for a landslide event, the landslide detection map, described in Section 3.3, is refined as follows: only landslide candidates located at $\overline{\delta }$ ≥ 20° are selected as final landslides (Equation (6)).

$$\{\begin{array}{cc}{\overline{NDVI}}_{t_{pre}}>+0.1\mathsf{\Lambda }{\overline{H}}_{t_{post}}\le 0.8\mathsf{\Lambda }\overline{\delta }\ge 20°& \text{final landslide}\\ {\overline{NDVI}}_{t_{pre}}>+0.1\mathsf{\Lambda }{\overline{H}}_{t_{post}}\le 0.8\mathsf{\Lambda }\overline{\delta }<20°& \text{no landslide}\end{array}$$

(6)

To detect the foot of a landslide region-growing into neighboring areas of $\overline{H}$ ≤ 0.8 and $\overline{\delta }$ > 12° is executed. Next, all neighboring landslide objects are merged to a common object.

As the Yeager Airport landslide took place on an artificial slope, which was constructed in the year 2007 (cf. Section 2), the Shuttle Radar Topography Mission (SRTM) DEM [60] from the year 2000 is too old and could not be used for the slope analysis at this site. Therefore, we used a bi-static TanDEM-X [61] dataset acquired on 2 April 2014 to generate via SAR interferometric (InSAR) analysis an up-to-date DEM of 12 m spatial resolution. Therefore, a more accurate measurement of the pre-failure slope was obtained.

Analysis of optical imagery (Landsat, as well as GoogleEarth) showed that the second study site, i.e., the mining waste landslide in Russia, is a very dynamic area with lots of changes that occurred after the last TanDEM-X acquisition on 23 August 2012 over this area. Therefore, neither the SRTM DEM nor a TanDEM-X DEM could be used for slope analysis in this area.

Finally, to decrease the number of false classifications, all detected landslides smaller than a minimum mapping unit (MMU) of 30 m × 30 m are excluded.

4. Results and Discussion

4.1. The Yeager Airport Landslide, Charleston, West Virginia, USA

First, the methodology is applied to Test Case 1: the Yeager Airport landslide, located near Charleston, West Virginia, USA. Figure 4, Figure 5 and Figure 6 show the results of the developed landslide detection procedure based on pre-event optical Landsat-8 imagery and post-event polarimetric VHR TerraSAR-X imagery acquired on 25 March 2015 (after the first failure of the slope) and 16 April 2015 (after the second slope failure). Figure 5 and Figure 6 demonstrate that the landslide is very well detected by the classification. In the SAR image acquired on 25 March 2015, one can also see a small over classification east of the landslide. At the SAR image acquired after the second slope failure the landslide detection procedure was able to detect the main part of the landslide. However, the scarp area of the landslide is not detected. Here, the second slope failure extended the scarp area of the landslide. At this very steep part, geometric distortions such as layover and foreshortening occurred, changing the SAR backscattering values, which influence the landslide detection procedure.

Table 2. Classification accuracies for the Yeager Airport landslide (12 March and 13 April 2015). Overall accuracy (OA), user’s (UA), and producer’s accuracy (PA).

Date	OA [%]	PA Landslide [%]	UA Landslide [%]	PA Other [%]	UA Other [%]	KHAT
25 March 2015	99.9	87.0	67.4	99.9	100.0	+0.759
16 April 2015	99.9	64.3	66.9	99.9	99.9	+0.655

shows the accuracy values of the methodology applied to the Yeager Airport landslide at the two stages of the landslide event. We applied an area-based accuracy assessment. The classification results are compared to a polygon of the landslide derived by visual interpretation of the SAR image and manual digitization. The landslide boundary is clearly visible in the SAR image due to its different roughness compared to the surroundings. Furthermore, photographs taken in the field or by airplanes (Figure 1), as well as the drone video (Supplementary Materials S2), were used to refine the reference polygon of the landslide.

It is important to note that not only the relatively small areas shown in Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9 were considered for the validation procedure, but the entire TerraSAR-X HS scene (5 km azimuth × 10 km ground range). The values of the confusion matrix are described as follows: the overall accuracy (OA), ranging from 0%–100%, describes the ratio of correctly classified area units (e.g., pixels) to the total number of pixels of the satellite scene. The producer’s accuracy (PA), ranging from 0%–100%, gives the percentage of reference data detected by the classification, while the user’s accuracy (also ranging from 0%–100%) describes the percentage of the classification matching with the reference data. Finally, the KHAT statistics gives information about the strength of the correlation between the classification result and the reference data. The KHAT coefficient ranges from 0, representing a completely random match between classification result and reference data, to +1, representing no random match between classification result and the reference [62].

The OA, as well as the values UA and PA, of the class other are very high (99.9%). This is due to a two-class problem, i.e., landslide and other (areas not affected by landslides), with the percentage area of the latter being much higher than the percentage area coverage of landslide. Consequently, the interesting parameters of the accuracy assessment are the UA and PA values of the landslides class.

The good matching between the landslide classification result and the reference data, shown in Figure 5 for the SAR image acquired after the first slope failure, is reflected by the high UA and PA values of the landslide class of ca. 64% and 87%, respectively (cf. Table 2). Contrary to this, the classification results of the image acquired after the second slope failure clearly show lower values for the landslide class PA while the landslide class UA is stable. The KHAT coefficient of the two SAR acquisitions shows a similar behavior.

4.2. Mining Waste Landslide near Bolshaya Talda, Kemerovo Oblast, Russia

Second, the methodology is applied to Test Case 2: the mining waste landslide near Bolshaya Talda, Kemerovo Oblast, Russia. Figure 7 shows a false color composite pre-event Landsat-8 imagery acquired on 14 September 2014 of the mining waste site.

The landslide was covered twice by dual-pol TerraSAR-X HS imagery: On 26 April 2015 and 14 August 2015. Figure 8 and Figure 9 show the corresponding classification results of the developed landslide detection procedure as well as the reference data derived by visual interpretation, i.e., visual change detection between the pre-event optical HR and post-event SAR image data, and manual digitization.

The figures clearly show that the developed methodology is only able to detect parts of the landslide. UA values of the landslide class are very high, with ca. 90%, i.e., most of the classification results are correct. However, the low PA values (landslide class) show that only ca. 50% of the real landslide area is detected (

Table 3. Classification accuracies for the Russian mining waste landslide (01 April 2015). Overall accuracy (OA), user’s (UA), and producer’s accuracy (PA).

Date	OA [%]	PA Landslide [%]	UA Landslide [%]	PA Other [%]	UA Other [%]	KHAT
26 April 2015	96.8	48.2	89.6	99.7	97.0	+0.612
14 August 2015	96.8	49.7	90.0	99.7	97.0	+0.625

). The reason for this is that the methodology assumes that the landslide area was covered by vegetation before the landslide event. However, this is only true for the area at the foot of the landslide, where the landslide ran over an area outside of the mining waste area (formerly vegetated area). The landslide moved from east to west. The mining waste area itself is free of vegetation and, therefore, shows low NDVI values in the pre-event optical imagery (Figure 7). Therefore, the developed methodology is not able to detect the part of the landslide which occurred inside the mining waste area. Furthermore, the entropy H, alone, would only detect the entire mining waste area and not only the landslide area, as the entire mining waste area is characterized by low H values, i.e., the area is free of vegetation (bare soil).

When considering only the area outside of the original mining site as affected landslide area, the PA of the landslide class increases to 83%–90%, while the UA of the landslide class slightly decreases (

Table 4. Classification accuracies for the Russian mining waste landslide (1 April 2015). Only landslide material outside the original mining site is treated as a landslide. Overall accuracy (OA), user’s (UA), and producer’s accuracy (PA).

Date	OA [%]	PA Landslide [%]	UA Landslide [%]	PA Other [%]	UA Other [%]	KHAT
26 April 2015	98.8	83.1	76.4	99.3	99.5	+0.790
14 August 2015	99.1	90.0	81.6	99.4	99.7	+0.850

).

The developed landslide detection methodology based on post-event polarimetric VHR SAR imagery showed promising accuracy values. However, the limitation of this methodology is that only landslides at slopes previously covered by vegetation can be detected.

4.3. General Discussion

Other landslide detection methodologies reported in the literature are based on change detection approaches using pre- and post-event PolSAR imagery [16,29,35,37,39,63]. As the side-looking geometry of SAR systems causes geometric distortions (e.g., shadowing, foreshortening, etc.) identical acquisition geometries for pre- and post-disaster imagery are required for change detection applications. However, as mentioned in Section 1, VHR archive SAR imagery acquired shortly before a disaster event, in particularly at the same image geometry as the next possible post-event acquisition, is in most cases not available. The advantage of the methodology presented in this article is that in addition to freely available and systematically acquired pre-event HR optical imagery only post-event VHR PolSAR is required, which can be programmed before the next overpass of the satellite, independent of any geometrical restrictions by a pre-event SAR imagery. Nevertheless, to guarantee useful image acquisition geometries, i.e., to decrease the influence of layover and shadowing effects, the terrain and the slope’s orientation in space should be considered [64]. Also, the procedure proposed in the current article is based on change detection. Contrary to the aforementioned studies, we use an optical HR imagery, which is systematically acquired and freely available for the entire Earth’s landmass (cf. Section 3.1).

Compared to other published procedures, which apply a pixel-based classification of the PolSAR imagery, e.g., [16,29,35,37], the proposed OBIA methodology enables an incorporation of optical pre-event imagery and also the additional use of DEM data, which strongly increases the classification accuracy.

According to Shimada et al. [17] and Dabbiru et al. [65], L-band is better suited for landslide detection than X-band, as due to its longer wavelength, L-band is characterized by a better penetration through the forest canopy. On the other hand, the shorter wavelength of X-band enables the acquisition of higher spatial resolution imagery compared to L-band. Consequently, smaller and thinner landslides are only detectable by X-band. In summary, one has to find a compromise between the degree of penetration through vegetation (i.e., detection of landslide mass under tree canopy) and spatial resolution (i.e., detection of smaller landslides, not covered by vegetation)—and data availability. L-band imagery is currently only provided by ALOS-2/PALSAR-2. However, its imagery is acquired at a pre-planned acquisition schedule. Contrary to this, X-band satellite missions, such as TerraSAR-X, acquire imagery after an on demand tasking, enabling a very flexible acquisition and fast reaction in case of a disaster, e.g., a landslide event. Therefore, this article presented a methodology based on VHR X-band PolSAR imagery.

The most important information for the landslide identification is the polarimetric entropy H derived from the post-event SAR imagery. As shown in Table 1, the proposed methodology is based on TerraSAR-X HS imagery of 1.2 m (range) × 2.2 m (azimuth) spatial resolution. Based on our experience we can report that landslides with a minimum size of 30 m × 30 m are detectable using only VHR PolSAR imagery. The pre-failure information on vegetation cover can be derived from 10 m to 30 m spatial resolution optical imagery (e.g., Sentinel-2 or Landsat-8). The aforementioned MMU of 30 m × 30 m was also chosen to consider the spatial resolution of Landsat-8. Since very small landslides as mentioned above cover only 1 to 9 pixels in the optical imagery, it is very difficult to identify such landslides in the optical imagery itself. However, this does not influence the results of the methodology, as the optical data is only used to derive NDVI information for selecting pre-failure vegetated slopes. The segmentation of the OBIA approach is based on the higher spatial resolution PolSAR imagery (cf. Section 3.3).

The proposed landslide detection methodology is based on change detection between pre-event optical and post-event PolSAR imagery. The change we focus on is the removal of the vegetation cover by the landslide. However, as slow-moving slides and slide-earth flows often still preserve entire portions of undisturbed vegetation cover, the proposed methodology is limited to the detection of rapid and faster movements [4], i.e., to the detection of landslides where the vegetation cover was removed.

Czuchlewski et al. [21] used VHR L-band airborne SAR imagery of full polarization to map landslides after the 1999 Chi-Chi earthquake in Taiwan. Despite of their promising results, a worldwide applicability of airborne PolSAR sensors is not feasible. Only satellite-based remote sensing enables fast response to disasters with global applicability.

Quad-pol imagery provides more information on the backscattering characteristics than dual-pol imagery. In general, the former achieves higher classification accuracies [28]. The methodology proposed in the current article is based on TerraSAR-X imagery. This SAR sensor provides only dual-pol imagery in operational mode (quad-pol data are only available in an experimental mode). Nevertheless, Section 4.1 and Section 4.2 showed that high accuracy values can be achieved even with dual-pol imagery.

5. Conclusions

This article presented a fast and transferable methodology for landslide detection. The procedure combines post-event Very High Resolution (VHR) Polarimetric Synthetic Aperture Radar (PolSAR) imagery of TerraSAR-X with pre-event multispectral imagery of Landsat-8 or Sentinel-2. First, vegetated slopes are selected in the pre-event imagery using the Normalized Difference Vegetation Index (NDVI) and a Digital Elevation Model (DEM). Second, based on the post-event VHR PolSAR imagery, areas of low polarimetric entropy H, derived by the entropy/alpha (H/α) polarimetric decomposition, are extracted. Low H values represent areas free of vegetation, e.g., areas covered by landslide material, whereas high H values are characteristic for densely vegetated areas, such as forests. Landslides are detected by change detection between the pre-event optical and the post-event PolSAR imagery. More precisely, possible landslides are areas, which are characterized by a high NDVI value in the pre-event imagery and a low H value in the post-event acquisition.

The developed landslide detection procedure is characterized by the following advantages:

• The utilization of SAR imagery allows fast response in a crisis situation due to the day/night availability and almost complete weather independency of the SAR system. As heavy rain events are an important trigger for landslides, optical sensors, relying on a cloud-free sky to be able to provide a useful imagery are, in many cases, not suited.
• The presented methodology requires only freely-available and systematically-acquired pre-event optical high resolution imagery and post-event VHR PolSAR imagery. Other landslide mapping procedures, which are based on change detection using SAR imagery, require pre- and post-event VHR SAR imagery. However, the VHR archive SAR imagery acquired shortly before a landslide event are, in most cases, not available. This is especially true when a certain imaging geometry is required determined by the next possible SAR acquisition over the crisis area. Modern VHR SAR missions, such as COSMO-SkyMed, TerraSAR-X, or RADARSAT-2 do not systematically cover the entire Earth’s landmass.
• The methodology proposed in this article is also based on change detection. However, high-resolution optical imagery of Landsat-8 or Sentinel-2 is used as pre-event information. As these imagery are freely available and systematically acquired on the entire Earth’s landmass at high repetition rates (cf. Section 3.1), it is guaranteed that useful, i.e., cloud-free, pre-event imagery is available for the entire Earth’s landmass. In the ideal case, the optical imagery is acquired shortly before the landslide event. However, in cases where cloud coverage is too high, cloud-free optical imagery acquired at the same season one year before could be used.

The methodology was successfully applied to two landslide case studies of different characteristics: a rotational slide near Charleston, West Virginia, USA, which occurred on 12 March 2015, and a mining waste earthflow near Bolshaya Talda, Russia, which occurred on 1 April 2015.

In the future, the developed methodology will be applied and tested on further upcoming landslide events, also including applicability tests during rapid mapping activities of DLR’s Center for Satellite Based Crisis Information (ZKI). In addition to the Landsat-8 data utilized in this study, imagery of the recently launched Sentinel-2 will also be employed as pre-event information.