Shallow landslide susceptibility map at a regional scale (Asturias, NW Spain). A heuristic-driven approach

ABSTRACT

The study area (Asturias, NW Spain) is a mountainous region extending over 10,000 km², where shallow landslides triggered after heavy rainfall episodes are very frequent and pose a major threat and hazard to infrastructure and human populations, causing frequenteconomic losses and damage. In this regard, Shallow Landslide Susceptibility Mapping (SLSM) represents a very powerful tool for managers and agencies dealing with landslide hazards and land planning. Here we produced an SLSM based on basic geological/geomorphological maps and a heuristic approach. The regional-scale model built here has a cell resolution of 50 m and combines bedrock geology, and surface deposits mapped at a1:25k scale together with a digital slope model. The resulting map was compared with two local inventories of landslides, giving a goodness of fit of 80 and 77.5%. To better understand the results, the erroneous data have been reviewed individually and the causes of error were analysed.

KEYWORDS:

• Rainfall-induced landslides
• susceptibility map
• heuristic models
• geological/geomorphological mapping

1. Introduction

Landslides induced by prolonged and heavy rainfall are frequent and widespread in mountainous regions, posing significant hazards that can cause damage to human populations and result in economic losses. As with many other natural and geological hazards, appropriate land planning will prevent these damages and losses (Ruff & Czurda, 2008). Consequently, the development of regional-scale, landslide-susceptibility maps constitutes a very useful tool for the administrations and agencies in charge of land planning (Fell et al., 2008; Gazibara et al., 2023; van Westen et al., 2006).

The landslide susceptibility map done by Brabb et al. (1972) broke new ground in the assessment of landslide susceptibility, by proposing a methodology that ‘portrays only those geologic factors that bear directly on the susceptibility of specific areas to landsliding’. Since then, geological and geomorphological maps have typically been applied to develop derived maps aimed at defining and delineating areas prone to landslide instability. Since the 1990s, digital terrain models (DTM), thematic maps and spatial analysis using GIS software have been the main tools for building Susceptibility Map of Shallow Landslides (SLSM), and several workflows and approaches have been developed (Carrara et al., 1991, 1995; Chung et al., 1995; Chung & Fabbri, 2008; Gazibara et al., 2023; Goetz et al., 2011; Guzzetti et al., 1999; Huang et al., 2020; Menéndez-Duarte et al., 2003; Nandi & Shakoor, 2009; Pourghasemi et al., 2014; Pradhan et al., 2019; Remondo et al., 2003; Salciarini et al., 2012; Soeters & van Westen, 1996; Sujatha & Rajamanickam, 2011; van Den Eeckhaut et al., 2010; Vergari et al., 2011). Used properly, these maps constitute an essential guidance document for urban planners, helping them to optimise the effort invested in designing corrective measures.

In many reviews of landslide-susceptibility, hazard, and vulnerability assessments, e.g. van Westen et al. (2006) or Reichenbach et al. (2018), four main approaches have been identified that can be followed in the development of landslide susceptibility maps: i. probabilistic models based on landslide inventories; ii. heuristic models; iii. statistical (bivariate or multivariate statistics) models; and iv. deterministic models. Probabilistic and statistical models are based on the analysis of numerical relationships between instability factors and the distribution of landslides. Deterministic models are physical models that simulate the physical laws governing stability/instability conditions. In the heuristic or expert-driven approach, investigators rank the known instability factors according to their assumed or expected importance in causing landslides.

A heuristic approach was used in this study; expert criteria are used to infer the type and degree of hazard for each area, using either a direct mapping or indirect qualitative map combination (Barredo et al., 2000; Castellanos-Abella & van Westen, 2008; van Westen et al., 2006). Heuristic models combining geology, geomorphology, and relief information for landslide-susceptibility assessment have been developed, for instance, by Barredo et al. (2000), Dymond et al. (2006), Ruff and Czurda (2008), Demoulin and Chung (2007), Wang and Peng (2009) or van Den Eeckhaut et al. (2010). An important advantage of heuristic models is that they do not require an extensive inventory of landslides for their construction. In this regard, we cannot ignore the fact that the acquisition of landslide inventories is costly and subject to important biases (Bucci et al., 2021; Guzzetti et al., 2012; Murillo-García & Alcántara-Ayala, 2017), which in turn represents a critical advantage compared to probabilistic or statistical models (Guzzetti et al., 2012; Smith et al., 2021; van Westen et al., 2006). Nevertheless, landslide-susceptibility models are conditioned by the data available for their validation (Guzzetti et al., 2006), especially in the case of heuristic models (i.e. Ruff & Czurda, 2008), since they rely heavily on expert criteria. If some landslide inventory were available, then they could be used for comparison with heuristic models and to assess their performance (Guzzetti et al., 2006; Menéndez-Duarte et al., 2003; Remondo et al., 2003; Rossi et al., 2010).

In this research, we wanted to test whether (or not) susceptibility maps constructed using a limited number of variables and applying some common criteria over a broad area could be successful in predicting rainfall-induced landslides at a regional scale. To do so, we combined geological, geomorphological and slope maps to develop a regional-scale shallow landslide susceptibility map (SLSM). The map was validated using an inventory of shallow landslides (Menéndez Duarte, 1994) and some recent shallow landslides triggered by heavy rainstorms. Finally, the errors found were discussed in order to better understand the strengths and weaknesses of the map.

2. Study area

The region of Asturias is located in the northwest of Spain and covers an area of approximately 10,000 km². It is bordered by the Cantabrian Sea to the north and the Cantabrian Mountains to the south (Figure 1). The entire region is characterised by a mountainous relief, including hundreds of peaks higher than 2,000 m, well scattered along the approximately E-W summit axis of the mountain range. The average slope of the hillsides often exceeds 20 degrees, and rivers are deeply incised into the bedrock and typically flow through steep valleys. Only along the coastline, the relief is somewhat smoother due to the exhumation of a paleo-cut wave platform (Alvarez-Marron et al., 2008; López-Fernández et al., 2020).

Figure 1. Location of study area (Asturias region in the northwest of Spain).

Location map of the study area in NW Spain (Iberian Peninsula) and location of Spain with respect to neighboring countries in SW of the European continent.

The geology of the region is largely known and mapped, since it is one of the best-preserved foreland thrust and fold belt of the European Variscides (i.e. Alonso et al., 2009; Julivert, 1971; Marcos & Pulgar, 1982; Pérez Estaún & Bastida, 1990; Pérez-Estaún et al., 1988). The dominant basement consists of sedimentary Palaeozoic rocks, alternating siliciclastic (slate, sandstone, and quartzite) and calcareous (mostly limestones) formations. Moreover, in the east-central sector, low deformed Mesozoic rocks of varied lithologies are disposed unconformable and overlying the deformed Paleozoic basement. The region underwent significant deformation during the Variscan Orogeny. However, the modern relief is the topographic expression of the convergent regime experienced by the Iberian and Eurasian plates during the Cenozoic (Alvarez-Marron et al., 1997; Fernández et al., 2021; Pulgar et al., 1996).

The deep incision of the drainage network conditioned the present-day relief. As a result, slope dynamics are mainly controlled by mass-wasting activity and locally by torrential dynamics. Large deep mass movements have been identified and mapped in many areas, most of which are now dormant (Jiménez Sánchez, 2002; Menéndez-Duarte & Marquínez, 1996; Rodríguez-Rodríguez et al., 2018) and only occasionally are active processes (Cuervas-Mons et al., 2021; Domínguez-Cuesta et al., 2022). Furthermore, rock falls and, to a lesser extent, large-scale rock-avalanches are frequent on rocky and steep slopes (Fernández et al., 2021; Marquínez et al., 2003; Menéndez-Duarte & Marquínez, 2002; Rodríguez-Rodríguez et al., 2023).

Most of the slopes are mantled by regolith, which is composed of several lithological types and textures. On these mantled slopes, present-day mass wasting processes are dominated by soil creep and small and shallow mass movements or lanslides (mudflows, debris flows, and small slides) with a maximum size of tens of meters (Figure 2). These shallow landslides, locally called ‘argayos’, are usually associated with episodes of heavy rainfall (Domínguez Cuesta et al., 1999; Domínguez-Cuesta et al., 2021; Valenzuela et al., 2018a, 2018b, 2019). These processes are frequent and often disrupt transport infrastructures, with critical economic consequences and occasionally loss of life. In this research, the susceptibility map has been created for these small rainfall-induced landslides.

Figure 2. Some examples of rainfall-induced landslides cases in the study area. (A) Shallow landslide triggered on colluviums above slate, and (B) Shallow landslide triggered on old mass movement deposit.

Two pictures of shallow landslides show the upper scar zone, with a sub-circular shape, and the destabilized material. These mass movents are classified as earthflows or mudflows in accord with their morphology.

3. Materials and methods

3.1. Landslide susceptibility model

Landslide susceptibility is assessed and mapped by combining bedrock geology, geomorphologic (surface deposits) and slope maps. A qualitative weighting of these variables was applied to produce the GIS-based SLSM (Figure 3).

Figure 3. Workflow followed to develop the shallow landslide susceptibility map. All the cartographic information has been managed in matrix format and combined according to the operations indicated in the flowchart.

A visual depiction of the workflow showing: the primary layers of the GIS analysis: geology, surface formations, rocky areas and digital elevation model, the derived layers: slopes and weighting of geological factors and the spatial analysis functions followed to calculate the numerical value in each cell of the final map.

Slope information was derived from a 50-m digital elevation model (DEM). For bedrock geology and surface deposits (colluvial, alluvial, residual), geological and geomorphological maps at a 1:25k scale were obtained from the CTAPA (Environmental Thematic Cartography of the Principality of Asturias https://sigvisor.asturias.es/SITPA/). In this digital repository, the bedrock geological map is provided in shapefile format (ESRI vectorial format) and was elaborated based on the MAGNA project of the Geological and Mining Institute of Spain (IGME). The geomorphological map is also stored, as a separate data layer, in shapefile format. This map was produced by CTAPA on the basis of detailed fieldwork. Both cartographies have been used as independent layers, resulting in two different inputs into the model (bedrock geology and surface deposits), each with its own weighting (see below). These two GIS layers are converted to a 50 m cell size matrix format with the same resolution as the DEM. A ‘maximum area’ criterion was followed when converting from vector to matrix format; that said, when more than one single feature fell within the same cell, the value assigned was that of the feature occupying the largest area.

3.1.1. Bedrock geology map

About 150 geological formations of sedimentary rocks are mapped in 1:25k scale geological maps of Asturias. These formations were grouped into 13 lithological classes and assigned a weighting between 1 and 10, with the highest value (10) corresponding to the formations most prone to instability and the lowest to the most competent lithologies. Table 1 summarizes the criteria that were considered to assign these weightings. The aim was to construct a lithological map based on engineering properties rather than stratigraphic classification (Brabb et al., 1972; van Westen et al., 2008). The results of the 13 classes employed and the weight assigned to each one are given in Table 2. A numerical model (in matrix format) was derived from the relative weight of bedrock geology on landslide susceptibility, covering the entire study area. This derived map is the ‘Layer BG’ in Figure 3.

Table 1. Criteria used to assign weights to bedrock geology units according to their lithology and meaning – favorable or unfavorable – in the development of shallow landslides.

Criteria	Meaning
Presence of unconsolidated material	Higher propensity to produce landslides
Presence of silt, clay or slates	Higher propensity to produce landslides
Layering of different lithologies	Higher propensity to produce landslides
Unit age (Paleozoic or post-Paleozoic)	In the studied region, post-Paleozoic units suffer more instability processes
Average slope of the unit for the entire region	The lithologies more unstable have a lower average slope
Most common lithologies in each unit	It weighs the presence of more or less resistant lithologies
Strength values available for some units and lithological types*	It quantifies the lithological strength
Regional knowledge and field experience	Expert criteria completing the weight obtained by the other criteria

*Gutiérrez Claverol & Torres Alonso, 1995

Table 2. Weights are assigned to bedrock geology units, based on the criteria shown in Table I. A higher weight indicates higher favorability for the development of instabilities and involves a greater value in the susceptibility map.

Bedrock geology units	Weight
Mesozoic formations of sandstones, sands, silts and clays	10
Mesozoic formations of sandstones, sands, silts and clays with conglomerates or with limestones	10
Mesozoic limestones and marls with clays and silts	7
Mesozoic or Tertiary unconsolidated conglomerates	7
Devonian and Carboniferous formations of sandstones, clays and limestones	5
Small outcrops of igneous rocks	4
Mesozoic formations of limestones, dolomites and marls	4
Palaeozoic formations of sandstones with clays (central and eastern Asturias)	3
Palaeozoic formations of clays (western Asturias)	3
Palaeozoic formations of laminated limestone (sometimes interbedded with other lithologies)	3
Palaeozoic or Mesozoic formations of consolidated conglomerates	3
Palaeozoic formations of quartzites and sandstones	2
Palaeozoic formations of massive limestone	1

3.1.2. Map of surface deposits

The map of surface deposits available in the CTAPA repository is very detailed and contains more than 72 different types, classified according to genetic and lithological criteria. These 72 types have been simplified to distinguish only three classes based on the presence and abundance of the sandy-clay matrix. This very simplified map of surface deposits may be compared to an ‘engineering soil map’ (van Westen et al., 2008) describing all loose materials on top of the bedrock and classifying them according to their geotechnical properties. The derived map is the ‘Layer SF’ in Figure 3, whose values ⁣⁣vary according to the dominant lithology of the surface deposit (Table 3).

Table 3. Weight is assigned to cells with surficial deposits according to their sandy-clay matrix content.

Lithology of surficial deposits	Weight
Surficial deposits with abundant sandy-clay matrix	5
Surficial deposits with sandy-clay matrix	3
Surficial deposits without clay matrix	0

The combination between bedrock geology (layer BG) and surface deposits (layer SF) was made by summing the values of these two raster; in both cases, a higher value represents a cell most prone to instability, obtaining a relative susceptibility map with a value between 1 and 15. Moreover, an extra value – two extra points – was added to those cells that fall over deposits of ancient mass movements. The reason for this weighting is that field observations and previous work (Domínguez Cuesta et al., 2007; Menéndez-Duarte et al., 2011; Nyssen et al., 2003) suggest that a fundamental factor in the development of slope instability is the presence of a significant thickness of loose material covering the rock substrate. Consequently, a new matrix layer was created (Layer MM in Figure 3) with a value of 2 in cells corresponding to mass movements and zero in the rest. This model was added to the previous creating a new raster containing values ranging from 1 to 17 (‘Sumatory of Geological Factors’ in Figure 3).

3.1.3. Digital slope model

The initially continuous slope values were then grouped into five intervals and a weighting score was assigned to each interval (Table 4; Layer SL of Figure 3). This weighting multiplied the values of the matrix model resulting from the combination of geological and geomorphological layers, and a relative susceptibility map was obtained with values ranging from 0.5 to 42.5 (‘Geological Factors and Slopes’ in Figure 3).

Table 4. Slope intervals were used to build the susceptibility map and weights were given (multiplication factor) to each of these intervals.

Degree slope	Multiplication factor
0–5	0.5
5 - 15	1
15 - 25	1.5
25 - 30	2
> 30	2.5

3.1.4. Rocky areas and flat deposits

Finally, rocky areas, talus scree deposits and flat deposits (areas mapped with alluvial deposits, estuaries, and other types of deposits typical of plains) were excluded. Although rocky areas and talus scree deposits are affected by other instability processes, they are not prone to the development of shallow landslides due to the lack of regolith or soil cover. Flat deposits are not subject to instability processes and do not form part of slopes.

The rocky areas and flat deposits were extracted from the CTAPA maps and converted to matrix format (Layers RS, TS and FD in Figure 3) by assigning a zero value to the cells corresponding to them, regardless of the value that the cell may have due to the combination of the other variables.

3.1.5. Shallow landslide susceptibility map (SLSM)

The Asturias shallow landslide susceptibility map is a matrix GIS layer with continuous values from 0 to 42.5. The values displayed in the raster map give specific weight to geology, geomorphology and slope, but exclude the cells corresponding to rocky areas, talus scree deposits and flat deposits. Such continuous values are discretised into five susceptibility levels: Very Low susceptibility (cells with 0-value corresponding to rocky areas and flat deposits); and Low, Medium, High, and Very High susceptibility, which have been defined after reclassifying the values ranging from 0.5 to 42.5 into four equal-area intervals. According to the specific weight of bedrock geology (BG), surface formations (SF), mass movements (MM) and slope (SL) and discounting rocky slopes (RS), talus scree (TS) and flat areas (FA), susceptibility model of shallow landslides might be expressed by the equation: SLSM = ((BG+SF+MM) x SL) x (RS,TS,FD).

3.1.6. Validation of the susceptibility map

The SLSM was validated using compiled shallow landslide data from two sources: i) shallow landslides reported in the local press and triggered after two heavy rainstorms episodes (June 2010 and November 2003) that were analysed in the field (Figures 4–6); and ii) a cartographic inventory of 93 scars and shallow landslide deposits mapped in a 200 km² zone of the study area (Figures 7 and 8).

Figure 4. Rainfall-induced landslide triggered by November 2003 rainfall: (A) Buscabreiro and (B) Zelúan. Map fragments show the location of the landslide on the susceptibility map. In both cases, the location of the rainfall-induced landslide on cells classes as Very High Susceptibility is related to the presence of surficial deposits (colluvium deposits).

Four images are included: two fragments of the susceptibility map, from two different areas, with the positions of the landslide starting points overlaid onto them, and two pictures of the slides showing the rupture zones and the unstable debris.

Figure 5. Rainfall-induced landslide triggered by November 2003 rainfall (Cabañaquinta) (A) Location of the landslide on the bedrock geology map, and (B) location of the landslide on the susceptibility model; the landslide falls on a cell classed as Low Susceptibility due to the proximity of a limestone level (see text).

Two images are included: one showing a fragment of the geologic map and the other showing an fragment of the susceptibility map. These maps were overlaid with the locations of shallow landslide initiation points.

Figure 6. Rainfall-induced landslide triggered by June 2010 rainfall (Viescas) (A) Location of the landslides on the bedrock geology map, (B) location of the landslides on the susceptibility model, (C) unconsolidated conglomerates and sands in which the landslides were triggered (probably tertiary conglomerates and sands, see text) and (D) photograph showing one of the landslides.

It is composed of four images: a fragment of the bedrock geology map, a fragment of the susceptibility map, a photograph showing the type of material involved in the slide and a photograph of one of the landslides displaying one of the landslides identified in the field. On the two maps the shallow landslide initiation points have been superimposed.

In the first case, rainfall-induced landslides reported in the regional press (a total of 10 localities) were visited. Seven landslides occurred immediately after the heavy rains of 15^th-17^th November 2003 and three others after the heavy rains of 10^th-11^th of June 2010. All of these landslides were georeferenced, photographically reported and analysed during the fieldwork.

In the second case, the cartography of 93 shallow landslides, mapped in detail in the context of a geomorphological study of the Somiedo Natural Park (Menéndez Duarte, 1994), was compared with the SLSM. The material involved in these landslides and their respective morphologies suggests that they were also triggered after heavy rainfall episodes, even though the exact date of these landslides is not known.

All landslides were cross-checked with the SLSM and classified according to their respective susceptibility level. Although the accuracy of the map is high in both cases (see results), some errors can also be observed. In order to better understand the results, the erroneous data were examined individually and the causes of the errors were analysed.

4. Results

4.1. The shallow landslides susceptibility map (SLSM)

When viewed globally, the susceptibility map shows two main areas of high susceptibility: (1) the north-central zone of the region, resulting from the predominance of poorly consolidated Mesozoic and Cenozoic rock formations (sands, clays and marls); and (2) the south-central zone, where carboniferous shale and sandstone predominate in combination with high average slopes and abundant old mass movements. In the rest of the region, the most susceptible areas basically follow the main traits of the geological structure, in good agreement with the alternation of more or less resistant materials. On a regional scale, the zones of rocky areas and rock screes (in some cases quite extensive) are also well represented on the map, as are the alluvial plains associated with the major rivers. Finally, a strip of low susceptibility parallel to the coastline is associated with the paleo-cut wave platform. By zooming in on the map, we can appreciate that the shallow landslide susceptibility is mainly controlled by the presence/absence and type of surface deposits, as well as by the variability of local slopes.

4.2. Comparison of susceptibility map with shallow landslides triggered by heavy rainfall episodes

The comparison of the available field information on shallow landslides triggered by heavy rains (November 2003 and June 2010) with the SLSM yields the results shown in Table 5. Seven of the ten documented landslides fall on cells classified as Very High susceptibility. Concerning the remaining three, one occurred in unclassified areas of the map (urban areas), while the other two occurred in zones classified as low and medium susceptibility.

Table 5. Summary and brief description of rainfall-induced landslides visited after two periods of heavy rain in the study area.

	Locality	Susceptibility class	Brief description
November 2003	Buscabreiro (Tineo)	Very high	Shallow landslide triggered on colluvium over slates (demolition of a house and two people dead)
	Serandias (Boal)	Very high	Many shallow landslides triggered on colluviums over slates (demolition of a garage and agricultural losses)
	Cabañaquinta (Aller)	Low	Shallow landslide triggered on carboniferous lutites and sandstones (demolition of a house)
	Oyanco (Aller)	Very high	Shallow landslide triggered on carboniferous lutites and sandstones (a broken road)
	La Felguera (Langreo)	Very high	Shallow landslide triggered on carboniferous lutites and sandstones (a damaged home)
	Carceda (Cudillero)	Very high	Shallow landslide triggered on colluvium over sandstones (a broken road)
	Tuiza (Lena)	Very high	Shallow landslide triggered on colluvium over sandstones and slates (a broken road)
June 2010	Raíces (Castrillón)	Unclassified (urban soil)	Shallow landslide triggered on an old mass movement deposit (a damaged home)
	Viescas (Illas)	Low and Medium	Many shallow landslides triggered on not consolidated sands and conglometates (a broken road and agricultural losses)
	Zelúan (Gozón)	Very high	Many shallow landslides triggered on colluviums above limestones and marlstones (agricultural losses and a home in danger)

Several of the satisfactorily assigned cases corresponded to the presence of surface deposits. For instance, this is the case of the Buscabreiro and Zeluán landslides (Figure 4), which fell over pixels classified as ‘Very high susceptibility’ in the SLSM and which correlate well with the presence of colluvium deposits. Moreover, in Table 5, we can notice that the lithological groups with the second highest registration of shallow landslides are carboniferous shales and sandstone, which are particularly frequent in the south-central zone of the region.

However, two of the landslides reported in the press took place in areas classified as Low (Cabañaquinta), and Low to Medium susceptibility (Viescas) (Figures 5 and 6). In the case of Cabañaquinta (Figure 5), the instability was triggered over shales, but only at the contact with limestones. Therefore, it is very likely the SLSM misclassified this cell as Low susceptibility due to the contact with limestones and pixel inaccuracies following the conversion of the geological map from vector to matrix format. Thus, the error can be considered somehow as a ‘boundary error’ associated with the vector conversion to matrix format.

In the case of Viescas (Figure 6), three landslides occurred over model cells classified as High, Medium and Low susceptibility. According to the available geological map, the landslides occurred over a Devonian formation composed of shales and dolomites. However, during the fieldwork, we observed that the landslides actually occurred over unconsolidated conglomerates and colluvium derived from these conglomerates. These conglomerates are named as Jurassic conglomerates on the original geological map, but: (1) the areal delimitation of these conglomerates in the map is smaller than the actual extent observed in the field; and (2) Jurassic conglomerates in this region are normally consolidated material, but the conglomerates observed in the field are unconsolidated materials; we believe they are probably younger, most likely similar to the Tertiary conglomerates and sands found in the vicinity. Thus, the errors detected appear to be associated with some deficiencies and inaccuracies in the original bedrock geology map.

4.3. Comparison of susceptibility map with the shallow landslide inventory

The comparison of the compiled cartographic inventory of shallow landslides with the susceptibility map yields the results shown in Table 6: 77.5% of the previously mapped shallow landslides coincide with areas classified as Very High and High susceptibility (more than 60% as Very High) in the map. We find these results very encouraging.

Table 6. Shallow landslides in each class of susceptibility. The total number of landslides is 93 (cartographic inventory in an area of 200 km²).

Num. of shallow mass movements	Percentage	Susceptibility class
58	62.4	Very high
14	15.1	High
8	8.6	Medium
4	4.3	Low
9	9.6	Flat deposits and rocky areas

Nevertheless, they remain some poorly classified landslides. The highest percentage of them corresponds to landslides located over ‘non-susceptible’ areas: flat deposits or rocky areas. More specifically, of the nine shallow landslides that fall over these cells (Table 6), five are over flat areas and four over rocky areas. A careful analysis of these cases of misclassification suggests that they might be considered as `boundary errors’ created after converting the input thematic layers from vector to matrix format. Two examples of these poorly classified landslides are shown in Figure 7 (landslides over cells assigned to flat surficial deposits) and Figure 8 (landslides over cells assigned to rocky areas).

Figure 7. (A) Shallow landslides, belonging to the landslide cartographic inventory, superimposed on the geomorphologycal map (surficial deposits and rocky areas map), and (B) the same landslides superimposed on the Susceptibility Model (matrix format 50×50 m cell). In figure A, it can be observed that some landslides are located over deposits of mass movement, very close to the contact between these deposits and alluvial plain. In the conversion to matrix format, the corresponding pixels have been assigned to the alluvial plain (flat surface deposits).

It consists of two images: the first is a fragment of the map of surface formations and rocky areas and the second is a fragment of the susceptibility map. On the two maps the landslide initiation points have been superimposed.

Figure 8. (A) Shallow landslides, belonging to the landslide cartographic inventory, superimposed on geomorphologycal map (surficial formation and rocky areas mapping), and (B) the same landslides superimposed on the Susceptibility Model (matrix format 50×50 m cell). Two shallow landslides were mapped next to the rocky areas (one of them on colluvium deposits) but in the matrix model, they are in a cell assigned to a rocky area.

5. Discussion

The SLSM built here gives an overall goodness of fit of 77.5% (Table 6), which, according to Guzzetti et al. (2006), corresponds to a model that provides results between acceptable and very satisfactory. Indeed, our SLSM has successfully classified nine of the twelve recent cases of shallow post-rainfall landslides reported in the press (Table 6), which means a goodness of fit close to 80%. In addition, as we pointed out above, most of the errors or misclassifications detected cannot be attributed to the criteria considered in the map construction (Table 1), nor to the weighting of their specific coefficients (Tables 2–4), but rather to inaccuracies related to the scale of the model or to an error in the original cartographic information.

Concerning error analysis, the main source of error is related to the generalisation of spatial information when converting vector data to a 50-metre cell size matrix format (Figures 5, 7, 8). Some of these errors could be rectified by using higher-resolution digital models of the study area. However, there are also questions about the use of a model constructed at a regional scale for the study of small-scale processes individually. The model presented here can be considered as a medium scale (1:100k to 1:25k). Following Cascini (2008), it should be used on an area of 1,000 to 10,000 km², to develop landslide inventory and susceptibility mapping for regional development or very large-scale engineering projects and for preliminary level hazard zoning of local areas. A similar idea is expressed by Fell et al. (2008), who discussed the use of maps for site-specific decisions.

Another type of error was found in the quality of the available thematic maps. This is the case of Viescas (Figure 6), where the bedrock geology map contains errors (incorrect cartography and incorrect lithological assignment) and there were some misclassifications. According to Guzzetti et al. (2006), determining the errors associated with the geomorphologic, geological and other thematic information is not a trivial task.

Finally, the assignment of incorrect weights to some of the lithological classes is one of the problems commonly described for heuristic approaches (van Den Eeckhaut et al., 2010). This problem may also affect the model presented here, at least partially, in the case of Cabañaquinta (Figure 5). There is a boundary error in this locality, but the adjacent pixels do not belong to the Very high Susceptibility category either, although these are very frequent in the area. There is probably an incorrect interpretation of the relationship between slope angle and bedrock lithology.

Nevertheless, these misclassifications and errors could be expected a priori, as any prediction of landslide susceptibility has some level of uncertainty. Sources of uncertainty include (1) errors and incompleteness in the landslide and thematic information available for the analysis; (2) an imperfect understanding of landslide processes and their geographical and temporal evolution; (3) limitations in the techniques used to determine the susceptibility; and (4) the inherent natural variability of landslide phenomena (Ardizzone et al., 2002; Carrara et al., 1992; Guzzetti et al., 2006). The important point is that although our map is still subject to some of these errors, we considered that overall our map has a high degree of accuracy. We believe that the workflow followed here provides a reliable alternative to deterministic or statistical/probability susceptibility mapping, which requires extensive detailed and high-resolution data, which are usually not available for the whole territory over large areas.

The susceptibility map presented here is based on very few variables and a model obtained using a reduced number of explanatory variables. This model is less expensive and, thus, superior to a model that uses a very large number of variables (Guzzetti et al., 2006). Usually, from of all the possible factors controlling landslides, only a few can be cost-effectively acquired over extensive areas.

Since the early days of landslide susceptibility mapping, lithology and slope have been two common factors used in models, as seen in pioneers such as Brabb et al. (1972) and recent works (Clerici et al., 2002; Domínguez Cuesta et al., 2007; Kamp et al., 2008; Pradhan et al., 2019; van Den Eeckhaut et al., 2010). In this case, mapping the bedrock geology and geomorphology provides a good understanding of the causal factors of landslides in the study area. In addition to bedrock geology, we gave strong importance to the presence/absence of surface deposits (colluvial, alluvial, etc.) and the availability of detailed geomorphological mapping has been indispensable for the construction of the model. This conclusion was also supported by van Westen et al. (2003, 2006, 2008), Guzzetti et al. (2012) or Reichenbach et al. (2018), who considered that detailed geomorphological mapping is essential to obtain reliable landslide hazard maps for large areas.

6. Conclusions

• - The results obtained in this work support the hypothesis that heuristic models, accomplished with a limited number of variables and using common criteria for a large area, are useful tools for assessing the level of hazard associated with rainfall-induced landslides and for producing medium-scale landslide susceptibility maps.

• - In our opinion, the main advantage of the methodology used here is that it can be applied to large areas, in the absence of region-wide landslide inventory maps. These inventories are necessary for the development of statistical and probabilistic susceptibility maps, and it is not easy to obtain them of the right quality for large areas.

• - The main limitations of the heuristic-driven models are probably related to their strong dependence on the quality of the expert criteria, i.e. the quality of the resulting susceptibility map is strongly dependent on the experience of the earth-scientist. This makes validation of the final map essential, although validation is also required for all other possible approaches to landslide susceptibility modelling.

• - Finally, we would like to emphasize that the high level of accuracy of our map is mainly due to the previous existence of detailed geological and geomorphological maps of good quality. This conclusion highlights the importance of having cartographic databases of high quality and up-to-date and of using these databases in territorial management strategies.

Software

The processing of the geographic data (vector thematic layers and digital elevation model), the spatial analysis (data reclassification, vector to raster conversion and layer overlapping) and the final output of the Shallow Landslide Susceptibility Map were carried out in ESRI’s ArcGIS 10.2.2.

Disclosure statement

No potential conflict of interest was reported by the authors.

Data availability statement

The data that support the findings of this study are available from the corresponding author, R. Menéndez-Duarte, upon reasonable request.

Additional information

Funding

This work was supported by 1) Project CN-01-087-B3 ‘Mapping Hazards in the Principality of Asturias’ supported by the General Direction of Planning and Urban Development (Principality of Asturias) and 2) Project PDI2021-126357NB-100 Spanish Ministry of Science and Innovation.