Mapping landslide phenomena in landlocked developing countries by means of satellite remote sensing data: the case of Dilijan (Armenia) area

ABSTRACT Landslide detection and mapping are essential issues for reducing impact of such natural disasters, and for improving the future built-up expansion and planning strategies, especially in developing countries where a reasonable land-use design is an important concern for sustainable growth and environmental management. Armenia is a landlocked country and its urban development is strongly tied to the improvement of infrastructures, which must takes into account the environmental setting and the slope instability of the area, in order to identify risks and possible damages to settlements and economic activities. The use of satellite-based Earth Observation data has advanced significantly in the last decade and has turned out to be very useful for measuring and monitoring slow-moving surface deformation phenomena with millimetric precision. In this framework, this study aims at providing a remote sensing-based Landslide Inventory Map (LIM) and a Landslide Susceptibility Map (LSM) over Dilijan (Armenia) area, performed within the Secondary Cities Urban Development in Armenia project. In particular, LIM and LSM in the study area were produced by using ground deformation measurements derived from satellite Synthetic Aperture Radar (SAR) data, acquired by ALOS and ENVISAT sensors from 2003 up to 2010, and integrated with photo-interpretation of recent optical images and morphological analysis of Digital Elevation Model (DEM). Given the extensive presence of vegetation in the area of interest, satellite SAR images were processed to produce both SqueeSAR™ and Temporary Coherent Scatterers data, which are PSI (Persistent Scatterer Interferometry) data conceived as evolution of PSInSAR™ approach and particularly suited for non-urban and rural areas characterized by low density of coherent terrain benchmarks over time. Landslide mapping produced through this work identifies the most hazardous landslide-affected and landslide-prone areas around Dilijan city, and can be used for further estimating environmental risks for urban infrastructure development in the area.


Introduction
Landslide hazard is one of the major natural disaster that causes damages and losses worldwide (IGOS 2004). Identifying and mapping these phenomena are thus essential for reducing their impact CONTACT Silvia Bianchini silvia.bianchini@unifi.it and for improving the future built-up expansion and planning strategies, especially in those countries where economic and urban growth is developing (Alc antara-Ayala 2002; El-Masri & Tipple 2002;Cascini et al. 2005;Guinau et al. 2005;Xiao 2011). The use of satellite-based Earth Observation (EO) data has advanced significantly in the last decade for a broad range of environmental management and situation awareness applications. In particular, satellites can capture the extent and the movement of large-scale phenomena such as landslides. Consequently, satellite optical images and Synthetic Aperture Radar (SAR) interferometric data, which can measure ground movements with millimetric accuracy across time, can be properly used for supporting landslide mapping activities, due to their non-invasiveness, long-term continuity and wide area covering (Cigna et al. 2010;Bellotti et al. 2014;Ciampalini et al. 2014).
This study focuses on landslide detection and mapping over Dilijan area in Armenia by means of visual interpretation of IKONOS orthophoto and thematic layers, and of Persistent Scatterer Interferometry (PSI), an advanced multi-temporal interferometric technique that exploits satellite SAR data to precisely measure slow-moving surface deformation phenomena (Ferretti et al. 2001;Crosetto et al. 2005).
This work has been carried out in the framework of the Secondary Cities Urban Development in Armenia (SCUDA) project, which is an initiative established to create urban investment and to support secondary cities development in Armenia, and relies on remote sensing EO information.
Armenia is a landlocked country, so its economic development strongly depends on the improvement of infrastructures that can easy the exchanges with other countries. In the early 1990's, relocations of commercial or residential areas and fast growth of some cities, i.e. Dilijan city, have produced important changes in the industrial and urban sectors throughout the country. Therefore, Armenian economic development is tied to infrastructure construction and urbanization, which need to follow reasonable land-use planning criteria for sustainable development and environmental management (Knuth 2006).
Within this framework, this study aims at analysing the slope stability in the wide area of the city of Dilijan through the creation of a Landslide Inventory Map (LIM) and a Landslide Susceptibility Map (LSM) of the area of interest, using remote sensing data such as ground deformation measurements derived from satellite SAR data and photo-interpretation of optical image and morphological analysis of Digital Elevation Model (DEM). Given the extensive presence of vegetation in the study area, satellite SAR images, acquired by ALOS (L-band) and ENVISAT (C-band) sensors from 2003 up to 2010, were processed to produce both SqueeSAR TM and Temporary Coherent Scatterers (TCS) data, which are PSI data conceived as evolution of PSInSAR TM approach (Ferretti et al. 2001). SqueeSAR TM and TCS are particularly suited for non-urban and rural areas, characterized by low density of coherent terrain targets, since the phase information for deriving ground movement is extracted not only from point-wise deterministic objects, but also from distributed scatterers represented by areas of moderate coherence (for SqueeSAR TM approach) and from benchmarks with variable signal-to-noise ratio across time (for TCS approach).
On one hand, LIM is a first mapping step and it addresses the spatial distribution of both past and current landslides by detecting their extension, boundaries, type, velocity and eventually their state of activity (Wieczorek 1984;Mantovani et al. 1996;Guzzetti et al. 2012). On the other hand, LSM is a second further step and expresses the probability of spatial occurrence of slope failures given a set of geo-environmental conditions, by providing a measurement of how prone a given area is to landsliding (Brabb 1984;Vandine et al. 2004;Guzzetti et al. 2005).
LIMs date back to the 1970s (Carrara & Merenda 1976) and nowadays are performed using a purely convention approach (i.e. stereoscopic interpretation of aerial photography) and modern remote sensing technologies. Exhaustive reviews on the application of remote sensing techniques for landslide mapping are given by Mantovani et al. (1996), Metternicht et al. (2005), and more recently by Guzzetti et al. (2012). LSI based on remote-sensing techniques, as the one compiled for the city of Dilijan, mainly relies on the visual interpretation of satellite radar interferometric scenes and optical very-high resolution satellite imagery.
Interferometric approaches are widely exploited primarily to measure ground deformations at specific points and to construct time series of displacement. Through the so-called 'radar interpretation' (first conceived by Farina et al. (2006), (2008)) a geomorphological meaning is assigned to the ground displacements measurements coming from the interferometric analysis. Moreover, the combined use of multi-sensors InSAR data allows to achieve better results by exploiting potentials of different microwave bands (e.g. in C-and L-band) (Bianchini et al. 2013;Herrera et al. 2013;Garc ıa-Davalillo et al. 2014).
Landslide susceptibility is the likelihood of a landslide occurring in an area, on the basis of local terrain conditions (Brabb 1984) and represents the degree to which an area can be affected by future slope failure, i.e. an evaluation of 'where' landslides are likely to occur (Guzzetti et al., 1999(Guzzetti et al., , 2005(Guzzetti et al., , 2006. From a mathematic point of view, landslide susceptibility is the probability of spatial occurrence of slope failures, given a set of geo-environmental conditions assumed as the main triggering factors for landslide occurrences in the area of interest (Chung & Fabbri 1999;Guzzetti et al. 2005Guzzetti et al. , 2006. Several works devoted to concepts, principles, techniques and methods for LSMs have been published in the last 30 years (e.g. Carrara 1983;Brabb 1984;van Westen 1994;Aleotti & Chowdhury 1999;Guzzetti et al. 1999;Catani et al. 2013).
The production of a LSM can be obtained through Geographic Information Systems (GIS) using several different methods, which are classified into four groups (Guzzetti et al. 1999;Chacon et al., 2006), i.e. deterministic, heuristic, statistical and landslide inventory-based probabilistic approaches.
To produce the LSM for the Dilijan area, a simple implementation of Random Forest (RF) was used, following the approach proposed by Catani et al. (2013).
The final outcome of the presented work is a remote sensing-based landslide mapping that identifies the most hazardous landslide-affected and landslide-prone areas around Dilijan city, and that can be useful for further estimating environmental risks for urban infrastructure development.

Study area
The Republic of Armenia is located NE of the Armenian Highland, surrounding the Ararat mountains, and covers an area of about 30,000 km 2 . The area investigated in this work extends up approximately 110 km 2 around Dilijan city, which is placed in the Tavush Province in the western part of the Republic of Armenia ( Figure 1).
Dilijan is a spa resort town, situated on naturally occurring mineral springs in the valley of Aghstev River at an elevation of 1250À1500 m a.s.l. The city area is surrounded by deep woody canyons, steep slopes and forests. The climate is characterized by mild summer and mild winter, and precipitation are mostly observed in spring months.
From a geo-structural point of view, Armenia is located in the seismically active zone that stretches from the Alps through the Caucasus and Central Asia to the Russian Federation, along with Turkey and other earthquake-endangered countries, i.e. Georgia and Azerbaijan (Aleksanyan et al. 2012). The area of Dilijan town is characterized by the presence of significant faults NWÀSE oriented and other minor faults.
The geological lithotypes that outcrop in the study area are volcanogenic and marine-sedimentary rocks, mainly made of PaleoceneÀEocene breccias and tuff sandstones, Holocene lacustrinealluvial sandy or clayey deposits, Miocene limestones and organogenous limestones (Kharzyan 2005).
Armenia is a typically mountainous country with steep topography, which contributes to cause many ground movements, mainly landslides.
Landslides in Armenia occur almost everywhere, but they are most widespread in the northern and southern regions, where folded block-mountains exist (Boynagryan 2009).
The formation and development of landslides in Armenia are caused by the geostructural-morphological conditions of the territory, worsened by anthropogenic factors (Balyan & Markosyan 2005).
The main triggering causes are the relief fragmentation and the presence of steep and convex slopes, the variety of rocks with a different degree of claying and weathering, the presence of large amounts of underground water and outflow of hydrothermal solutions to the surface, high seismic conditions due to the presence of numerous faults with different activity that generate dislocation of separate blocks with differential vertical displacements, and strong earthquakes occurrence (Balyan & Markosyan 2005;Boynagryan 2009).
The human activities that increase favorable conditions for mass movements over Armenian regions include trimming, cutting and overloading of slopes, terrain remoistening during water leakages from irrigation canals, water-supplies and excessive watering for agricultural purposes, slope vibrations during vehicle passages (Boynagryan 2009).
Among the naturally occurring landslides, the seismogenic landslides, which develop during and after earthquakes, as well as creep dislocations of boards of active faults, are very frequent in Armenia, followed by landslides connected with other causes (i.e. over-moistening of slopes during abundant precipitation, etc.).
Seismogenic landslides are characterized by small amplitudes of vertical displacement and great horizontal displacements that are represented by landslides-blocks and floods. Landslides occur in the most active seismogenic zones, while floods are rather related to hydrothermal areas where rocks are weakened, since they liquefy during strong earthquakes and begin to flow, producing long landslide bodies floods (Boynagryan 2009). Landslides in Armenia display various sizes, from small landslides on terraces or valleys and small surface dislocations and slope ruptures with slidings, overflows and slumps to large landslide rockfalls, landslide blocks and flood slides of huge size (Boynagryan 2009).
In the territory of Dilijan town, landslides are generally developed on the right slope of the Aghstev river valley. Overall, more than 160 landslide occurrences of different intensity and sizes can be distinguished (Boynagryan 2009;Aleksanyan et al. 2012). The generation of Dilijan landslides is conditioned by several factors, such as the geomorphological and geologic setting, the intensive weathering and bentonitization of sedimentary rocks, their strong fracturing and fissuring.
In particular, ground instability is triggered by the wide spread of subsoil waters and their outflow, and by seismic activeness determined by differentiated movements of neotectonic blocks along the main faults. Another triggering factor comes from the increasing overload in rock weight on slopes, as a result of works during road laying and buildings construction, excessive watering of farmlands, as well as rock moistening by atmospheric precipitations that make erosion and progress of gravitational processes more intensive (Aleksanyan et al. 2012).
The sliding of masses became more active after the 1988 Spitak earthquake and especially after USSR dissolution when the financing of landslide mitigation and prevention measures has stopped (Balyan & Markosyan 2005;Boynagryan 2009;Aleksanyan et al. 2012).

Methodology
The methodological procedure for landslide mapping in Dilijan area consisted in two working steps.
The first one is the compilation of a new LIM of the study area. The LIM was produced by the combination of various input data layers with SAR data processed through PSI techniques, i.e. SqueeSAR TM and TCS algorithms, using the so-called photo-and radar-interpretation procedures (Farina et al. 2006(Farina et al. , 2008. This analysis consists of the integration of PSI measures with optical data (e.g. ortophotos, optical satellite images) and other geo-thematic data (e.g. topographic, geological maps, DEM, etc.). In particular, the procedure combines the photo-interpretation of the optical images with the 'radar-interpretation' of PSI data (Farina et al. 2008;Cigna et al. 2011), allowing the spatial extension of the point-wise ground motion measures provided by satellite data and the identification of geomorphologic evidences related to landslide geometry, provided by the further input layers. On one hand, the photo-interpretation, as well as the morphological analysis of DEM and thematic layers, permitted to recognize landslide evidences in natural environments such as vegetated and rural areas, but not in urbanized areas where landslide-induced geomorphologic features are instead difficult to be identified, due to dense urban fabric. On the other hand, radar-interpretation is very suitable for mapping ground displacements in urbanized areas where much more radar benchmarks (e.g. buildings, other reflective structures and urban infrastructure) are retrieved; whereas, in vegetated areas, the vegetation coverage prevent a high density of PS targets to be identified (Bianchini et al. 2012).
The second step is the creation of LSM, in order to assess the location where a mass movement is more likely to occur. The LSM was compiled by using as input data the thematic layer maps (i.e. land use, DEM and its derivatives) and the LIM.
To produce the LSM for the Dilijan area, a simple implementation of RF algorithm developed in Matlab environment was used. It is a machine learning algorithm developed by Breiman (2001) that carries out a multivariate classification. The methodology applied to LSM relies on Catani et al. (2013). A simple implementation of the RF is performed, in order to produce an ensemble of LSMs for a set of different model settings, input data types and scales. RF is a combination of Bayesian trees that relates a set of predictors to the actual landslide occurrence. Many widely acknowledged landslide predisposing factors are taken into account as mainly related to the lithology, the land use, the geomorphology, the structural and anthropogenic constraints. In addition, for each factor, we also include in the predictors set a measure of the standard deviation (for numerical variables) or the variety (for categorical ones) over the map unit. As in other systems, the use of RF enables us to estimate the relative importance of the single input parameters and to select the optimal configuration of the classification model. The model is initially applied using the complete set of input variables, then an iterative process is implemented and progressively smaller subsets of the parameters are considered.

Input data
The different types of input data that are needed to perform the landslide mapping procedures through the above-mentioned integrated analysis can be grouped in two categories: thematic layers and SAR data.

Thematic data
Thematic layers include thematic maps such as topographic, geomorphologic, geological and land use data and optical images. In particular, for the Dilijan area, we collected the following data.
IKONOS optical image acquired in 2009, with a pixel resolution of 0.5 m. Shuttle Radar Topography Mission (SRTM) DEM data available at 3 arc-seconds (about 90 m pixel spacing). Digital Elevation Model, generated by both COSMO-SkyMed ascending and descending interferometric pairs acquired in Stripmap mode. After SRTM component removal and under the hypothesis of zero displacement, the interferometric pairs have been filtered using the Goldstein approach for noise reduction (Goldstein & Werner 1998). All DInSAR pairs have been spatially unwrapped using the Ghiglia-Romero approach (Ghiglia & Romero 1994). Dates, temporal and perpendicular baselines (Bt and Bn) of the used COSMO-SkyMed interferometric pairs are reported in Table 1.
Ascending and descending geometries were analysed and used separately. For each geometry, every DInSAR pair, after the unwrapping step, has been converted to meters using the phase to height conversion factors, in order to retrieve the SRTM correction field in meters. On those areas where all pairs show good coherence, a combination of the results was made on the basis of a weighted average on Bn values. The final product has been geocoded and summed to the SRTM, in order to get the corrected SRTM values; empty values are present in case of non-coherent pixels or not visible areas. The ground resolution of the derived DEM cell is 15 £ 15 m. The deviation standard value on elevation measurement is 7 m and it has been retrieved from the statistics of the atmospheric noise power and from the normal baselines of the differential interferograms used in the analysis.
Topographic and cartographic maps, consisting in contour isolines derived by DEM surface and in a cadastral map at 1:2000 scale. Land use map, consisting in a Baseline Urban Classification Map. This map derives from a Pl eiades orthoimage acquired in 2014 with bundle resolution of 50 cm for Panchromatic and 2 m for 4-bands multispectral combined in a GIS with other information sources, i.e. administrative boundaries, socio-economic and environmental data of the study area and a SRTM DEM base map. The land coverage map provides data on urban settlement (highly dense, medium dense or discontinuous urban fabric and road/rail network) and natural and seminatural areas (e.g. green areas and shallow water, irrigated and semi-irrigated field zones, open spaces with no or sparse vegetation etc.). The broadest level of categorization distinguishes five classes, i.e. urban, agricultural, forest, water, wetlands.
Furthermore, some bibliographic information (i.e. Boynagryan 2009) about pre-existing LIM on the Dilijan area was also used to compile an investigation of landslide phenomena of the study area.

SAR data
In this work, we exploit satellite SAR images processed by means of two advanced multi-temporal InSAR techniques, i.e. SqueeSAR TM algorithm and TCS method.
SqueeSAR TM algorithm is a second generation PSInSAR TM analysis, developed in 2010 exploiting both 'point-wise' PS (Persistent Scatterers) and 'spatially distributed scatterers' DS (Ferretti et al. 2001(Ferretti et al. , 2011. PS generally correspond to man-made objects (e.g. buildings, linear structures and rocky areas), while DS are areas corresponding to rangeland, pastures, bare earth, debris fields, scattered outcrops or dry soils. These targets do not produce the same high signal-to-noise ratios of PS, but are nonetheless distinguishable from the background noise and their reflected radar signals are less strong, but still statistically consistent.
The SqueeSAR TM algorithm was developed to process the signals reflected from these low-reflectivity homogeneous areas, but it also incorporates PSInSAR TM , hence no information is lost and movement measurement accuracy is improved (Ferretti et al. 2011). Therefore, a higher number of points are retrieved at an increased confidence of ground motion, by identifying and 'squeezing' all possible ground target information with acceptable coherent levels for estimated optimum phase values for the PSI analysis (Ferretti et al. 2011) SqueeSAR TM also produces improvements in the quality of the displacement time series and has been successfully used in some works for analysing ground motion (e.g. Lagios et al. 2013;Raspini et al. 2013;Notti et al. 2014).
TCS is a further evolution of SqueeSAR TM . On one hand for SqueeSAR TM algorithm, for each DS the so-called coherence matrix is estimated and a good phase stability in all SAR images is required, in order to retrieve a full time-series of displacement values for each measurement point. On the other hand, for TCS approach, targets whose signal-to-noise ratio values vary dramatically over time, are typically discarded, because for these points it is not possible to extract time series of motion. In TCS approach, by exploiting a maximum likelihood estimation, it is possible to extract information on DS exhibiting coherence not necessarily in all the interferograms but only in few ones. Thus, TCS approach permits to retrieve average ground velocity and topographic elevation, achieving a significant increase in the density of measurement points (Ferretti et al. 2012).
For Dilijan study area, we used 27 ENVISAT ASAR (Advanced SAR) scenes acquired in C-band along descending orbits in 2003À2010, and 14 ALOS images acquired in L-band along ascending geometry in 2007À2010 (Table 2).
These ENVISAT and ALOS SAR images were processed by means of both SqueeSAR TM and TCS techniques (Figure 2). SqueeSAR TM processing provided 2,650 ENVISAT radar targets and 17,476 ALOS ones. For each of them, geographic coordinates (latitude, longitude and elevation) with meter precision, average Line Of Sight (LOS) displacement rate over the whole acquisition span, as well as displacement time series (i.e. LOS displacement at each acquisition) with millimetre precision are available.
TCS technique produced raster maps of ENVISAT and ALOS terrain motion data that allowed a better overview of ground movements on the study area, including the higher data coverage also on rural and forest areas where mere point-like radar benchmarks were poorly identified due to vegetation.

Landslide inventory map (LIM)
For the Dilijan area, a LIM has been performed at basin scale and includes landslide detection and mapping. Detection, which means the identification of those features related to topographic surface movements, and mapping were mainly based on visual interpretation of remote sensing data, i.e. thematic layers and SAR data (Figure 3). The pre-existing landslide inventory of the study area included 40 main events, classified in seismogenetic landslides, and in differently triggered phenomena further classified in active and dormant landslides (Boynagryan 2009). Within this previous inventory, the most widespread mapped phenomena are located on the southern slope of the Dilijan valley. Radar and photo-interpretation, as well as analysis of DTM and derived map allow newly detecting and mapping a large amount of instability phenomena affecting the study area (Figure 3). The detailed analysis allows enlarging the boundary of most of the already mapped phenomena in the study area and identifying new ones. A total of 204 landslides have been mapped, covering an area of about 24 km 2 .
In particular, visual interpretation of IKONOS optical image and other auxiliary data (i.e. contour slopes derived from DEM) permitted to detect morphological evidences of landslide phenomena, especially flows and superficial erosion (Figure 3(a)). PSI data permitted to identify the main slides and estimating their mean ground velocity (Figure 3(b)).
As in the previous inventory, the most important phenomena have been detected on the southern slope of the Dilijan valley. Here several widespread landslides have been recognized. They can be probably classified as complex landslides affected by partial reactivations. Furthermore, in this area, other several phenomena have been mapped. They generated in the upper part of the relief where the vegetation cover is absent or very reduced. These events are channelized in the right tributaries of the Aghstev river. Areas with superficial and channelized erosion affecting the surface deposits have been classified as 'superficial erosion'. These areas are characterized by lack of soil and vegetation cover, with an ephemeral drainage.
The state of activity of mapped phenomena can be considered as dormant or inactive for almost all the phenomena by merely evaluating the PSI LOS velocities. The stability threshold range of PSI rates for distinguishing active/inactive landslides was chosen as §1.5 and §4 mm/yr, respectively, considering ENVISAT and ALOS data, fitting the standard deviations of each PSI populations (Bianchini et al. 2013).
Considering the radar data, the active and most dangerous landslide is the one located in the Mets-Tala district, southward Aghstev river (Figure 3).
Another relevant and active landslide is located in the Golovino area (Figure 4), already detected within the pre-existing LIM. This sismogenetic landslide is affected by a partial reactivation in its central part probably related to a more recent superficial complex phenomenon. Also westward of Mets-Tala district, on the left hydrographic side of Aghstev river (Figure 3), an active complex landslide has been detected and mapped. It is characterized by a velocity of ¡11 mm/yr, but its location does not involve urban areas.
The northern slope of the valley is characterized by a lower landslide hazard degree. The mapped phenomena involve smaller areas with respect to the southern flank and the most common phenomena are here related to superficial erosion processes.
Landslides were classified with respect to the prevalent type of movement (Cruden & Varnes 1996), the estimated depth, and the related age. Landslide types were determined on the basis of the local morphological characteristics that are detectable on the available thematic layers.
Pre-existing landslides are represented by flows, slides, and complex landslides. For a small number of landslides (about 15%), the determination of landslide typology was not possible due to lack of information on geologic material or absence of any evidences from visual interpretation of optical images, and so these phenomena have been classified as 'Not determined'.
Flows are rapid movements of soil and debris along open-slope or existing channels. They occur isolated or in groups and are common where debris and loose material are abundant along the slopes, i.e. in the lower sectors of the slope profile and within the valleys. Superficial erosion and flow have been distinguished on the basis of the existence (for the latter) of the source area, the transport area (channel), and the deposition (accumulation) area. In the source areas of flows, minor debris slides, rills, and gullies are present locally. Generally, the transport channels coincide with drainage lines.
Slides and complex/composite landslides are abundant, and concentrated in groups in the area of Dilijan. Small slides are likely shallow translational movements located inside other landslides and along undisturbed slopes. Largest complex/composite landslides are deep-seated, rotational and translational slides with a well-defined depletion zone characterized by a concave profile with multiple vertical escarpments, and a distinct bulging deposit characterized by an irregular or convex profile. Deep-seated landslides involve large volumes of material that affect the local morphology and geological structure. The complex/composite landslides are characterized by the combination of two or more types of movement with a spatial and temporal predominance of one of them.

Landslide susceptibility map (LSM)
For Dilijan area, an LSM was created by considering numerous parameters, in order to avoid the subjectivity in the choice of explanatory variables: elevation, curvature (general, planar, profile), slope, aspect, flow accumulation and land use. The source data are represented by a cell of 81 £ 81 m from SRTM DEM data that are characterized by a suitable cell resolution and the most continuous elevation surface for the study area extension. For each parameter the standard deviation for numeric variables and the variety for the categorical ones were also calculated, considering the variability of the eight neighbouring cells (Catani et al. 2013). The training area used to calibrate the model is 10% of the total landslide areas, chosen with a random sampling. The model performs also a feature selection, choosing the variables set that involves the lowest value of Misclassification Probability (MP). The variables used to perform the susceptibility map are then: aspect, aspect variety, DEM, DEM standard deviation, flow accumulation, standard deviation of logarithm of flow accumulation and slope. This configuration is characterized by the lowest MP value (0.15), computed on the test set.
The obtained results are presented in Figure 5: they highlight a large area, along the right bank of the Aghstev River, characterized by very high susceptibility. Other very-high or high-susceptibility zones are presented in the Golovino area and in the lower part of the valleys of the right tributaries of the Aghstev River.
The misclassification to realize a LSM, whatever adopted, is a relevant problem which may led to a decrease of the accuracy of obtained results. The quality evaluation of model accuracy is performed by analysing the agreement between the model results and the observed data that comprises the presence/absence of landslides within a certain terrain unit used for the analysis.
Evaluating the performance of landslide susceptibility models is needed to ensure their reliable application to risk management and land-use planning. The most common errors in the susceptibility models can be summarized in two different groups (Frattini et al. 2010): (1) Error Type I (false positive): terrains not affected by landslides are classified as unstable. This means that their economic value incorrectly decreases because they are limited in their use. In this case the misclassification leads to a social cost equal to the loss of economic value of the terrain. The cost of a terrain usually depends on different environmental (slope gradient, aspect, altitude, etc.) and social aspects (distance from urban areas, roads, etc.). A false positive error implies a reduction of the economic value of terrains and thus corresponds to a private economic loss, for instance in case of residential zoning on sale; (2) Error Type II (false negative): terrains affected by ground deformation are classified as stable and incorrectly used without restrictions. In this case, the cost of the misclassification is higher because they correspond to the loss of element at risk (life, buildings, etc.) that can be involved in a landslide event. The cost is strictly related to the economic value of the element at risk and to their vulnerability. False negative errors are higher in term of social loss with respect to the false positive ones and they can be potentially higher also in term economic losses in case of a landslide event.
The performance of the model of the Dilijan area was evaluated by performing a Receiver Operating Characteristic (ROC) analysis over the whole landslide database. The built ROC curve is reported in Figure 5. The Area Under the ROC Curve (AUC) can be used as a metric to assess the overall quality of a model: the larger the area, the best the performance of the model. A ROC curve is better than another if it is closer to the upper left corner. The AUC value of 0.81 highlights a good result for the LSMs of the Dilijan area.

Discussion
Landslide maps, dealing with the spatial distribution of the main mass movements (LIM) and with the terrain propensity to produce these phenomena (LSM), were created for Dilijan study area by using remote sensing data.
The effectiveness of remote sensing techniques is particularly relevant for landslide mapping in wide and inaccessible regions, for which conventional analysis would not provide comparable potential in terms of timely update, systematic coverage and precision.
On Dilijan study area, the newly generated LIM permitted to map many additional phenomena with respect to the previous landslide inventory already existing in literature (Boynagryan 2009), achieving a more recent and accurate overview of the instability of the area up to December 2010.
The contributions of the analysis rely on analysis of DEM and thematic layers, and on the radarand photo-interpretation procedure, and consisted in the detection of geomorphologic evidences otherwise not emerged from conventional analyses or bibliographic studies, in the verification or modification of already mapped landslide boundaries, in the detection of new phenomena, in the assessment of the main displacement direction and representative ground deformation rate for each landslide, as well as in the definition of the state of activity and the type of phenomenon, if possible.
A total amount of 204 landslides have been mapped. In particular, the most common detected phenomena are related to superficial erosional processes and counted for 75 mapped areas. Translational slides were also detected (52 slides, including 26 areas for displaced material and 28 for main scarps); 25 phenomena were classified as complex landslides (among which 14 refer to displaced material and 17 to main scarps) and 17 as flows.
The presence of highly vegetated areas frequently led to a lack of PSI measurements due to the absence of good radar reflectors.
SAR data only provide information on boundaries and/or activity to 30 landslides (15% of the total LIM), while remaining 177 phenomena were mainly mapped throughout morphological analysis and interpretation of optical image and DEM-derivated layers. Thus, in the Dilijan area, the use of photo-interpretation and DEM morphological analysis turned out to be a useful contribution for landslide mapping, in particular in forested mountainous zones, where a low density of radar targets was identified.
It is important to underline that the use of TCS algorithm with respect to SqueeSAR TM permitted to increase the radar measurements coverage in vegetated zones, which are the most part of the study area, since the information is extracted on those DS exhibiting coherence only in few interferograms and not necessarily in all of them. If we consider SqueeSAR data, a total of 33 landslides (16% of the total) do not include any radar benchmarks at all within their boundaries. This percentage is reduced to 7% (13 landslides) if we consider TCS data.
One of the main cons of TCS is that without prior information it is not possible to reconstruct the displacement time series as for SqueeSAR data, so time series are not available for TCS data. Another drawback is that the precision on estimated elevation values and average displacement rate can be very different within the radar TCS data (Ferretti et al. 2012). The best results can be obtained by a joint use of the two algorithms in order to gain an enhanced insight into ground deformation phenomena: TCS approach permits to better identify the footprint of unstable areas, while displacement time-series provided by SqueeSAR allow to monitor the temporal evolution of the deformation scenario.
It is worth to highlight that, despite the chosen algorithm to process SAR scenes, ALOS data show a better performance on the study area than ENVISAT thanks to a significant higher amount of PS radar benchmarks. This occurs because L-band data allow a more feasible investigation on rural and natural areas than C-band, due to its higher penetration coefficient on the ground and the consequent lower volumetric and temporal decorrelation in vegetated zones (Strozzi et al. 2005;Bianchini et al. 2013;Garc ıa-Davalillo, et al. 2014). In particular, ALOS data permitted to detect and analyse 25 landslides (12% of the whole inventory) that included a sufficient number of PSI data (at least 3 radar targets, Bianchini et al. 2012) with high mean yearly velocity (LOS velocity > 4 mm/yr, towards or away from the satellite). ENVISAT data provided information to only 5 phenomena, which are the ones characterized by a mean yearly velocity higher than the stability threshold ( §1.5 mm/yr). As a result, a total amount of 27 landslides (13%) were classified as active by considering the average yearly velocity measured by radar data. Among these phenomena, 13 were related to landslides (complex, translational slides and flows), 8 were determined as superficial erosion and creep, while 9 were 'Not determined' as it was difficult to distinguish the type of terrain movement.
The most dangerous phenomena triggering urban settlements turned out to be the active slides located in the Mets-Tala district and in the Golovino area that revealed the highest motion rates (>10 mm/yr).
LSM was also performed on Dilijan area, starting from thematic input layers and previously produced LIM.
The LSM provided for Dilijan area achieved a good level of reliability, since it was evaluated by building a ROC curve, and obtaining a good AUC value of 0.81, which can be assumed as an indicator of the overall quality of a used model to provide the landslide susceptibility.
Despite the many advantages, LIM produced by means of remote-sensing techniques suffer some limitations. Firstly, remote-sensing approaches for landslide studies are capable of observing only the behaviour of the surface of the ground and therefore they can evaluate only superficial deformations, while deep displacements cannot be observed, for instance the sliding surface cannot be identified. Another critical issue is that the produced LIM requires support of field survey to validate obtained results. Consequently the compiled inventory map may be subjective products, whose quality depends on the skill and the experience of the investigators, the complexity of the study area, and the reliability of the available information, including the aerial photographs used to identify the landslides (Guzzetti et al. 2000;Malamud et al. 2004;Galli et al. 2008).
Overall, the obtained landslide mapping results highlight the presence of a large area, along the right bank of the Aghstev River, affected by the occurrence of a series of complex landslide, with partial reactivation. This area is also characterized by very high landslide susceptibility. Other landslides have been mapped in the Golovino area and along the valleys of the right tributaries of the Aghstev River. The northern slope of the valley of the Aghstev River where few minor landslides have been mapped, is characterized by a lower landslide hazard degree.

Conclusions
Geospatial information derived from EO satellite data provide added value for detecting natural hazard and for supporting decision-making activities for economic and territorial expansion of landlocked developing countries.
This work led to the newly generated LIM of Dilijan area in Armenia by means of ground deformation measurements derived from SAR data, integrated with optical photo-interpretation and morphological analysis of DEM, and to a LSM delineating hazard areas, using as input data DEM, thematic maps and the remote sensing-based LIM. In particular, this work exploited SqueeSAR and TCS data acquired by ENVISAT and ALOS satellites over Dilijan area in 2003À2010.
The LIM allowed to detect 204 landslides, most of them not previously mapped. Most phenomena were classified as surface erosional processes, even if complex slides and slides are well represented. A good correspondence was found with the LSM produced by a RF implementation. The most critical area was detected along the right bank of the Aghstev River, southward Dilijan town. This area reveals to be affected by high ground movements, corresponding to the occurrences of a series of complex landslides partially reactivated, and turns out to be also characterized by very high landslide susceptibility.