1. Introduction

Argentina shows a complex scenario as regards natural hazards, given by its extension and by the diversity and heterogeneity of its environments (Fernandez Bussy et al., 2010). Landslides are common in the Sierras Subandinas (NW of Argentina); several events have occurred during the past 50 years resulting in significant damage to people and property. The geographic and geologic context of this region together with anthropic activity have favored the occurrence of mass wasting processes and subsequently placed additional strains on the most vulnerable sectors of society.

On 9 February 2009, intense rains triggered landslides in the Tartagal River basin that coalesced and evolved into a debris flow that travelled along the Tartagal River (Fig. 1) and resulted in the destruction of the main bridge of the city of Tartagal and severe flooding. Officially, two people killed, 1,000 evacuated and severe material damages were reported (Baumann et al., 2009, 2011; Brea et al., 2013). Events of this type occurred in the past but none of the magnitude of the February 2009 disaster (Adler, 2009). For instance, destructive flooding occurred in Tartagal town in March-April 2006, following a very rainy summer season. Several landslides and debris flows had taken place within the upper and middle Tartagal River sub-basins that evolved into hyperconcentrated flows that reached the city. High discharge levels resulted in the removal of the channel coating and subsequent deep, lateral and backward erosion that eventually produced the collapse of residential buildings that were lying in the margins. Approximately 130 people had to be evacuated (Adler, 2009; Brea et al., 2013). Within that period in 2006, floods resulting from landslides upstream affected the city of Tartagal more than once and rainfall-triggered landslides took place also in other basins, such as in “Quebrada Yacuy” (Baumann et al., 2009, 2011).

 

Fig. 1. The Tartagal River basin and sub-basins over a color composite that combines the bands centered at 1.63, 0.83 and 0.55 µm in the red, green and blue channels, respectively, of a SPOT-4 image of 15-February-2009. Violet pixels in the upper and middle basin represent the areas affected by landslides and other land cover types that were masked in the classification (see section 3.2.). Yellow pixels represent a sample of the events occurred in 2006, which we used for partial validation of the susceptibility zonation (see Section 3.4.). The inset shows the location of Tartagal River basin within Argentina.

The aim of this study was to map the 2009 events in the upper and middle Tartagal River basin and to elaborate a landslide susceptibility map for the area (Fig. 1). For the landslide inventory map we used the best satellite dataset available acquired just after the disaster and for the elaboration of the susceptibility map, we implemented a heuristic model by using as input a set of conditioning factors that included slope, curvature, lithology and the inventory map of the February 2009 events. We did not consider precipitation as a conditioning factor by using parameters such as accumulated rainfall because data at the appropriate spatial scale was unavailable and the spatial resolution of meteorological models would not have been suitable given the size of the basin (see section 2). As for the choice of the type of susceptibility model selected, different approaches to address this problem are available. These can be broadly grouped in (i) landslide inventory analysis, (ii) statistically based models, (iii) geo-technical or physically-based models and (iv) heuristic modelling (e.g., Guzzetti et al., 1999; Barredo et al., 2000; Rossi and Reichenbach, 2016; Zieher et al., 2017; Lara et al., 2018). Given the limited previous work in the area and in the region (cf. Hermanns et al., 2012), the data, and the resources available, a heuristic approach was the best option for this study. This is not about an improved methodology for landslide susceptibility zonation but an application example, where we have prioritized the generation of preliminary intermediate and final products in order to provide a baseline for future studies and hazard mitigation work.

2. Site description

2.1. The study area

The Tartagal River basin is located in San Martin Department, Salta province, NW of Argentina. The river length is 19 km up to the bridge in the town of Tartagal. The basin covers almost 90 km² and extends along an altitudinal range between 400 and 1,200 meters above sea level (m a.s.l.). The mean elevation is 730 m a.s.l. and the slope ranges between 0 and 84° around a mean of 12° . The upper and middle sub-basins occupy 62 and 21 km², respectively, and altogether they represent more than 90% of total basin’s surface. The town of Tartagal is within the lower river basin at 500 m a.s.l. at the foot of the Sierras Subandinas (Fig. 1) and it had a population in 2010 of almost 80,000 inhabitants (http://www.indec.gob.ar) (Last visit 18/12/2020).

2.2. Climate, vegetation and geology

Rainfall occurring during the summer, i.e., December-March, constitutes the main trigger of landslides in the area. Winds arriving from the East almost dry gain altitude -due to the presence of the Sub Andean hills- and discharge as precipitation all the remaining moisture brought from the Atlantic Ocean. As a result, a sub-tropical mountain forest called "Yungas" covers these hills. This forest extends between 22° S and 28° S latitude, and between 400 and 3,000 m a.s.l. The annual rainfall range is between 900 and 2,000 mm (Grau, 2005) and it exhibits a monsoon regime. The wet season concentrates about 80% of the annual precipitation namely between December and March (Volante et al., 2005).

An asymmetric anticline that heads along a sub-meridian and southwards characterizes the hills, where the Tartagal River basin is located (Mingramm et al., 1979) (Fig. 2). At the transition between the upper and middle sub-basins the nucleus is exposed and shows Upper-Paleozoic sedimentites (Las Peñas and San Telmo formations (Padula and Reyes, 1958)). In the eastern and western flanks, the outcrops belong to Tranquitas Formation, which is from the lower/middle Miocene age. This formation includes sandstones, calcareous sandstones and gray-greenish partially competent pelites. The upper river sub-basin is at the western flank of the anticline, where Tertiary deposits outcrop; they consist of sandstones and highly erodible and reddish pelites (Hernández et al., 1996). On the other hand, conglomeratic and median sandstones compose the Tertiary outcrops towards the east (Gebhard et al., 1974). The areas with low slopes are conformed by Paleozoic rocks within deep valleys characterized by straight sections. Conversely, Tertiary sedimentary rocks form asymmetric homoclinal ridges and hogbacks with secondary rivers flowing parallel to the rock bedding. The ridges and hogbacks slopes, formed by highly erodible sandstones and pelites are prone to shallow landslides (i.e., the 2009 events). However, most of the slope-failures in February 2009 took place on soils and colluvial deposits, with some of them also affecting the underlying sedimentary rocks. Bed scour in rivers located at the foot of slopes has also an important role in the induction of slope failure. In the upper basin, where soft rocks (Tertiary) are present, the Tartagal River valley is broad and winding, whereas in the middle basin where rocks are hard (Paleozoic), the valley is straight and narrow. Finally, the town of Tartagal is located at the apex of the alluvial fan of the river of the same name (Fig. 2).

 

Fig. 2. Geology of the study area; original scale: 1:250,000 (after Bonorino et al., 2001). Silty sandstones and subordinated pelites are associated to the Sub Andean Tertiary unit, brittle sandstones to Tranquitas Formation, pelites and subordinated sandstones to San Telmo Formation and competent sandstones to Las Peñas Formation (see Fig. 5).

The two main types of soils in the Tartagal River basin are from the orders mollisols and entisols (Nadir and Chafatinos, 1994). Entisols develop from fluvial or sandy parent material and distribute along riverbeds and steep slopes. They lack manifestation of soil forming processes and distinct genetic horizons are absent. Mollisols include clays in their profile, a dark colored and well-structured horizon with high base saturation, i.e., a mollic horizon, and a horizon with moderate to high organic matter content (Soil Survey Staff, 1996). Mollisols develop along gentle slopes and from clayey parent material, i.e., Tertiary sandstones and pelites.

2.3. The February 2009 disaster

On 9 February 2009, intense rainfall triggered landslides that evolved into a debris flow, which resulted in severe flooding, two people killed and severe damage to people and property in the town of Tartagal, Salta, Argentina (Fig. 1). The violet pixels in figure 1, represent the 2009 landslides, while yellow pixels correspond to a sample of the landslides that contributed towards the 2006 floods in Tartagal and served for the preliminary validation of the landslide susceptibility zonation map presented by this study (see section 3.4.). Field surveys by the National Geologic and Mining Service of Argentina (SEGEMAR) followed the February 2009 disaster at Tartagal and focused on two areas covered by rocks of the Sub Andean Tertiary unit and Tranquitas Formation (Fig. 2). Most of the landslides that were associated with the former were shallow, whereas in the latter, the majority were larger events that probably resulted from the reactivation of previous historical processes. The two photos at the top of figure 3 correspond to a landslide and a debris flow associated with the Sub Andean Tertiary rocks, while the other two at the bottom were landslides surveyed in a zone associated with rocks of the Tranquitas Formation (Fig. 2).

 

Fig. 3. Photos at top left and top right are a landslide and a debris flow in two sites called Agu. 1 and Agu.2, respectively, in the sector Quebrada del Aguay located towards the W-SW of the basin (S1 in Fig. 4) in the Sub Andean Tertiary (see Fig. 2). The two photos at the bottom are both landslides in two sites called Yar. 1 and Yar. 3 in a zone called Yariguarenda, which lies towards the E-NE (S2 in Fig. 4) within Tranquitas Formation (see Fig. 2).

On February 9, the rainfall gauge located in the village measured 21 mm but local eyewitnesses mentioned that in the higher Tartagal River basin the precipitation might have been higher (Baumann et al., 2009). Debris flows broke the principal railway bridge of the city, the water pipes, sewers and lifelines. Eventually, the deposits of debris flows occupied a significant portion of the urban area. As a result, a large number of persons were deprived of fundamental services for an extended period. The flow reached heights of about 2 m above the ground in the city and its measured discharge rate and velocity were of 1,250 m³ s^-1 and 7,5 m s^-1, respectively (Baumann et al., 2009, 2011).

The rise and development of the town of Tartagal occurred during the 20^th century and there are no reports of landslide related disasters of the magnitude of the one occurred in February 2009 (Adler, 2009). The main cause was initially attributed to the removal of sub-tropical rainforest. However, although selective logging and petroleum extraction has taken place in the upper and middle sub-basins since the foundation of the town of Tartagal in 1924, there is no proof of massive deforestation in the Tartagal River basin. Nevertheless, a substantial increase in agricultural activities (i.e., soybean) has taken place since the 90’s, which led to deforestation and required the channelization and artificial enlargement of the Tartagal River for more than 30 km into the Chaco plain (Sanchez, 2008). The latter may have produced a change in the river base level, which resulted in erosion and channel deepening upstream and in turn in the increased susceptibility to landslides. Thus, the February 2009 disaster could be considered as the result of the combination of a set of unfavorable natural conditions and anthropic actions (Baumann et al., 2009, 2011). Finally, the historical record of El Niño events of the Ministry of Environment of Perú (http://www.minam.gob.pe/fenomenodelnino/el-nino-en-el-peru-y-sus-caracteristicas/registro-historico-de-el-nino/) (Last visit 18/12/2020) indicates two low and two high intensity events occurred during the first decade of the 21^st century. The low intensity were in the summer of the years 2004/2005 and 2006/2007, whereas the high intensity ones, were in the years 2002/2003 and 2009/2010. Thus, there is no apparent relationship of the disaster at Tartagal with an El Niño event. In this study, we intended to retrieve and analyze the most important conditioning factors in order to obtain a preliminary appraisal of the susceptibility to landslides of the river basin.

3. Data sources and methods

3.1. Datasets

We relied on an image acquired on 15 February 2009 by the instrument High Resolution Visible and Infrared (HRVIR) on board the satellite SPOT-4 for the elaboration of the landslide inventory map. These are multi-spectral data of 10 m spatial resolution with four bands centered at the green (0.50 to 0.59 μm), red (0.61 to 0.68 μm), near infrared (0.78 to 0.89 μm) and short wave infrared (1.58 to 1.75 μm) portions of the electromagnetic spectrum. It was the best dataset available acquired just after the event; higher spatial resolution data was unavailable. We applied geometric, radiometric and atmospheric corrections in order to obtain reflectance at the ground surface. We also developed a DEM of 12.5 m spatial resolution by using digital photogrammetry with a stereo-pair of satellite images captured by the sensor ALOS-PRISM on 27 April 2008. This digital topography served for the derivation of topographical variables used in the analysis of conditioning factors. Finally, we used also vector data of the Tartagal River basin’s borders and rivers at a scale 1:50,000 and of the geology at scale 1:250,000 (Fig. 2).

3.2. Landslide inventory map of the 2009 event

The identification of the landslides events of 2009, which led to the elaboration of the inventory map, relied on the classification of the SPOT data by using a standard approach, i.e., the maximum likelihood classifier. In order to train the algorithm, we defined regions of interest (ROIs) based on visual interpretation of color composites. We aimed to detect shallow landslide source areas and debris-flow channels and levees; nevertheless, different types of processes (e.g., landslides, debris flows, etc.) could not be discriminated. We were also unable to differentiate source, transport or deposition areas. The near infrared band was critical in the identification of landslides, since vegetation covers most of the area. The separability of the ROIs was good; both Jeffries-Matusita and Transformed Divergence measures were higher than 1.9. We masked out from the classification water bodies, roads, urban areas, riverbeds and the lower Tartagal River basin. We excluded the lower basin because it included completely modified forest because of anthropic activity and mass movements within it were the result of fluvial erosion and not of slope instability. The output classification map resulted from the implementation of a probability threshold; we worked with a range of thresholds in order to optimize the results.

We used a set of six sites that served as ground truth taken from the field survey that followed the February 2009 events by SEGEMAR. Two of them (Agu. 1 and Agu. 2) were within a sector called Quebrada del Aguay (Q. Aguay) located towards the W-SW of the basin (S₁) (Figs. 3, 4), whereas the other four sites (i.e., Yar. 1 to Yar. 4) were in a zone called Yariguarenda, which lies towards the E-NE of the basin (S₂)(Fig. 4). They were all landslides except for Agu. 2, which was a debris flow. We used these sites to carry out a preliminary and qualitative evaluation of the accuracy of the inventory map of the 2009 landslides events.

 

Fig. 4. Output classification maps obtained with probability thresholds of 0.97 (top left) and 0.87 (top right). The black rectangles indicate the validation zones: S1. Quebrada del Aguay with two sites, i.e., Agu.1 and Agu.2, and S2. Yariguarenda, with four sites (Yar.1 to Yar.4)(Fig. 3).

Finally, the volume of removed material is an important parameter for modelling purposes and in turn for hazard assessment and the development of hazard mitigation measured. Thus, we used the area covered by the 2009 events based on the classification of the satellite data and an average scarp thickness of 0.75 m to estimate the volume of material removed; field observations and geologic and geomorphological analysis indicate that the scarps in the landslide source areas were roughly between 0.5 and 1.0 m thick. Thus, we estimated a volume range and we compared the result with the volume of the debris flow deposit based on observations in the deposition/inundation area in the town. It is important to bear in mind that the volume of the material removed is important but at an order of magnitude scale. In these respects, we required to assume that the landslides events that eventually turned into a debris flow were contemporaneous. Besides, the volume estimate represents the maximum amount of material that collapsed within the perimeters of the upper and middle Tartagal River basins (Fig. 1). While this assumption may lead to an overestimate, we justify this, in part, by considering this estimate as the worst-case scenario for the debris flow that affected the town of Tartagal.

3.3. Analysis of conditioning factors

We used a set of input conditioning factors for the heuristic model that included lithology, slope, curvature and the presence/absence of landslides obtained from the satellite image classification. Many other factors may have conditioned the occurrence of the landslides at Tartagal River basin but it would have been impossible to incorporate them all (cf. Moreiras, 2005).

We evaluated the influence of each factor on the susceptibility separately by defining categories for each one and by ranking them into three discrete levels of relative susceptibility: (i) low, (ii) moderate and (iii) high. To determine whether a low, moderate or high susceptibility would correspond to each category, we used the surface and proportion of landslides mapped within each factor category, out of the total area identified by the classification, and we relied on our knowledge of the study area and the literature (e.g., Ellen et al., 1997; Corominas et al., 2003).

In the case of lithology, the map was a vector layer based on the Geological Sheet 2363-I from SEGEMAR originally at a scale 1:250,000 (Bonorino et al., 2001). For slope and curvature, we used the 12.5 m ALOS-PRISM DEM developed for this study. The relationship between slope and the occurrence of landslides is widely documented; the higher the slope angle, the larger the tensile stresses of the soils and unconsolidated materials (e.g., Guzzetti et al., 2006; Moreiras, 2005, 2009; Martha et al., 2010). As regards curvature, we derived and evaluated the planform curvature, which is normal to the maximum slope’s direction; it is tangential to a contour line. Positive and negative values represent convex and concave surfaces, respectively. Concave plan surfaces favor the convergence of flows and in turn erosion. Thus, the more negative the curvature, the more prone to failure is the surface (Vasantha Kumar et al., 2007; Lee and Evangelista, 2008). Finally, the output classification map itself, i.e., the 2009 landslides inventory map, constituted also an input factor for the heuristic model. To justify this, we relied on the fact that future events are more likely to occur in locations where events took place in the past (e.g., Carrara et al., 1995; Moreiras, 2005, 2009; Ruff and Czurda, 2008; Muñiz and Hernández, 2012). We bear in mind, however, that not all landslides are recurrent in the same area and that the recurrence depends on the types of slope failure. We must also consider the bias imposed by the incorporation of the 2009 landslides inventory map as a conditioning factor into the heuristic model on the resulting susceptibility map. This is however, the best dataset available; a more comprehensive analysis during a larger period at the Tartagal River basin’s scale was not possible.

3.4. Heuristic modelling

First, we carried out a quick evaluation of the independent contribution of each conditioning factor to the landslide susceptibility. This required us to assign a numerical value of one, two and three to low, moderate and high susceptibility levels, respectively and we mapped the susceptibility according to each conditioning factor.

The heuristic model requires on the one hand, a weight for each conditioning factor and the sum of all is one while on the other; each factor has a susceptibility score, which in the case of this study can assume values of one, two or three, as explained above. The final implementation of the model consists of a weighted sum of the susceptibility scores of the conditioning factors as follows:

where S is the susceptibility score or value assigned to each conditioning factor (CF_i) group or category and w is the weight of each conditioning factor. Thus, given the values, that S can assume and that the sum of all weights is one (see above and Table 1), the value of the final susceptibility index (SI) varies between one and three. It must be born in mind that despite the index is a continuous variable it can take a finite number of values.

In order to determine the weights (Table 1) we relied on our expert knowledge (e.g., Castellanos and van Westen, 2008). The landslide areas detected by the classification (i.e., the inventory map of the 2009 events), contributed towards the final susceptibility with the highest weight, topographical factors with an intermediate weight of 0.36 (i.e., slope with 0.23 and curvature with 0.13) and lithology with the lowest contribution (Table 1). The landslide events that led to the disaster at Tartagal started with the infiltration of rainwater in soils, fact that decreased the slope stability leading to the occurrence of shallow landslides and subsequent transformation into debris flows, hence the increased relative importance of the lithology/geology. From the 36% contribution of topography, 63% corresponds to slope; this contrast to the weight of 0.7 decided by Castellanos and van Westen (2008).

In order to obtain the final susceptibility map, we added the four input susceptibility layers weighted by the values in Table 1 within a geographic information system (GIS). The pixel size of this map was 12.5 m, i.e., the same as the spatial resolution of the DEM and topographic variables. Thus, the vector layer with the lithology was converted to a raster with a pixel size of 12.5 m and the inventory map of the 2009 landslides resulted from the classification of SPOT data with a spatial resolution of 10 m was resampled to a pixel size of 12.5 m.

Finally, for verification purposes, we used a dataset that resulted from the interpretation of a true color composite of an image Quickbird-02 from 9-Sep-2006 available in Google Earth. This is a small sample of the events that occurred in early 2006, that covered 0.12 km² and is located in the middle basin of the Tartagal River (Fig. 1); an inventory map of the landslides occurred in early 2006 was unavailable. Therefore, this constitutes a preliminary and partial verification of the landslide susceptibility zonation.

4. Results

4.1. Satellite image classification

The classification of the SPOT data resulted in the February 2009 landslide inventory map. A probability threshold of 0.97 on the maximum likelihood classifier was the maximum possible in order to obtain an output classification map with a coherent pattern of the areas detected as landslides (post-classification clean-up processing included the application of a 3x3-window majority filter). This first output map resulted in an area covered by landslides of 2.9 km²; forested and unclassified pixels covered 68.3 and 7.4 km², respectively (Fig. 4, top left). A threshold 10% lower, i.e., 0.87, did not result in a map with increased noise and included areas of 8.0, 69.5 and 1.2 km² for landslides, forests and unclassified pixels, respectively (Fig. 4, top right). Thus, the introduction of a 10% uncertainty in the probability threshold almost did not change the area covered by forest but it did increase the surface of the detected landsides at the expense of the unclassified pixels. The lower threshold provided an output map that was aesthetically better. In terms of the affected area, the results were within the same order of magnitude.

In Quebrada del Aguay (Fig. 4) both sites surveyed matched locations of pixels classified as landslides at the border with pixels that were unclassified or identified as vegetation with either of the two probability thresholds. In Yariguarenda, Yar. 1 and Yar. 4 matched landslides detected pixels, Yar. 3 matched always unclassified pixels and Yar 2. coincided with landslides when we used the lower probability threshold; with a probability of 0.97 the site matched unclassified pixels. The classifier did not identify any of the six locations as vegetation (Figs. 3, 4).

Finally, based on the assumed scarp thickness of 0.5-1.0 m and the areas obtained with both probability thresholds, we obtained an average volume that ranged between 2.2×10⁶ and 6.0×10⁶ m³. By combining the range of probability thresholds and the range of scarp thickness, the volume varied between 1.5×10⁶ m³ and 8.0×10⁶ m³. Thus, the amount of material removed in February 2009 was within the order of the millions of cubic meters.

4.2. Analysis of conditioning factors

For the evaluation of the conditioning factors, we used the 2009 landslides inventory map based on the 0.87 probability threshold on the satellite image classification algorithm, since this classification left less pixels unclassified and did not result in increased noise (see Section 4.1).

4.2.1. Lithology

The Tartagal River basin features six lithological units that match the formations presented in Section 2 (Figs. 2, 5; Table 2). The largest unit consists of silty sandstones and subordinated pelites (55 km²); it represents about 60% of the whole basin and belongs to the Sub Andean Tertiary. This lithological group lies in the west and it covers most of the upper and about 10% of the middle basin. The classification identified 5.6 km² affected by landslides within this group. This surface represents almost 71% of the 8.0 km² detected by the classification. Thus, we assigned to this group a high susceptibility (Table 2).

Brittle sandstones belongs to Tranquitas Formation and cover 23 km²; it represents about 25% of the whole basin, almost 75% of the middle and 12% of the upper basin. The landslides detected covered 2.2 km² of this unit and represent close to 30% of the classification (Table 2). Thus, we assigned also to this unit a high susceptibility. The quaternary alluvial fan deposits are in the lower basin (they cover 90% of it) and were therefore excluded from the susceptibility zonation (we masked out the lower basin from the classification of SPOT data). The surface occupied by the fourth group (pelites and subordinated sandstones) is close to 3 km² and represents 3% of the basin; half of it is in the upper and half in the middle basin (both areas represent less than 10% of each sub-basin). From the remaining 2% of landslides (see above), 1.5% were detected within this lithology (Table 2). We therefore assigned to this lithology a moderate susceptibility. Finally, competent sandstones are present in the middle basin only, while the 0.4 km² covered by conglomerates and conglomeratic sandstones are mostly in the lower basin (Fig. 5, Table 2). The areas covered by landslides within these latter two lithologies were extremely low. We assigned to both groups a low susceptibility.

 

Fig. 5. Lithology map and landslides detected by the maximum likelihood classifier with a probability threshold of 0.87 (modified from Bonorino et al., 2001 and Baumann et al., 2014).

4.2.2. Slope

The mean slope of the whole basin is 12° (90% of the pixels have slopes that are lower than 23°). At the upper basin, the mean slope is also about 12° whereas the middle basin shows a mean slope that is higher, i.e., 15.5°. Finally, the lower basin has a mean slope of 4°. We reclassified the slopes into four categories. The first group Very Low Slopes includes values of 0° to 7°; these correspond to terraces and alluvial plains. Low slopes are within the range 7°-15° and involve secondary rivers and transverse and lateral valleys. Moderate slopes are between 15° and 30° and correspond to structures-aligned and cutting-hills or tilting valleys. Finally, the category called high slopes includes slopes steeper than 30°. These are related with erosion scarps, which are common in the middle river basin, where the oldest rocks are present (Table 3 and see Sections 2.2 and 4.2.1.).

Low (7.0° -15° ) and very low (0.0° -7.0° ) slopes represented about 40% and 30% of the whole basin, respectively. The remaining 30% consists of moderate slopes except for 3% occupied by steep slopes. Very low slopes cover about 80% of the lower basin, while slopes within the range of 7°-15° represent 15% of it. Within the middle basin, moderate slopes represent about 40% of its surface and low slopes followed with 33%, both with a surface of 7-8 km². Two thirds of high slopes are within this middle basin. The remaining are within the upper basin, where moderate slopes occupy 16 km² but low slopes 26 km², i.e., slopes are generally lower at the upper than at the middle basin (Table 3).

Landslides on slopes lower than 7° occupy 1.8 km², i.e., 20% from the 8 km² identified by the classification (Table 2). If we base the assignation of the susceptibility solely on the incidence of the events within the group, this slope category should have a moderate susceptibility. We decided, however, to assign a low susceptibility to this group because the classification overestimates the areas affected by landslides, since we were unable to discriminate among source, transport and deposition zones nor among different processes (i.e., landslides, debris flows, etc.)(see Section 3.2). Low slopes (the second group) register the highest area of landslides, i.e., they cover a surface of 3.1 km² and represent 40% of the total. For the same reason as in the case of very low slopes, we decided to assign a moderate instead of a high susceptibility to this group (Table 3). Between 15° and 30° slopes, the surface covered by landslides was similar as in the case of low slopes. Moderate and low slopes together gathered 75% of the landslides detected. We assigned a high susceptibility to this group. Finally, landslides on high slopes occupied 0.2 km² and a 3% of the detected areas (Table 3). We assigned a high susceptibility to this category, since it is reasonable to expect at such slopes landslides are more likely to occur. Figure 6 illustrates the slope reclassification together with the landslides pixels mapped by the maximum likelihood classifier with a threshold of 0.87.

 

Fig. 6. Slope ranges and the landslides detected by the image classification (pb=0.87).

4.2.3. Curvature

We grouped the pixels according to three categories of curvature. Concave and convex areas occupy 31 and 34 km², respectively, while flat surfaces encompassed 24 km² (Table 4). In turn, the former two, represent 35% and 38% of the whole basin, respectively, while flat zones the remaining 27%. Within the upper and middle basins (not shown in Table 4), we were able to observe that convex and concave zones occupy 22-24 and 8-9 km², respectively, while flat zones cover 15 and 4 km², respectively. The proportions of concave and convex areas are both around 40% out from the surface occupied by both sub-basins. Finally, at the lower basin, flat zones predominate; they represent 73% of its surface.

At concave locations, the classification detected more landslides, i.e., a surface that represents almost 50% of the whole area obtained from the algorithm. At flat surfaces, the number of landslide-detected pixels was somewhat lower than at convex areas. In turn, the former represents 21%, while the latter the remaining 33% of the landslide detected surfaces (Table 4). Again, we should consider that the classification lead to an overestimation of the areas covered by landslides, especially at convex locations. Thus, we have assigned low, moderate and high susceptibilities to convex, flat and concave areas, respectively. Figure 7 shows how the different curvature classes distribute spatially together with the landslides detected by the satellite image classification.

 

Fig. 7. Curvature ranges and the landslides detected by the image classification (pb=0.87).

4.2.4. The 2009 landslides inventory map

Finally, pixels identified as landslides in the inventory map were assigned a high susceptibility, whereas pixels, where landslides were undetected, were assigned a low susceptibility. From the 8.0 km² of the landslides detected, 87% of this surface took place in the upper and the remaining 13% in the middle basin.

4.2.5. Susceptibility analysis

Lithological units of high susceptibility occupy 95% of the upper and middle basins together. These are the first two lithological groups in Table 1, i.e., "Silty Sandstones and Subordinated Pelites" and "Brittle Sandstones" from the Sub Andean Tertiary and Tranquitas Formation, respectively (Table 2). In the case of slope, moderate susceptibility (7-15°) predominates with 33 km² and a fraction of 40%, while high susceptibility slopes cover about the same surface as the group of slopes of low susceptibility (Table 3). It is important to recall that slopes of the middle basin are higher and coexist with one of the high susceptibility lithological unit “Brittle Sandstones” from Tranquitas Formation (Bonorino et al., 2001). As regards curvature, concave surfaces of high susceptibility cover nearly 40% of the upper and middle basins, while convex areas (i.e., low susceptibility) represent also about 40% and together with moderate (flat) zones occupy the remaining 60% (Table 4). Finally, the 8-km² surface of landslides detected by the classification resulted in 10% and 90% of areas of low and high susceptibility, respectively. Figure 8 summarize the susceptibility given by each of the fours conditioning factors.

 

Fig. 8. Susceptibility according to A. Lithology, B. Slope, C. Curvature and D. Landslides.

4.3. Landslide susceptibility zonation

The weighted sum of the input layers resulted in a raster of the susceptibility, where the mean susceptibility index was 1.8 and its standard deviation 0.3. We used quantile classification (cf. Bathrellos et al., 2009) in order to split up the susceptibility index into four ranges, i.e., (i) Low: 1.00-1.68, (ii) Moderate: 1.68-1.92, (iii) High: 1.92-2.18 and (iv) Very High: 2.18-3.00 (Fig. 9).

The distribution of the susceptibility zones in the final map is quite even. While the two high susceptibility classes occupy 50% of the basins together, both moderate and high susceptibilities represent about 30% of the total area under study. The same is observed with the remaining 40%, where half are low and very high susceptibilities, respectively (Table 5).

As regards the preliminary validation carried out, the 2006 landslides identified were located within the middle basin (yellow pixels in figure 1 and black pixels in figure 9). About 98% of these are located within the high or the very high susceptibility areas; 53% overlaid the zones of very high susceptibility (Table 5).

 

Fig. 9. Landslide susceptibility zonation map.

5. Discussion

The identification of landslides by means of the application of a traditional image classification approach has proved to be rapid and practical given the resources available. The spatial resolution of the SPOT dataset constituted the main limitation, which prevents the identification of events smaller   than 100 m², the differentiation of processes (e.g., landslides, debris flows, rock falls, etc.) and the identification of source, transport and deposition areas. The resulting 2009 landslides inventory map, however, shows good correspondence with the field sites, especially when we used a probability threshold of 0.87; none of the six sites were identified as vegetation (Fig. 4). We were also able to observe that the sensitivity of the output map to the probability threshold imposed to the maximum likelihood classifier was not significant; the landslides detected areas and volumes remained within the same order of magnitude. The average amount of material removed from the upper and middle basins was within the order of the millions of cubic meters. This value matched well the estimate obtained for the volume of the February 2009 debris flow deposit in the town of Tartagal (Baumann et al., 2009, 2011). Besides, the impacts in the town are in accordance with the consequences postulated by Jakcob (2005) for such magnitude, i.e., destruction of towns, obliteration of valleys, river damming. Finally, this volume estimation is essential for modeling purposes and subsequent hazard assessment.

Zones of high and very high susceptibility covered half of the upper/middle Tartagal River basin. Moderate and high susceptibilities areas occupied surfaces that were very similar and the coverage of zones of low and very high susceptibility were almost identical (Table 5). Zones of high susceptibility were not larger probably because steep slopes (>30°) were underrepresented (Table 3) and besides, in conditioning factor groups where a moderate or low susceptibility was expected, the presence of landslides had been overestimated. The latter occurred especially in the case of slope and curvature, namely because different landslide processes nor source, transport and deposition areas could be discriminated. In the middle basin, where slopes are higher than at the upper basin and coexist with brittle sandstones (i.e., Tranquitas Formation) more landslides were expected (Bonorino et al., 2001). This did not occur, however, because most of the precipitation fell within the upper basin (Baumann et al., 2009, 2011). If only lithology is considered, 95% of the total area under study could be considered as of high susceptibility.

The decision for excluding rainfall as a conditioning factor into the heuristic model is based on the lack of appropriate data at a suitable spatial scale. Meteorological modelling of rainfall was not possible because of the coarse spatial resolution of models given the size of the basin under study. As far as the susceptibility analysis is concerned, we bear in mind that results are biased by the fact that the landslides identified by the classification of SPOT data, included also as input conditioning factor of the heuristic model, were triggered by a rainfall event.

Among other possible conditioning factors, we intended to consider the aspect. The upper Tartagal watershed is on the western flank of the anticline (see Section 2) and the regional inclination of sedimentary rock strata is westward. Therefore, it would have been reasonable to expect more landslides on the eastern erosion scarps and less on western dipping slopes. However, we were unable to observe this correlation between the slope orientation towards the east and the occurrence of landslides. As already expressed elsewhere in this paper, we did not map individual landslides events nor differentiated among different processes. This would be possible with a digital elevation model and satellite data of higher spatial resolution.

We also intended to link changes in vegetation cover with the February 2009 events through the analysis of a time series (1999-2011) of Normalized Difference Vegetation Index (NDVI) based on 30 m spatial resolution Landsat data. We observed a normal pattern of changes in vegetation coupled with the monthly variation in precipitation but we did not identify any correlation of such changes with the occurrence of the February 2009 event at Tartagal River basin and we were therefore unable to incorpórate changes in vegetation into the heuristic model. At the basin scale of this study, remotely sensed data of higher spatial resolution would have been more suitable for this purpose. It is important to point out that these observations do not imply that deforestation is unrelated to the occurrence of landslides. We must recall that the February 2009 events at Tartagal were the result of the action of both anthropic and natural factors combined (see Section 2.3).

Finally, there was good correspondence between the validation dataset and the susceptibility map (Table 5). The landslides that led to the floods in Tartagal in 2006 were restricted to the middle basin and were the results of precipitation that affected and triggered landslides in more than one basin (e.g., Quebrada Yacuy), in contrast to the February 2009 disaster that was restricted only to the Tartagal River basin.

6. Conclusions

In this work, we mapped the landslides that led to the debris flows that resulted in the February 2009 disaster in the town of Tartagal, Salta, Argentina, by using 10-m spatial resolution SPOT-4 data acquired just after the event and a standard image classification approach. The maximum likelihood classifier enabled the identification of an area of 8 km² affected by landslides and the output map showed a good correspondence with a small dataset collected in the field. The landslide surface detected by the standard image classification in combination with a minimum scarp thickness of 0.5 m resulted in a total volume of material of at least 2.2x10⁶ m³. The sensitivity of the classifier to a 10% change in the probability threshold was not significant; the resulting area and volume remained within the same order of magnitude, i.e., within the order of the millions of cubic meters. This estimate constitutes essential information for modelling purposes and hazard assessment.

We evaluated a set of conditioning factors for the occurrence of landslides including lithology, slope and curvature and observed that while the high susceptibility lithological units occupied 95% of the upper and middle catchments, slope and curvature contributed more to moderate and low susceptibilities, respectively. The reason for this is that steeper slopes were underrepresented within the river basin plus the fact that the landslides incidence within moderate to lower susceptibility areas was overestimated, because of our inability to identify source areas and discriminate different slope failure processes. Finally, the 8 km² surface affected by the 2009 events based on the classification of SPOT data, which was also an input for the heuristic model, constituted areas of high susceptibility that represented 10% of the river basin.

This study resulted in the first map of landslide susceptibility for the Tartagal upper and middle river basins through the implementation of a heuristic model that used as input lithology, slope, curvature and the 2009 landslide inventory. Zones of high and very high susceptibility occupied 25 and 17 km², respectively, and covered about 50% of the upper and middle basins together. From the remaining half, areas of low susceptibility covered an almost identical surface as very high susceptibility areas. The map showed good correspondence with the middle basin validation dataset based on the sample of landslide events of early 2006; 95% of it fell within the high and very high susceptibility areas.

The comprehensive analysis presented constitutes a starting point for more detailed future studies. Besides, the outcome data such as the volume estimation of the removed material and the digital topography are essential for further studies and hazard assessment.

Acknowledgements
We would like to thank the National Commission of Space Activities of Argentina (CONAE) and the Italian Space Agency (ASI) for the fellowship granted to P. Cardozo for the Master’s Program in Remote Sensing (CONAE and Universidad de Córdoba, Argentina) and for providing the facilities, satellite data and financial support for this study. We thank Dr. J. Viramonte and Dr. C. Peralta from the Universidad de Salta and CONICET, Argentina, for supporting the field survey of ground control points for the elaboration of the digital topography. We also thank F. Bucci from the Istituto di Ricerca per la Protezione Idrogeologica (CNIR-IRPI), Italia, for his collaboration in the identification of the 2006 landslide events. Finally, we are very grateful to D. Fernández, M. Lara, an anonymous reviewer and the editor of Andean Geology for their detailed and constructive comments and suggestions.

References

Adler, F. 2009. El desastre natural de Tartagal (provincia de Salta, Argentina). Facultad de Ciencias Exactas y Tecnología, Universidad Nacional de Tucumán, República Argentina. Revista de Ciencias Exactas e Ingeniería 31: 1-8.

Barredo, J.I.; Benavidez, A.; Hervhl, J.; van Westen, C.J. 2000. Comparing heuristic landslide hazard assessment techniques using GIS in the Tirajana basin, Gran Canaria Island, Spain. International Journal of Applied Earth Observation and Geoinformation 2 (1): 9-23.

Bathrellos, G.D.; Kalivas, D.P.; Skilodimou, H.D. 2009. GIS-based landslide susceptibility mapping models applied to natural and urban planning in Trikala, Central Greece. Estudios Geológicos 65 (1): 49-65.

Baumann, V.; Seggiaro, R.; Azcurra, D.; García. V. 2009. Inundación de detritos en la ciudad de Tartagal, provincia de Salta. Serie Contribuciones Técnicas: Peligrosidad Geológica N° 15. Servicio Geológico Minero Argentino (SEGEMAR): 17 p. Buenos Aires.

Baumann, V.; Seggiaro, R.; Azcurra, D.; García, F. 2011. Inundación de detritos en la ciudad de Tartagal. In Congreso Geológico Argentino, No. 18, Actas. Neuquén.

Baumann, V.; Villegas, D.; Seggiaro, R.; Azcurra, D. 2014. Carta de Peligrosidad Geológica 2363-I “Tartagal” Provincias de Salta. Servicio Geológico Minero Argentino (SEGEMAR) Boletin 352. Buenos Aires.

Bonorino, G.; Rivelli, F.; Bartoloni, M. 2001. Hoja geológica 2363-I, Tartagal. Provincia de Salta. Instituto de Geología y Recursos Minerales. Servicio Geológico Minero Argentino (SEGEMAR). Boletín: 331 p. Buenos Aires.

Brea, J.D.; Spalletti, P.; Amores, G. 2013. The 2006 and 2009 Landslide and Flood Events in the Tartagal River Basin (Argentina). In Landslide Science and Practice (Margottini, C.; Canuti, P.; Sassa, K.; editors). Springer: 341-347. Berlin, Heidelberg.

Carrara, A.; Cardinali, M.; Guzzetti, F.; Reichenbach, P. 1995. GIS Technology in Mapping Landslide Hazard. In Geographical Information Systems in Assessing Natural Hazards. Advances in Natural and Technological Hazards Research 5 (Carrara, A.; Guzzetti, F.; editors). Springer: 135-175. Dordrecht.

Castellanos, E.A.; Van Westen, C.J. 2008. Qualitative landslide susceptibility assessment by multicriteria analysis: A case study from San Antonio del Sur, Guantánamo, Cuba. Geomorphology 94: 453-466.

Corominas, J.; Copons, R.; Villaplana, J.M.; Altimir, J.; Amigó, J. 2003. Integrated Landslide Susceptibility Analysis and Hazard Assessment in the Principality of Andorra. Natural Hazards 30: 421-435.

Ellen, S.D.; Mark, R.K.; Wieczorek, G.F.; Ramsey, D.W.; May, T.E.; Graham, S.E.; Beukelman, G.S.; Barron, A.D.; Graymer, R.W. 1997. Map showing principal debris flow source areas in the San Francisco Bay region, California. US Geological Survey, Open File Report 97-745E: 8 p.

Fernández-Bussy, J.; Dall’Armellina, M.; López, P. 2010. El riesgo de desastres en la planificación del territorio-primer avance. Ministerio de Planificación federal, Inversión Pública y Servicios. Programa Nacional de Prevención y Reducción del Riesgo de Desastres y Desarrollo Territorial (PNUD-ARG 0/020), Programa de las Naciones Unidas para el desarrollo: 392 p. Argentina.

Gebhard, J.; Giudici, A.; Oliver Gascón, J. 1974. Geología de la comarca entre el río Juramento y arroyo Las Tortugas, provincias de Salta y Jujuy. Revista de la Asociación Geológica Argentina 39 (3): 359-375.

Grau, H.R. 2005. Dinámica de bosques en el gradiente altitudinal de las yungas argentinas. In Ecología y Manejo de los Bosques en Argentina (Goya, J.F.; Frangi, J.L.; Arturi, M.F.; editors). Editorial de la Universidad Nacional de la Plata (EDULP): 1-30.

Guzzetti, F.; Carrara, A.; Cardinalli, M.; Raichenbach, P. 1999. Landslide hazards evaluation: a review of current techniques and their applications in a multi-scale study, Central Italy. Geomorphology 31: 181-216.

Guzzetti, F.; Galli, M.; Reichenbach, P.; Ardizzone, F.; Cardinalli, M. 2006. Landslide hazard assessment in the Collazzone area, Umbria, Central Italy. Natural Hazards and Earth System Sciences 6: 115-131.

Hermanns, R.L.; Valderrama, P.; Fauque, L.; Penna, I.M.; Sepúlveda, S.; Moreiras, S.; Zavala, B. 2012. Landslides in the Andes and the need to communicate on an interandean level on landslide mapping and research. Revista de la Asociación Geológica Argentina 69 (3): 3121-327. doi: 10.5027/andgeoV45n3-3151.

Hernández, R.; Reynolds, J.; Di Salvo, A. 1996. Análisis tectosedimentario y ubicación geocronológica del Grupo Orán en el río Iruya. Boletín Informaciones Petroleras 45: 80-93. Buenos Aires.

Jakob, M. 2005. Debris flow probability and magnitude. In Debris flow hazards and related phenomena (Jakob, M.; Hungr, O.; editors). Springer-Praxis Books in Geophysical Sciences: 411-443. Chichester.

Lara, M.; Sepúlveda, S.; Celis, C.; Rebolledo, S.; Ceballos, P. 2018. Landslide susceptibility maps of Santiago city Andean foothills, Chile. Andean Geology 45 (3): 433-442. doi: 10.5027/andgeoV45n3-3151.

Lee, S.; Evangelista, D.G. 2008. Landslide Susceptibility Mapping using Probability and Statistics Models in Baguio City, Philippines. Department of Environment and Natural Resources, North Avenue, Diliman, Quezon City, Philippines: 12 p.

Martha, R.T.; Kerle, N.; Jetten, V.; van Westen, C.J.; Kumar, K.V. 2010. Characterizing spectral, spatial and morphometric properties of landslides for semiautomatic detection using object-oriented methods. Geomorphology 116: 24-36.

Mingramm, A.; Russo, A.; Pozzo, A.; Casau, L. 1979. Sierras Subandinas. Segundo simposio de Geología Regional Argentina. Academia Nacional de Ciencias 1: 95-138. Córdoba.

Moreiras, S.M. 2005. Landslide susceptibility zonation in the Río Mendoza Valley, Argentina. Geomorphology 66: 345-357.

Moreiras, S.M. 2009. Análisis Estadístico Probabilístico de las variables que condicionan la inestabilidad de las laderas en los valles de los río Las Cuevas y Mendoza. Revista de la Asociación Geológica Argentina 65 (4): 780-790.

Muñiz, J.; Hernández, V.M. 2012. Zonificación de procesos de remoción en masa en Puerto Vallarta, Jalisco, mediante combinación de análisis multicriterio y método heurístico. Revista Mexicana de Ciencias Geológicas 29: 103-114.

Nadir, A.; Chafatinos, T. 1994. Los suelos del NOA (Salta y Jujuy). Universidad Nacional de Salta. CIUNSa: 86 p.

Padula, E.; Reyes, F. 1958. Contribución al léxico estratigráfico de las Sierras Subandinas. Revista Técnica Yacimientos Petrolíferos Fiscales Bolivianos 1 (1): 9-70.

Rossi, M.; Reichenbach, P. 2016. LAND-SE: A software for statistically based landslide susceptibility zonation, version 1.0. Geoscientific Model Development 9: 3533-3543.

Ruff, M.; Czurda, K. 2008. Landslide susceptibility analysis with a heuristic approach in the Eastern Alps (Vorarlberg, Austria). Geomorphology 94: 314-324.

Sánchez, L.G. 2008. Evaluación Temporal de la Cobertura Vegetal y Estimación de la Pérdida de Suelo en la Cuenca del Río Tartagal, Depto. San Martín, Salta. Dissertation in Eng. in Natural Resources and the Environment (Unpublished). Universidad Nacional de Salta, Facultad de Natural Sciences: 163 p.

Soil Survey Staff. 1996. Keys to Soil Taxonomy, 760^th ed. United States Department of Agriculture (USDA)-Natural Resources Conservation Service: 644 p. Washington DC.

Vasantha Kumar, S.; Raja, N.; Prasad Babu, G. 2007. Extraction of topographic and morphometric features for landslide zonation- a case of study Mettupalayam Highway. In Proceed. Map World Forum. Hyderabad.

Volante, J.N.; Bianchi, A.R.; Paoli, H.; Noé, Y.; Elena, H.; Cabral, C. 2005. Análisis de la dinámica del uso del suelo agrícola del Noroeste Argentino mediante Teteledetección y Sistemas de Información Geográfica. Período: 2000-2005. Proyecto de Relevamiento de Cultivos del NOA-PRO.RE:NOA. Estación Experimental Agropecuaria Salta: 70 p.

Zieher, T.; Schneider-Muntau, B.; Mergili, M. 2017. Are real-world shallow landslides reproducible by physically-based models? Four test cases in the Laternser valley, Vorarlberg (Austria). Landslides: 1-15. doi: 10.1007/s10346-017-0840-9.