Introduction 

Forests provide timber and non-timber products and a wide range of services. Most forests are multifunctional, to a certain extent, simultaneously fulfilling ecological, economic and social functions. However, a clear functional specialization of forests also exists and can thus determine focused management strategies. This is the case of protection forests.

The protective role of forests is due to the ability of stands to control and reduce the occurrence and impact of natural hazards such as landslides, soil erosion, surface runoff, debris flow and avalanches ([12], [32]). Forest protection may be classified according to indirect or direct effects ([53]). Direct protection implies that the forest directly defends people, their goods and infrastructures against the impact of natural hazards. Indirect protection is related to other natural resources, such as protection of soil against erosion or protection of the quantity and quality of water resources. Forests offering direct protection also provide indirect protection, whereas the opposite is not always true. Berger & Rey ([9]) introduce also the terms “active protection” (which is possible in the departure zone of avalanches, floods or erosion) and “passive protection” in the transition and stopping zones. Protection forests require management decision other than, for example, those of productive forests aiming at preserving and enhancing the protective effectiveness of the stands ([11], [28]). However, in many cases, protective stands are also harvested for wood production and must be managed weighing the relative importance of each function ([20]). Thus assessing and mapping the degree of protection against a certain hazard guaranteed by a forest are key tools for planning silvicultural options at stand level. This paper had this goal referring to protection against landslides. We did not distinguish between direct and indirect protection, the scope of the work being to evaluate the protective potential of each stands apart from the object to protect.

The relationships between vegetation cover and mass movement are diverse and multifaceted, as discussed by Sidle et al. ([68]), Sidle & Ochiai ([69]) and Marston ([48]). Forest cover enhances slope stability through its effect on both hydrology and soil strength. Through interception, evapo-transpiration and enhancing soil permeability, forests improve the hydrological characteristics of the soil. Tree root systems enhance the shearing strength of the soil, enabling it to resist landslides and erosion. Among the preparatory causes, slope failures may take place under modification of forest stand conditions by human activities, such as road building and wood harvesting ([69]). Landslides which occur some years after deforestation are often attributed to decay and loss of strength in root networks of a former forest cover ([33], [1]).

A landslide can be defined as the movement of rock, earth or debris down a slope under the force of gravity ([24]). Despite this simple definition, a landslide is a complex phenomenon characterized by five fundamental mechanisms (fall, topple, slide, lateral spread, flow) and their combinations ([23]). Varnes ([75]) defines landslide hazard as the probability of a landslide occurring within a specified period and given area. This definition includes two features: geographical location (spatial probability) and the occurrence of the events (temporal probability). Spatial probability, as represented by landslide susceptibility, is the likelihood of a landslide occurrence in a given location (or land unit) conditioned to its characteristics ([19], [35]). Landslide susceptibility maps are very useful for planners for selecting suitable locations for developmental activities in mountain regions ([65]).

A variety of approaches, including heuristic, statistical and deterministic ones, were developed to assess landslide susceptibility. The heuristic approaches take into account the relation between the location of previous landslide events and geo-morphological and geo-lithological settings ([2], [5], [77], [10]). This approach strongly depends on the expert experience: each geo-morphological and geo-lithological factor is weighted according to its assumed or expected influence in causing mass movements. A final map is obtained by overlaying the thematic layers relative to each factor.

The deterministic approach is based on the hydrogeo-technical models of slope stability. A variety of frameworks are available in literature referring to different physical models. For instance, many authors have used the infinite slope stability paradigm ([52], [64]), hydrological models with slope stability component ([50], [42]), more complex dynamic landscape evolution models ([21]), and shallow subsurface flow and landsliding models ([73]) to assess landslide susceptibility. The primary drawback of deterministic approach is that it requires detailed information on many environmental conditions within a landscape at a high spatial resolution.

Statistical approaches to assess landslide susceptibility involve the determination of the combinations of variables that have led to landslides in the past. Thorough reviews of the approach are given by Guzzetti et al. ([37]) and Brenning ([14]). Both bivariate and multivariate statistical approaches were exploited. Bivariate analysis compares each predisposing factor with a landslide inventory map, ignoring any significant intercorrelation among the factors. Under bivariate statistical analysis, each causal factor layer is combined with the landslide distribution map and weighting values are calculated according to landslide densi- ties. One of the common methods applied is the weight-of-evidence model ([43], [63], [62], [59]).

Multivariate analysis builds statistical models by simultaneously considering the influences of all factors within a defined land unit on landslide susceptibility. For this analysis, several authors have successfully used logistic regression modeling ([4], [78], [46], [18], [27], [55], [61]).

This paper compares the results of bivariate and multivariate statistical approaches in order to introduce a framework for estimating the degree of protection offered by forest at the stand level as a decision-support tool for forest management and planning with a protective role against landslide. Assuming that the forest cover generally enhances slope stability, we propose an approach to asses forest protection by considering landslide susceptibility scenarios derived from various landslide predictive factors dataset.

The methodological approach, which is illustrated with references to a case study in the Italian Alps, has the following goals: (i) to evaluate the influence of different factors on the model prediction; (ii) to quantify, in relative terms, the protection against landslide occurrence provided by forests; (iii) to determine the spatial extent of forest that have a preeminent protective function.

  Material and methods 

Models established for distinct functional units of landslides involve different patterns for landslide susceptibility assessment ([76]). Magliuolo et al. ([47]) demonstrated that consideration of landslide detachment zones rather than landslide body (the transport zone and the deposition zone) for statistical analysis conveys more meaningful results. Therefore, only the area of landslide detachment zone was considered for the statistical analysis in the present study.

The bivariate (weight-of-evidence) and multivariate (logistic regression) approaches were exploited to examine the relationship between landslide occurrence and geo-environmental factors, whose maps are typically available for large areas. For each statistical approach, two models were established, with the second one including forest cover as a predicting factor. The predictive potentials of the established models were evaluated using the Receiver Operating Characteristics (ROC) curve ([8], [31], [71]).

After identifying the best statistical approaches, the landslides susceptibility maps produced by either models that exclude (model type 1) and the models including (model type 2) forest cover were compared. By subtracting the model type 1 from the model type 2 a coefficient map was obtained by scaling between 0 and 1 this difference. Our interpretation is that this coefficient quantifies the relative effect of the forest cover on slope stability.

Weight-of-Evidence (WofE)

The WofE technique is based on a statistical Bayesian bivariate approach. Bayes’ theorem expresses the conditional probability of finding supporting evidence s in the study area given the class B_i of the B predicting factor, as B_i variable represents the presence of a factor that predicts the potential evidence, B_i^ is the absence of a factor that predicts potential evidence, s is the presence of the event and s^ is the absence of the event. The above-mentioned conditional probability can be written as follows (eqn. 1):

where P(B_i|s) is the conditional probability of B_i given s; P(s) is the prior probability that s exists in the study area; P(B_i) is the prior probability of finding the class B_i in the study area.

Similarly, the conditional probability of finding s in areas not covered by class B_i of the factor B is as follows (eqn. 2):

where P(B_i^|s) is the conditional probability of not finding the class B_i given s; P(s) is the prior probability that s exists within the study area; P(B_i^) is the prior probability that class B_i does not exist within the study area.

These conditional probabilities are incorporated into two logarithmic equations to determine positive (eqn. 1) and negative (eqn. 2) weights, from which a weighted contrast (eqn. 3) can be determined for each predictive factor’s class, which reflects the overall spatial association between landslides and the factor ([13] - eqn. 3, eqn. 4, eqn. 5):

Landslide susceptibility is calculated from the combined probabilities associated with the different components of the model ([6] - eqn. 6):

where O_f=P_f/(1-P_f) is the prior odds of a landslide in the study area, and P_f=A_f/A_t where A_f is that portion of the study area affected by landslides, and A_t is the total study area. The WofE approach was here applied using the spatial data modeller compiled in the Python language with ArcGIS^® 9.3 software ([30]) by Sawatzky et al. ([66]).

Logistic Regression (LR)

Logistic regression is a common multivariate approach that: (i) can be used to predict a result measured by a binary variable, such as the absence or presence of landslides, based on a set of one or more independent variables (predictive factors); (ii) does not require the variables to be normally distributed; and (iii) the independent variables can be either nonlinear, continuous, categorical or a combination of both ([39], [45]).

The conditional probability that an event will occur is given as follows (eqn. 7):

and the logistic model can be expressed in its simplest form as (eqn. 8):

in which z can be expressed as a linear combination (eqn. 9):

where β₀ is the intercept of the model, n is the number of variables, X_i (i = 1, 2, …, n) are the independent variables, and β_i (i = 1, 2, …, n) are the values associated with each of the independent variables. P varies between 0 and 1 on an S-shaped curve as z varies from -∞ to +∞ ([39]).

The logistic multiple regression algorithm within the free statistical programming language “R” ([60]) was used in this study to obtain the coefficients for the logistic multiple regression model according to a forward stepwise method.

The results were then imported in the ArcGIS software for mapping landslide susceptibility with a probability of occurrence ranging from 0 to 1.

Case study

The Alpine region is characterized by a strong link between local communities and forests ([56]) and by particular attention to forest multifunctionality ([36]). This is particularly relevant in the Italian Alps where all forests have been affected by humans to some extent, either through direct periodic wood harvesting or by other forms of land-use (e.g., grazing). Nevertheless, in the last decades there has been a noticeable reduction of the anthropogenic disturbances and, as a consequence, many forest stands have developed naturally even if their composition and structure still reflect past human activity ([54]).

The study area is inside the Veneto region (Italy) within the geographic zone of the Dolomites (south-eastern European Alps), dominated by mountain ridges rising up to 3000 m a.s.l., with their uppermost portion consisting mainly of carbonate rocks (limestone and dolostone) that form isolated peaks or high standing plateaus. The case study refers to the Cadore-Longaronese-Zoldo administrative district, covering approximately 32 000 ha and contains 7 municipalities consisting of small and sparsely populated villages. The administrative district covers a latitude of 46° 11′ - 46° 26′ N, and a longitude of 12° 05′ - 12° 30′ E (Fig. 1).

Fig. 1 - Map of the case study area and landslide inventory map (a). The locations of the area within the Veneto region (b) and within Italy (c) are shown on the right.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Artificial surfaces such as urban and built up areas cover a small percentage of the study area (2%). Approximately, 10% is used for agriculture and consists mostly of meadows and alpine pastures; screes and rocks outcrops cover a similar surface (11%). The most dominant land cover type is forest (77%).

The survey protocol was based on the following steps:

1. delineation of the landslide detachment areas from the IFFI Project (Italian Landslides Inventory);
2. selection of the predictive factors and acquisition of the corresponding GIS layers;
3. WofE modeling and evaluation of the predictive factor classes via W⁺and C;
4. LR analysis by forward stepwise regression method;
5. validation and assessment of the model’s predictive potential using the ROC curve;
6. comparison of the modeling results, i.e., the susceptibility map generated using only geo-environmental variables as predictive factors (model type 1) vs. the susceptibility map generated including also forest cover as a predictive factor (model type 2).

Landslide database

The IFFI project constitutes the Italian national landslide inventory, accomplished from 1999 to 2004. In the 2005-2006 biennium the inventory was updated based on field surveys, collection of the historical and archive data, and aerial photo interpretation ([41]). Classification of landslides was made referring to Varnes ([75], modified by [23]), partly modified to meet the practical needs of surveying and mapping the landslides. Additional types of movement were introduced: sinkholes, deep seated gravitational slope deformation, areas affected by numerous rockfalls/topples, areas affected by numerous sinkholes and areas affected by numerous shallow landslides. The last three classes were added in order to classify those landslides which are limited in size, recurrent and referable to same type of movement, which affect large sectors of slopes in Italy ([41]).

A total of 120 active landslides were identified in the study area, covering close to 350 ha and accounting for 1.1% of the territory (Fig. 1). Although the exact date of occurrence for most of the considered landslide is unknown, they are all referred to the last two decades. According to the landslide classification by Varnes partly modified, 70 % of the inventoried landslides account for three types of mass movement: rapid flow, rotational/translational slides and complex slope movements (Tab. 1). Rock falls, that can occur from Dolomite peaks, were not considered in this study.

The landslide detachment zone is where failure occurs. Generally, landslide detachment surface is located in the upper part of steep slope within the boundaries of crown portion of the slide, beyond the main scrap and the active scrap itself (Fig. 2).

Fig. 2 - A rotational landslide within the case study area (Pontesei zone, municipality of Forno di Zoldo, province of Belluno, Italy). The detachment zone is delineated.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Landslide detachment identification and delineation were drawn from the IFFI inventory through an on-screen monoscopic-visual interpretation technique of digital aerial photos taken in 2006-2007. The detachment area of each landslide was delineated in a vector-type spatial database using the ARCGIS 9.3 software package. The geo-referenced point located at the highest point of the landslide crown was selected as starting/ending point to map each detachment area. When the landslide detachment zone was hidden or obscured by tree cover, the identification was carried out by considering the size of the features, the contrast between landslides and surrounding areas, and the morphological expression (Tab. 2).

Landslide predictive factors

To apply the statistical approaches, a spatial database was created. According to literature ([16], [70], [37]) and to the data available over the entire study area, six potential predictive factors were considered: aspect, slope, lithology, permeability, soil type, and forest cover. Influential climatic factors, such as rainfall, were not taken into account because their variation is negligible within this area.

The DEM resolution is the basis for the mapping unit in landslide susceptibility assessments ([4], [29]). Lee et al. ([44]) undertook a comparative study of landslide susceptibility assessments at different spatial resolutions: relationship between the landslide distribution and the predictive factors generated by DEM (constructed from a 1:5.000 scale topographic map) and vector data were evaluated with spatial resolutions of 5, 10, 30, 100 and 200 m. The authors determined that a 30 m pixel size was the maximum advised to obtain suitable predictive capacity. For this study, it was decided to use the best precision level of the available data (20 m) as the basic mapping unit. By using a small pixel size, we were reasonably confident that the sampled parameters are constant across the pixel size. Applying a cell size of 20 m resulted in 807 255 pixels, of which 4 014 were classified as landslide detachment area.

The spatial distribution of the potential predictive factors was crossed with those of the landslides and detachment areas to calculate the density values for each class of the predictive factors (Fig. 3). The density was calculated according to Ayalew & Yamagishi ([4]): the ratio between the area occupied by the landslide pixels in a given class of a certain factor and the total pixels in that class was first determined; this was repeated for each of the classes available for that factor; these ratios were then added, and each ratio was divided by the total sum to obtain the density with respect to each factor.

Fig. 3 - Density of all landslides (black) and of the landslide detachment zones (grey) with respect to the predictive factors: aspect (a), slope (b), lithology (c), permeability (d), soil types (e) and forest cover (f). Histograms showing landslide densities as described in the text. Sub-divisions of the X-axis represent classes of the six landslide-predictive parameters. The definitions of the acronyms for the lithology classes, permeability classes and soil types are reported in Appendix 1: Rock Units (Tab. S1), Perm Classes (Tab. S2), Soil Types (Tab. S3).

  Enlarge/Shrink   Download   Full Width  Open in Viewer

DEM-derived factors: aspect and slope

The aspect and slope layers were derived from DEM (resolution of 20×20 m) based on digital topographic maps at scale 1:5.000 of the study area revised in 2004. The DEM was created by using an interpolation method based upon the ANUDEM program developed by Hutchinson ([40]), which is included in the ArcGIS 9.3 software package.

These predictive factors were treated as discrete variables under the WofE analysis (categorized according to their different classes), while they were considered continuous variables as input for LR analysis.

The aspect indicates the influence of solar radiation, exposition to preferential winds and local weather conditions on landslide distributions. Aspect discontinuities may also control the occurrence of landslides ([72], [51], [38]). The aspect was classified according to direction: north (0°-22.5°; 337.5°-360°), northeast (22.5°-67.5°), east (67.5°-112.5°), southeast (112.5°-157.5°), south (157.5°-202.5°), southwest (202.5°-247.5°), west (247.5°-292.5°) and northwest (292.5°-337.5). No flat class was found in the study area. Most events occurred on southern-oriented slopes (19% - Fig. 3). For LR analysis the aspect was treated as a continuous variable after cosine transformation.

Slope gradient is one of the factors most closely associated with landslides ([37], [57], [51], [55]). The slope layer was derived from the DEM using the average maximum technique ([15]) included in ArcGIS 9.3. The slope value was calculated using the maximum change in elevation over the distance between each pixel and its eight neighbors. We defined nine categories of slope percentage slicing the values into equal intervals and computed the corresponding percentage of landslide occurrence (Fig. 3). In the study area, most landslides occur on slopes of 40 to 100%.

Other geo-environmental factors

Referring to the case study, three factors hypothesized to act as other geo-environmental predictors of landslide susceptibility were extracted from the regional spatial datasets available: lithology, permeability and soil types.

Different lithological types are affected by different landslide susceptibilities ([26], [78], [55]). The lithology layer, rasterized from a 1:250.000-scale geological map of the Veneto region by Dal Piaz ([25]), consists of 19 rock units (see Tab. S1 in Appendix 1 for a short description). Most landslides occur in arenitic rocks belonging to Auronzo, Fernazza and La Valle formations (Ar - Fig. 3).

It is widely recognized that geological parameters greatly influence the occurrence of landslides, because lithological and structural variations lead to a difference in strength and permeability of rocks and soils, runoff and groundwater movement. In many cases, landslides are due to large changes in the pore-water pressure due to the rise or fall of the groundwater level. The Department of Geology of the Veneto region compiled the permeability data through the digital integration of the geological map at a scale of 1:250.000 (see Tab. S2 in Appendix 1 for the classes of this predictive factor). The permeability map was prepared to show the relative permeability of lithological units. The landslide densities computed for rocks with low primary and secondary permeability were 31% and 19%, respectively. The landslide density was also high in rocks classified as having medium permeability (25%). When assessing the density of these permeability classes with respect to the detachment zones, the maximum density value was confirmed for the rocks with low primary permeability, while the density increased to 26% for rocks with low secondary permeability and decreased to 18% for rocks with medium primary permeability (Fig. 3).

Sidle et al. ([68]) observed that soil properties such as particle size and pore distribution in the soil matrix can significantly influence the slope stability. These properties control the rate of water movement and the capacity of the soil to hold water. In addition, finer soils tend to hold higher volumes of water under unsaturated conditions than coarse-textured soils. The digital soil information was obtained from a 1:250.000-scale soil map of the Veneto region ([3]). The soil types identified in the study area are reported in the Tab. S3 in Appendix 1. In the study area, the soil type mostly affected by landslides was the DB3 type (Fig. 3).

Forest cover

In order to clarify the effects of forest stands on the occurrence of landslides the local forest cover map, available at a scale of 1:10.000 ([3]), was used for the analysis. Each landslide detachment area was assigned in the field to either forest or open land and checked in the field. About 77% of the study area is covered by forest, of which approximately 0.5% is included within the landslide detachment areas. The forest stands are primarily dominated by Larix decidua L., Fagus sylvatica L. and Pinus sylvestris L.

The protective function of the forest against landslide is immediately apparent from Fig. 3: the density of landslides in non-forest areas is more than twice that in forest areas.

Preparation of the calibration and validation data, the ROC curve

The tested statistical approaches were applied by two stages: model calibration and model validation. Based on the random partition approach, the 90% of the original dataset was selected as the calibration subset used for building the model and the remaining 10% was designed as validation subset used for independent verification of the model’s predictive performance via the ROC curve. The technique using a ROC curve has been frequently adopted to evaluate and compare models ([78], [35], [55], [17]). ROC curves can be prepared: (i) using the same landslide information used to construct the classification model (in this case, the ROC curve measures the degree of model fit); and (ii) for independent landslide information, in which case they measure the classification prediction skill ([74]), as in the present case study. The curve plots the sensitivity (true-positive) vs. 1-specificity (false-positive) to indicate the ability of a model to correctly discriminate between the positive and negative observations in the validation space. High sensitivity often implies low specificity and vice versa. For example, a conservative model would have high sensitivity (a high number of correct predictions) but low specificity (a high number of false positives). In a standard ROC plot, the area under the curve (AUC) is a quantitative measure of the model performance: the higher is this value, the better is the model fit ([58]). The statistical range from 0.5 (random prediction, represented by the diagonal straight line) to 1 (perfect prediction) can be used to compare models. In contrast to success and prediction rate curves, ROC curves are not sensitive to prevalence, i.e., considerable difference between landslide free and landslide-affected mapping units ([8]): therefore, ROC curves are considered a more appropriate assessment and validation tool ([74]).

  Results 

The values of the weights and the contrasts and their standard deviations for each predictive factor under the WofE approach are summarized in Tab. 3. According to the weights assigned to the aspect classes, the highest positive correlation is observed for the south-western slopes. A steeper slope indicates a greater probability of landslide detachment: below 20%, the contrast is -2.224, indicating very low probability; from 40 to 60%, the contrast is 0.104 and increases in the three subsequent classes. Among the considered lithological classes, arenitic rocks belonging to Auronzo, Fernazza and La Valle formations (Ar) and calcareous and dolomitic rocks belonging to S. Vigilio group (Co) exhibit the greatest importance. In the permeability classes, the contrast value is higher when the permeability is low. The contrast is higher for some classes of Cambisols and is relatively lower for Leptosols, indicating an increase in the landslide occurrence in soil with medium-to-fine texture. Finally, the contrast value is much lower when forest cover is present (model type 2), indicating a statistically significant reduction in the occurrence of landslide detachments compared to the areas outside the forest (Tab. 4).

Tab. 5 summarizes the overall model statistics of LR. The McFadden pseudo R² is respectively 0.18 and 0.2. The McFadden pseudo R² measures the reduction in maximized log-likelihood. It is conceptually and mathematically close to the ordinary least squares R² ([49]). A McFadden pseudo R² of 1 indicates a perfect fit, whereas a pseudo R² equal to 0 indicates no relationship. When a pseudo R² is greater than 0.2, it shows a relatively good fit ([22], [4]).

The regression coefficients obtained by LR considering only the geo-environmental factors (model type 1) are provided in Tab. 6. The influence of the selected factors on landslide occurrence can be deduced from the sign of the estimated coefficients: a negative coefficient indicates that the corresponding factor is related to a safe area, while a positive coefficient signifies that the probability of landslide detachment is low. The slope is included with a positive coefficient, while the aspect is not statistically significant: landslide probability increases with increasing slope, and no significant additional interaction occurs with the aspect. Among the categorical variables, arenitic rocks belonging to Auronzo, Fernazza and La Valle formations (Ar) and to S. Cassiano formation (ACc) and calcareous and dolomitic rocks belonging to Serla Inf. Dolostone and Werfen formation (DC) have the strongest relationship with landslide occurrence, whereas the DB1 soil type prove to have a slight positive influence. Permeability is not a relevant predictor of landslide detachment in the case study: only rocks with medium secondary permeability are included in the model with slightly significant coefficient.

Forest cover has a high negative coefficient (model type 2) which provides further evidence that forest cover significantly improves slope stability (Tab. 7).

Although similarities exist with respect to predictor factors (e.g., the regression coefficient sign for the slope and the exclusion of aspect as a significant factor), differences arise between models of types 1 and 2. For example, the model type 2 replaces DA1 and DB3 soil classes included in the model type 1 with DB2 and DB4. For the model type 2, medium secondary permeability class has a more significant positive coefficient than model type 1. Among the lithological variable, arenitic rocks belonging to Auronzo, Fernazza and La Valle formations (Ar) and to S. Cassiano formation (ACc) and calcareous and dolomitic rocks belonging to Serla Inf. Dolostone and Werfen formation (DC) confirm the strongest relationship with landslide occurrence, but the comparison between models reveals difference in estimated coefficient values.

However, on the whole, the obtained evidence corroborates the known interrelationship between forest cover and geoenvironmental factors.

The results of the ROC analysis on the validation dataset are reported in Fig. 4. The LR approach provides results similar to those of the WofE approach: the AUC ranges from 0.757 for the WofE model type 1 to 0.768 for the WofE model type 2, while the AUC ranges from 0.734 for the LR model type 1 to 0.779 for the LR model type 2. LR approach shows a slightly better predictive capability than WofE when forest cover is included in the model.

Fig. 4 - ROC curves obtained using the studied model types with respect to the validation dataset.

  Enlarge/Shrink   Download   Full Width  Open in Viewer

Considering that LR approach shows better performance than the WofE one, only LR models were used to estimate the forest protective effect. We subtracted the modeled susceptibility value by the LR model type 1 (established without forest cover as predictor) from the value provided by the model type 2 (including forest cover), and re-scaled this difference between 0 and 1. The result is an empirical map which shows the ranking score of protection function associated to the forest stand (Fig. 5): values close to 0 indicate forest stands with an indirect protective function while, values close to 1 identify forest stands with a high potential protective function. The forest with high protective function are mainly located in the middle part of the slope.

Fig. 5 - Map of the relative protective function of the forest stands in the study area (LR approach).

  Enlarge/Shrink   Download   Full Width  Open in Viewer

  Discussion 

In this study, we proposed an approach to assess the protective function of forests through the comparison of landslide susceptibility assessment with and without the current forest cover, with reference to a case study in the Italian Alps. The WofE and LR approaches were comparatively applied to generate multifactor empirical maps of landslide detachment susceptibility.

Both the WofE and the LR methods pinpointed the same predictors of the spatial distribution of landslides. Among the considered factors, those related to slope, soil type and forest cover prove to be the most influential. The probability of occurrence is highest when slope is between 40% and 100%. The soil type affects the occurrence of landslide detachment when soil particles are large. The inclusion of forest cover significantly improves accuracy in the prediction of landslide occurrence: in this case study, the density of landslide detachment areas outside forests is more than twice than that within forests.

Established statistical models proved to be characterized by fair predictive potential, as shown by AUC higher than 0.7 for the validation dataset. The LR appears to give more reliable results than WofE, as it is more sensitive to intermingling between predictive factors.

Overall, the comparison of LR model type 1 with LR model type 2 yielded a high significance of forest cover. The reason is related to the reinforcing effect on slope stability of the diverse root systems of the trees and the associated forest vegetation, as stressed by various authors ([68], [69], [48]). This effect can be also described as an increase in the angle of internal friction or as a virtual increase in soil density ([34]).

Assessing the protective effect of mountain forests on slope stability using landslide susceptibility maps contributes to support resource managers in evaluating forest management practices and setting up suitable measures in forest planning. However, the method may have some limitations that should be stressed in order to correctly interpret the obtained map of forest protective function. The first limitation relates to the scale of analysis. This work was based on factor data available for the study area at different scales. As a result, the accuracy of boundaries in protective function within the map may be potentially misleading at a certain extent. The second key limitation is that the work evaluates forest cover only, and does not consider forest stand characteristics and management practices, since no data are available for the entire study area. On the other hand, it is worth to notice that research on stand characteristics (e.g., stand structure, crown cover, tree density and composition) that optimize the protective functions of forests for specific hazards ([7]) is only at the beginning of its route: stand over-maturity and over-stocking may determine the death and fall rate of large trees or the windthrow rate over large areas with a consequent prolonged loss in protection ([53], [12]), but root cohesion of old-growth forests is greater than that of post-harvest, second-growth stands ([67]).

  Conclusion 

The main objective of this study was to develop a multifactor empirical model to assess the relationship between forest cover and landslide detachment at stand scale. Multifactor methodologies typically estimate landslide susceptibility using the relationship between past landslide patterns and observed site characteristics.

The difference between models with and without the inclusion of forest cover as predictive factor has been exploited to spatially assess the potential role of forest stands in controlling landslide detachment. No comparable results or proposals can be found in the literature relative to this issue, despite a wealth of general references regarding landslide susceptibility assessment. By the empirical approach proposed here, we are able to extract the effect of forest cover on the probability of landslides.

To simultaneously sustain several different functions, large areas of forest should be managed as multi-purpose forests. Under such perspective, the spatial assessment of protective functions allows foresters to make conscious management decisions by comparing the relative importance of this distinctive function within the territorial context. Ultimately, the modeling approach proposed here can be used to pinpoint forest areas with high degree of protection against landslide occurrence on a stand scale, in which targeted policies and careful planning could be important in maintaining their effectiveness for providing protection under a multifunctional sustainable management framework.

  Acknowledgements 

This work was prepared in the Forestry Inventory and Information Systems Laboratory, University of Tuscia, under a contact between Department for Innovation in Biological, Agro-food and Forest systems (DIBAF) and Department of Agriculture, Forests, Nature and Energy (DAFNE).

CD carried out the design of the study and performed the statistical analysis, MS participated in the geological mapping, GS participated in the DEM pre-processing, PL and CP conceived the study, participated in its design and coordination, and helped to draft the manuscript.

  References

  Supplementary Material

Tab. S1 - Lithological classes.
Tab. S2 - Permeability classes.
Tab. S3 - Soil types.