Interactive comment on “ Combination of empirically-based and physically-based methods to assess shallow slides susceptibility at the basin scale ” by Sérgio C

The paper deals with an interesting topic, zoning of shallow landslide susceptibility at basin scale. The original contribution of the paper is the combined use of two different methodologies (in the paper called empirically-based and physically-based methods) to evaluate and map landslide susceptibility over a catchment. The paper is well structured. The proposed procedure, which is based on a “set of integration rules defined by the cross-tabulation of the susceptibility classes of both maps and analysis of the corresponding contingency tables”, is clearly described in the paper. The application of the procedure in the test site (a study area in Portugal) effectively demonstrates the effectiveness of the proposal. The methods chosen in the application (i.e. the bivariate statistical information value method and the infinite slope method) adequately serve

the conditions that cause instability (predisposing factors) can be identified, registered and used to build predictive models; (iii) the occurrence of landslide can be spatially inferred.Within this conceptual scheme, it is assumed that future landslides are more likely to occur in areas where geologic and geomorphologic conditions are similar to those that originate the slope instability in the past (Guzzetti et al., 1999).This conceptual scheme has been extended to different methods of landslide susceptibility assessment, regardless of their nature (Varnes et al., 1984;Hutchinson, 1995;Aleotti and Chowdhury, 1999;5 Carrara et al., 1999;Fell et al., 2008b).This is nonetheless surprising since the conceptual model is perfectly applied to any empirically-based method used to assess landslide susceptibility, but the same is not true for the physically-based methods.
Indeed, the latter methods are based on physical laws and soil mechanics principles, being the slope understood as a system where shear stress and shear strength are continually in opposition.That is, unlike what happens with statistical methods, deterministic methods are applicable not accounting the landslide inventory, which, however, is still essential to validate the 10 obtained landslide susceptibility results.
The comparison between different methods to assess landslide susceptibility is not a new research topic when performed exclusively between different empirically-based statistical methods (Gorsevski et al., 2003;Süzen and Doyuran, 2004;Brenning, 2005;Davis et al., 2006;Lee et al., 2007;Felicísimo et al., 2013;Bui et al., 2016) or between different physicallybased methods (Zizioli et al., 2013;Formetta et al., 2014;Pradham and Kim, 2015;Teixeira et al., 2015).Regarding the 15 comparison of the predictive capacity between empirically-based and physical-based methods, a few number of works exist (Crosta et al., 2006;Carrara et al., 2008;Frattini et al., 2008;Yilmaz and Keskin, 2009;Cervi et al., 2010;Goetz et al., 2011) and from those only a limited number of studies have combined the results obtained with empirically-based and physicallybased approaches (Chang and Chiang, 2009;Goetz et al., 2011).According to Zizioli et al. (2013) the different methods used to assess shallow slides susceptibility are not mutually exclusive.The latter authors pointed out that the use of different 20 strategies to assess landslide susceptibility and the comparison of their predictive capacity can help to: (i) enhance the quality and reliability of each method; (ii) highlight and identify the most important factors affecting the slope instability system; (iii) neglect less influential aspects to simplify the models; and (iv) select the most appropriate methodology to achieve a specified goal.
In this study we aim to verify two hypotheses: (i) although conceptually and methodologically distinct, the statistic and 25 deterministic methods generate similar results for shallow landslides susceptibility regarding the model's predictive capacity and spatial agreement; and (ii) the combination of the shallow landslides susceptibility maps obtained with empirically-based and physically-based methods, for the same study area, generate a more reliable susceptibility map for shallow slides occurrence.

Study area 30
The study area comprises the two small catchments of Monfalim and Louriceira (13.9 km 2 ), which are located 25 km NNW of Lisbon, Portugal (Fig. 1).The elevation ranges from 442 m at the West to 134 m in the northeast sector of the study area, near the confluence of both Monfalim and Louriceira rivers with the Grande da Pipa River (GPR), which is affluent of the Tagus River.
The lithological units are mainly constituted by sedimentary rocks dated from the Kimmeridgian to the Lower Thitonian (Upper Jurassic).Alluvium deposits of Holocene age and a complex of dikes and volcanic masses are also present, both covering only 1.1 % of the study area.The detailed lithological units map for the study area shown in Fig. 1 was constructed 5 based on existing official geological maps (Zbyszewski and Assunção, 1965;INETI, 2005), but also on interpretation of aerial photographs and validation of lithological units limits through field work.Therefore, it was possible to identify eight lithological units enumerated here according the age criteria, being progressive older: (i) alluvium; (ii) Arranhó formation (limestones and marls); (iii) Sobral formation (sandstones and limestones); (iv) Sobral formation (mudstones and marls); (v) Amaral formation (limestones); (vi) Amaral formation (marls); (vii) Abadia formation (mudstones and marls).The 10 lithological unit (viii) is constituted by dykes and volcanic masses (basalt, teschenite, dolerite and weathered rocks).
The study area suffered since the Miocene a wide curvature angle tectonic rebound (Zbyszewski and Assunção, 1965) and the layers dip typically to SE/SW.This structural setting together with the alternation of soft rocks such as marls, clays and mudstones with more resistant rocks as the limestones allowed the development of cuesta-like landforms resulting from differential erosional processes (Ferreira, 1984;Ferreira et al., 1987;Zêzere, 1991).Therefore, in the study area gentle 15 reverse slopes are found over the lithologic units of Sobral and Arranhó formations, whereas abrupt cutting slopes are present along the Amaral limestones lithological unit that outcrops over the erosive depression mainly excavated in the Abadia marls and mudstones formation (Ferreira, 1984).The slopes within the study area are typically moderate: 78.1 % of the total area is within the slope range from 5º to 20º.The gentle slopes (0º -5º) represent only 12.9 % and the steepest slopes (> 20º) occur only in 9 % of the study area . 20 Landslides in the study area have been triggered by rainfall (Zêzere et al., 1999(Zêzere et al., , 2005(Zêzere et al., , 2015;;Zêzere and Rodrigues, 2002;Oliveira, 2012).The climate is Mediterranean and the Mean Annual Precipitation (MAP) is 730 mm (at São Julião do Tojal gauge located 20 km south from the study area) (Zêzere et al., 2015).Shallow slides have been triggered mainly by intense short duration rainfall episodes, typically not exceeding 1 to 15 days (Zêzere and Trigo, 2011;Zêzere et al., 2015).These rainfall events generate the accretion of pore water pressure and the reduction of the soil shear strength, including the loss of 25 cohesion on fine sediments, which promotes the failure along the superficial soil formations or along the contact between the soil and the impermeable bedrock (Trigo et al., 2005).

Methods and data
The methodological procedures for assessing shallow slides susceptibility based on the application and combination of empirically-based and physically-based approaches are summarized in Fig. 2. Two commonly used methods were chosen: 30 the bivariate statistical Information Value (IV) method (Yin and Yan, 1988) and the Infinite Slope method (IS) (Sharma, 2002) sustained on the calculation of the Factor of Safety (FS).Both methods are in line with the experts panel Nat. Hazards Earth Syst. Sci. Discuss., doi:10.5194/nhess-2016-381, 2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.recommendations to assess landslide susceptibility (Cascini, 2008;Fell et al., 2008aFell et al., , 2008b;;Corominas et al., 2014) and have been applied successfully in similar geological and geomorphological context in the region north of Lisbon (Zêzere, 2002;Pimenta, 2011;Guillard and Zêzere, 2012;Oliveira et al., 2015).For shallow slides susceptibility modelling, the dependent variables (shallow slides modelling and validation groups), the independent dataset of variables used as predisposing factors, and the maps representing geotechnical and hydrological parameters were rasterized using a pixel of 5 5 x 5 m.

Landslide inventory
The landslide inventory is used twofold in this study: (i) to establish the statistical relationships between shallow slides and the data-set of environmental factors assumed as shallow slides predisposing factors in the empirically-based approach; and (ii) to validate the shallow slides susceptibility models obtained with both empirically-based and physically-based models.10 The landslide inventory of the study area (Fig. 1) includes 111 shallow slides (translational slides and rotational slides with high curvature angle of the slip surface) that were classified according to Cruden and Varnes (1996) classification.The depth of the slip surface is typically less than 1.5 m.The shallow slides inventory was extracted from (Oliveira, 2012) and was based on interpretation of aerial photographs (1983, 1989) and orthophotomaps (2003, 2004, 2007), as well as on extensive field work made during the period 2006-2010. 15 The inventory of shallow slides was further subjected to a partition based on a temporal criterion (Fig. 1, Table 1).The landslide training group includes the shallow slides that occurred until the end of 1983 (51 cases, 0.027 km 2 , and 0.19 % of the study area).The landslide validation group includes all landslides occurred after 1983 until the end of 2010 (60 cases, 0.03 km 2 , 0.22 % of the study area.The training group was used to weight classes of shallow slides predisposing factors in the statistical model using the IV method, and was also used to calibrate the shear strength parameters (cohesion and friction 20 angle) of the lithological formations in the IS model.The validation group was used for the independent validation of both empirically-based and physically-based shallow slides susceptibility models.

The Information Value method
The Information Value (IV) (Yin and Yan, 1988) was used to compute the susceptibility score for each class of each variable 25 considered as landslide predisposing factor based on the log normalization of the ratio between the conditional probability to find a shallow slide in a certain class of a predisposing factor and the a priori probability to find a shallow slide in the study area, following the Eq. ( 1).
where: Lxi is the Information Value of class xi belonging to an independent variable (predisposing factor); Si is the number of pixels with shallow slides belonging to the training group and the presence of the variable class xi; Ni is the number of pixels with variable class Xi; S is the total number pixels with shallow slides belonging to the training group; and N is the total number of pixels of the study area.Due to the logarithmic normalization li is not calculated when Si = 0.In those cases li was determined as the lowest value considering the complete data set of predisposing factors.Final IV scores (Lxi) for 5 each terrain unit (j) was obtained using Eq. ( 2).
where: m is the total number of variable classes; and Xij is either 0 if the variable class is not present in the pixel j, or 1 if the 10 variable class is present.

Landslide predisposing factors
In this work, we selected as independent variables seven landslide predisposing factors (Fig. 3, Fig. 1 and Table 4 for the description of classes) that have been used with success in previous studies in the region north of Lisbon (e.g., (Oliveira et al., 2015): lithology, slope angle, slope aspect, slope curvature, topographic position index (TPI), slope over area ratio and 15 land use.
The lithologic map includes 8 classes that were already described (cf.Sect. 2. Study area).The Land use map was obtained from the official map representing the land use observed in 1990.Although does not match to the current land use in the study area, this is the one that best fits the time span of shallow landslides included in the present landslide inventory and the temporal land use frame closer to the age of landslides training group.The remaining variables (slope, aspect, curvature, 20 topographic position index and slope over area ratio) were derived from a Digital Elevation Model based on elevation data interpolated from a topographic contours map (equidistance 10 m).For the curvature map, a DEM generalization based on a 50 m pixel size grid was considered to calculate the profile of the slopes, which prove to best fit the morphology of slopes in the study area (Oliveira et al., 2015).The Topographic Position Index (TPI) was calculated based on the Facet Corridor Designer tool for ArcGIS (Jenness et al., 2011).As the index is heavily dependent on the scale (Piacentini et al., 2015) an 25 interactively neighbourhood radius of 25 meters for index calculation proved to be the most appropriate to the work reference scale.The Slope Over Area Ratio (SOAR) was used to express the importance of the topography in hydrological processes by the relation between the slope and the contribution area (Sørensen et al., 2006), which allow to infer the areas prone to surface saturation (Fonseca, 2005).The calculation of the SOAR was made using the TauDEM 5.2 (Terrain Analysis Using Digital Elevation Models) tool and considering the algorithm D8 (O'Callagham and Mark, 1984) to 30 minimize the dispersion of accumulation flow.

The Infinite Slope method (IS)
The most popular formulations of the Infinite Slope method consider a subsurface flow/water table level parallel to the topographic surface, introducing that way, the saturated soil thickness factor.In this context, the development of a steadystate hydrological model in static conditions can be related to the ratio between the thickness of saturated soil and the 5 thickness of the potentially unstable soil as provided in the formulation of SHALSTAB model (Dietrich and Montegomery, 1998).The FS for each terrain unit (pixel) was thus calculated based on the Infinite Slope method, incorporating a soil thickness model and an hydrologic model for the study area, following Eq.( 3) (Sharma, 2002):

) 10
Where: c′ is the effective cohesion (kN/m²); h is the potentially unstable soil depth; β is the slope of the terrain unit; m is the equation component of the hydrological model, considered as the ratio between the saturated soil depth and the potentially unstable soil depth; ϕ′ is the internal friction angle (°); γ is the specific soil weight (kN/m³); γ is the saturated soil weight (kN/m³) and γ is the submerged soil weight (kN/m³).The FS values can be interpreted in two ways.In the more restrict 15 sense it is assumed that all terrain units with FS values ≤ 1 are unstable.A more broad interpretation turns possible to compare FS results with results obtained from the empirically-based approach; this is, to consider that each terrain unit within a study area could be tier according FS values being more susceptible the terrain unit as lower the FS value.
The development of the IS model was supported by the following parameters: (i) topographical variables (slope and catchment area), (ii) soil thickness, (iii) hydrologic parameters (hydraulic conductivity, soil transmissivity and daily rainfall 20 threshold), (iv) geotechnical parameters (natural, saturated and submerged specific soil weights; cohesion; and internal friction angle).Most geotechnical parameters were deduced from references with regional validity that were summarized by (Pimenta, 2011).

Soil thickness model
The depth of the soil potentially unstable is a critical parameter that strongly influences the stability of slopes.The soil depth 25 model for the study area was obtained following Eq.( 4), as proposed by (Catani et al., 2010): Where: h is the soil thickness, K c is a constant calibration parameter, C is an index based on the slope profile curvature, η is 30 the relative soil depth dependent on the topographic position; ψ -1 is the critical slope angle associated to landslide Nat. Hazards Earth Syst. Sci. Discuss., doi:10.5194/nhess-2016-381, 2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.
occurrence.The three parameters C, η and ψ-1 were expressed by linear normalization into a dimensionless index with values ranging between 0 and 1.The constant K c was estimated independently for each lithological unit based on trial and error estimation to fit as much as possible the soil thickness values obtained by Eq. 4 to the soil thickness values measured in 110 sampling field points.These sampling field measurements were spatial distributed in order to guarantee a reasonable number of soil thickness measurements in each lithological unit but also along different geomorphological units (interfluve 5 areas, slopes, valley floors), although dependent of the existence of slope cuts where the soil depth was measured.To accept the K c constant calibration for any lithological unit, the differences between the maximum estimated soil thickness and the maximum soil thickness measured in the field should not exceed 1 m.Table 2 summarizes the K c constant calibration values obtained for each lithological unit in the study area.Soil profiles were not found in LU1, LU3 and LU8 during the field work.In the case of LU 3, we adopted a K c value equal to the one estimated for the other lithologic unit belonging to the 10 Sobral formation (LU 4, K c = 3.6).In the case of alluvium (LU1) and complex of dikes and masses (LU9) we adopted a K c = 2.9, which is the arithmetic mean of all K c values obtained for lithological units where it was possible to measure soil thickness during field sampling.Fig. 4 shows the final soil thickness map of the study area.

Hydrological model
The adopted hydrological model is based on the equation provided by SHALSTAB (Dietrich and Montegomery, 1998), 15 supported by (O' Loughlin, 1986) model.According to Sharma (2002), the hydrologic model corresponds to the ratio between the thickness of saturated soil and the thickness of the potentially unstable soil, according to Eq. ( 5).
Where: h/z is the ratio between the thickness of the saturated soil above the impermeable layer and the thickness of the potentially unstable soil; Q is the effective precipitation (m/day); T is the transmissivity of the soil (m²/day); a is the upstream contribution area (m2); b is the cell length (m); and β is the slope gradient (°).The increase of the hydrologic ratio (Q/T) indicates that soil saturation will be faster and more extensive.The topographic ratio (a/(b * Sinβ)) describes the topography effect on runoff (Dietrich and Montegomery, 1998;Montgomery et al., 1998).The transmissivity of the soil was 25 estimated using Eq. ( 6) (Lencastre and Franco, 2006): Where: T is the soil transmissivity (m²/day); k is the saturated hydraulic conductivity (m/day); and z is the soil thickness (m).30 As the hydraulic conductivity based on field measurements was not available for the study area, this parameter was estimated for the different soil types existing in the study area based on the work developed by (Rawls et al., 1982) The effective precipitation was estimated based on the Eq. ( 7) proposed by Trigo et al. (2005) which defines the rainfall threshold for triggering translational and rotational landslides in the region north of Lisbon, which includes the study area. 10 Where: Cr is the rainfall threshold that is associated to landslides occurrence (mm), and D is the number of consecutive rainfall days.
Most landslide events occur in the study area during the Winter season, so we assume that effect of evapotranspiration can 15 be neglected and the effective precipitation can be assumed equal to total precipitation, namely for short rainfall periods.
Using Eq. ( 7) we obtained a critical daily rainfall for failure of 114.4 mm (0.1144 m).The rainfall concentrated in a single day is a feasible scenario for triggering of shallow landslide events, as it happened in the Lisbon Region in 1967 and 1983 (Zêzere et al., 2005(Zêzere et al., , 2015)).
The hydraulic conductivity was estimated based on the critical precipitation for failure and the soil texture.In the study area 20 k ranges from 5.05 m/day in the luvisols with dominantly sandy texture to 0.0144 m/day in vertisols with dominantly clayey texture.The computed transmissivity ranges between 0 and 13.45 m 2 /day (Fig. 5A).The final hydrological model is shown in Fig. 5B.

Geotechnical parameters of superficial soils
All geotechnical parameters mentioned in this section, related to soil weight (Υm, Υsat Υsub) cohesion (c') and friction 25 angle (ϕ'), were based on literature and were defined for the superficial soils above the bedrock within each lithological unit.
The strength parameters of the lithological units obtained in laboratory with direct shear tests (Pimenta, 2011) proved to be too high to explain the observed slope instability.Therefore, the optimal combinations of cohesion and effective internal 30 friction angle values for each lithological unit were defined iteratively through back analysis.Different combinations of cohesion and effective internal friction angles were tested in the Infinite Slope method and validated with the landslide training group (landslide area) using as reference the maximum and minimum friction angles suggested by (Geotechdata,  (Chung and Fabbri, 2003) is expressed, by the ratio between the percentage of landslide area predicted in the class and the percentage of the class area in the study area.For LU2 and LU5 it was not 5 possible to comply with the criterion (i), but the corresponding critical pair cohesion / internal friction angle were selected respecting criterion (ii).In addition, strength parameters of LU1 and LU8 were not possible to estimate with this method, due to the absence of landslides in these lithological units.In these cases, the cohesion and effective internal friction angle were derived directly from (Pimenta, 2011), which gathered information from technical reports, geotechnical laboratory tests and standard values reported in the literature (Baptista, 2004;Cernica, 1995;Fernandes, 1994;Jeremias, 2000;Vallejo et al., 10 2002).Table 3 summarizes the geotechnical parameters of the lithological units used for implementing the physically-based model.

Validation, comparison and combination of shallow slides susceptibility models
The validation of susceptibility maps produced by empirically-based and physically-based models was made independently using the landslide validation group.ROC (Receiver Operating Characteristic) curves were computed and the corresponding 15 of susceptibility shows a clear contrast between north/northeast sectors of the study area in which the susceptibility is predominantly classified as low to very low, whereas in the centre/south part of the study area the susceptibility to shallow slides is typically higher.This contrast is mainly justified by the lithological differentiation.In fact, in the northern part of the study area the LU7 (Abadia formation: marls and clays) and LU5 (Amaral formation: limestones) are found, which apparently have a low predisposition to shallow slide occurrence (Table 4).By opposition, lithological units more prone to 5 slope instability (LU2 -Arranhó formation: limestones and marls; and LU3 -Sobral formation: sandstones and limestones) outcrop in the centre and south part of the study area.In addition, the slope angle tends to be higher in the latter part of the study area, thus contributing for the higher landslide susceptibility.
The ROC curve of the landslide susceptibility model is shown in Fig. 7.The model predictive capacity is reasonable/good, as expressed by the AUC ROC of 0.75.10

Physically based landslide susceptibility assessment
The shallow slides susceptibility map computed with the IS method is shown in Fig. 8A.The susceptibility class with FS ≤1 (Very high susceptibility) covers 17.9 % of the total study area and validates 53.4 % of the shallow slides belonging to the landslide validation group, which explains the higher effective ratio (2.98) of this susceptibility class (Table 5).By comparison with the IV susceptibility map, it is evident the increment of area classified with very high/high susceptibility in 15 the north sector of the study area where LU7 outcrops, whereas the spatial expression of the two highest landslide susceptibility classes decreases in the southwest/south sector were the LU2 outcrops.The ROC curve of the model based on the landslide validation group is shown in Fig. 7.The ROC curve is detached to the upper left corner of ROC space, which confirms the best predictive capacity of the IS susceptibility map when compared with the IV susceptibility map.The AUC of 0.81 also supports the better predictive capacity of the IS model.20 As it was already mentioned, shallow landslides in the study area have been triggered by rainfall, typically during intense short duration (1 -15 days) rainfall events (Zêzere et al., 2005(Zêzere et al., , 2015;;Zêzere and Trigo, 2011).Additionally, it is known by extensive field work in the study area (Oliveira, 2012) a total absence of instability signs during the summer, reflecting the dryness that characterizes this season.Therefore, it can be assumed a typical situation of superficial absence of water in the soil during summer, i.e., m = 0. Assuming this situation, an additional physically-based shallow slides susceptibility map 25 was prepared considering no water in the soil (m = 0).Figure 8B shows the results of modelling.Given the assumed boundary conditions, it was expectable that model do not generate FS ≤ 1.However, Fig. 8B shows a small fraction of the study area classified with Very high susceptibility (FS ≤ 1, 2.25 % of study area) in a condition of absence of water into the soil, which is interpreted as an error of the IS model.It is worth mentioning that most of the model errors occur over the LU2 (Arranhó formation) indicating that corresponding resistance parameters (cohesion, internal friction angle) may be 30 underestimated.

Comparison of landslide susceptibility models
The comparison of the susceptibility maps produced with IV and IS methods demonstrates that spatially, the susceptibility ranking differs substantially depending on the method used.Indeed, the Kappa coefficient is only 0.23, which means that spatial correlation is moderate, although the reasonable/good predictive capacity of both models attested by the AUC ROC (Fig. 7). 5 The two highest susceptibility classes in the IV landslide susceptibility map spread over 34.1 % of the total study area and the corresponding percentage of predicted shallow slides approaches 69.4 %.The performance of the predictive model is less interesting for the intermediate susceptibility classes (moderate and low), in particular for the low susceptibility class that includes a relevant portion (15.7 %) of shallow slides belonging to the landslide validation group.The IS landslide susceptibility model reveals a better predictive capacity, which is attested by the presence of 83.1 % of the landslide 10 validation group within the two highest susceptibility classes.
The effective ratios calculated for landslide susceptibility classes of both models are summarized in Table 5.The effective ratios for the IS model are higher for the Very high and High susceptibility classes and lower for the Low and Very low susceptibility classes, which indicate a better predictive capacity when compared with the IV model.
The spatial comparison of the two susceptibility maps is shown in Fig. 9.The value zero means the spatial 15 agreement between landslide susceptibility classes, whereas the other values mean disagreement.Negative values indicate that landslide susceptibility obtained with IV is lower when compared with the map obtained with IS, with the difference increasing from -1 to -4.For example, a grid cell with a score -4 means this terrain unit was classified as very high susceptibility in the IS susceptibility map and as very low susceptibility in the IV susceptibility map.Positive values indicate the opposite relationship between map classes.The perfect spatial agreement between susceptibility classes in both maps 20 occurs in 39.9 % of the study area (Table 6).However, adding the minimum mismatch classification (-1 and +1 in Fig. 9) the previous feature rises to 73 % of the total study area.The major discrepancy between both susceptibility maps (-4, -3, 3 and 4 in Fig. 9) occurs along 10.5 % of the study area, namely where the Abadia formation (LU7) and the Arranhó formation (LU2) outcrop.In the north part of the study area where the LU7 is present, the landslide susceptibility obtained with the IV method is lower when compared with the IS method, whereas the opposite occurs in the central and southern part of the 25 study area where the LU2 is present.These results can be interpreted according to particular specifications associated with the physically-based and empirically-based methods.The resistance parameters estimated for the superficial soil over LU7 (c '= 2 kPa, φ' = 19 °) are higher than those estimated for LU2 (c '= 0.5 kPa φ '= 17 °).However, the landslide susceptibility computed using the IS tends to be higher over LU7, which is related to the soil water content and eventually to the presence of thicker soils, 30 particularly along the lower part of slopes where topographic conditions are more prone to soil saturation.On the other hand, the empirically-based approach generated IV scores of 0.494 and -0.857, respectively for LU2 and LU7.The positive IV score for LU2 clearly indicates a higher chance for shallow slides occurrence.We admit that shallow slides inventory may be Nat. Hazards Earth Syst. Sci. Discuss., doi:10.5194/nhess-2016-381, 2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.incomplete in the area corresponding to LU7, which could justify the negative IV score.Indeed, the LU7 clays and marls are associated with gentle slopes and are characterized by intense agricultural use; thus, the footprint of small shallow slides is easily erased on the landscape, as the "original" slope profile is recovered for agricultural activities.On the contrary, the LU2 is constituted by sequences of marl and limestone layers, which induce larger topographic irregularities and less productive soils on steep to moderate slopes.These geological and geomorphological conditions favoured along time a land 5 use mainly associated to forest and annual crop cultures.In this context, the landslide footprint over slopes tends to last more in time, which justifies a more complete shallow slides inventory, and consequently, the higher IV score.

Combination of landslide susceptibility models
The results of the cross-tabulation between landslide susceptibility classes of both susceptibility maps (empirically-based and physically-based) are summarized in a contingency table (Table 6).The distribution of shallow slides belonging to the 10 validation group on the same contingency table is summarized in Table 7. Table 6 also shows the considered combinations within the contingency table to classify the final landslide susceptibility map resulting from the integration of empiricallybased and physically-based predictive models, where the colours (red, orange, yellow, light green, green and grey) represent the final susceptibility classes (Very high, High, Moderate, Low, Very low, and uncertain, respectively).The corresponding final shallow slides susceptibility map is shown in Fig. 10 and information about final landslide susceptibility classes is 15 detailed in Table 8.
The Very high susceptibility class covers 16.4 % of the study area and includes 55.6 % of the shallow slides validation group.Similarly, the High susceptibility class covers 14.3 % of the study area and includes 18.6 % of the shallow slides.In opposition, , the Very low and Low susceptibility classes cover 33.4 % and 10.6 % of the study area, respectively, and include only a small fraction of the landslide validation group (1.4 % each class).20 Terrain units classified as Very high or High susceptibility by one method and simultaneously as Very low or Low susceptibility by the other method were considered as uncertain regarding susceptibility to shallow slides occurrence in the final map.The 'grey' class, although classified as Uncertain, is potentially High or Very high landslide susceptible and covers 16.3 % of the study area and includes 16.0 % of the shallow slides belonging to the validation group.However, the distribution of landslide validation group in the Uncertain susceptibility class is different in the upper right corner and in the 25 lower left corner of the contingency table (see Tables 6 and 7).Terrain units classified as Very high or High susceptibility by the IS susceptibility map and as Very low or Low susceptibility by the VI method (upper right corner in Tables 6 and 7) includes 14.7 % of shallow slides belonging to the validation group, whereas terrain units with inverse classification (lower left corner in Tables 6 and 7) only contain 1.2 % of the shallow slides validation group.These values, once more, reflect the higher quality of the physically-based susceptibility model in comparison with the empirically-based model.30 The predictive quality of susceptibility classes constituting the final landslide susceptibility map is demonstrated by the estimated effective ratios (Table 8).The effective ratio of the Very high susceptibility class (3.39) is higher than those obtained for the equivalent susceptibility class with the empirically-based and physically-based methods (cf.Table 5).In Nat. Hazards Earth Syst. Sci. Discuss., doi:10.5194/nhess-2016-381, 2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.addition, effective ratios corresponding to the Very low and Low susceptibility classes (0.04 and 0.12, respectively) are lower than those obtained with empirically-based and physically-based methods (cf.Table 5), which indicates a better predictive performance.Moreover, the effective ratio is higher for the Uncertain class than for the Moderate class (Table 8), which is consistent with the potential for high or very high susceptibility considered for the Uncertain class.

Conclusion 5
Empirically-based and physically-based methods used to assess landslide susceptibility at the basin scale are conceptually distinct as the former are based on weighting environment predisposing factors, whereas the latter are supported by the computation of shearing and resistance forces along potential slip surfaces.The existence of a landslide inventory is crucial to weight predictive variables within empirically-based methods, which is not the case of physicallybased methods that can be computed independently on the landslide inventory.Both types of methods have advantages and 10 drawbacks.The major constrains associated to empirically-based approaches have been summarized in previous works (Corominas et al., 2014;Fell et al., 2008a) and result from: (i) the difficulty of establishing causal (cause-effect) relationships between variables; (ii) problems arising from self-correlation between variables; (iii) the typically not normal statistical distribution of predictor variables; (iv) the limitations related to the quality of data, in particular the completion of the landslide inventory; and (v) the difficulty in transferring the results from the study area to other areas, even with similar 15 characteristics.In the case of physically-based methods, the major constrains were listed as follow (Corominas et al., 2014;Fell et al., 2008a): (i) the high level of generalization and/or simplification regarding the spatial distribution of geotechnical or hydrological parameters; (ii) the feasibility of model application is limited to areas with relatively homogeneous ground conditions (e.g., geology and geomorphology); (iii) the uncertainties about the depth of the soil and of the slip surface; and (iv) the difficulties in predicting groundwater pore pressures and their relationship with rainfall.20 In this work we intent to test two hypotheses: (i) although conceptually distinct, empirically-based and physicallybased methods generate similar results concerning susceptibility to shallow slide occurrence; and (ii) a reliable landslide susceptibility map can be obtained for a single study area by combining two landslide susceptibility models (empiricallybased vs physically-based).
To achieve the proposed objectives the Information Value method and the Infinite Slope method were chosen to 25 build two landslide susceptibility maps.A shallow slides inventory was separated into two independent landslide groups adopting a temporal criterion.The training group was used twofold to define the statistical relationships between landslides and the dataset of variables assumed as landslide predisposing factors by the IV method, and to calibrate the resistance parameters (cohesion and internal friction angle) within the IS method.The landslide validation group was used to independently validate both susceptibility maps.30 When analysed independently, both methods generated good predictive results, although the physically-based model revealed to be more effective to spatially predict shallow landslides, which is attested by the AUC ROC and the Nat. Hazards Earth Syst. Sci. Discuss., doi:10.5194/nhess-2016-381, 2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.
, whichNat.Hazards Earth Syst.Sci.Discuss., doi:10.5194/nhess-2016-381,2016   Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.summarized the typical hydraulic conductivities for different soil types starting from the respective textural properties.The national digital soil map at 1: 25,000 scale was used to extract the clay, silt + sand, and coarse sand fractions for the different soils types present in the study area.The soil taxonomy of the US Department of Agriculture was used to distinguish between soil types, through the Soil Texture Triangle Bulk Density.Rocky outcrops and urban areas were assigned with a -1 value, thus corresponding to 0 (absence of water) in the hydrological model.The castanozems soils were also assigned with 5 a -1 value because the typical pedological stage of castanozem soils within the study area is a stony soil phase.At final, 55 types of soils were identified, in addition to social areas and rocky outcrops.
Nat. Hazards Earth Syst.Sci.Discuss., doi:10.5194/nhess-2016-381,2016   Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.2013).Critical pairs of cohesion and internal friction angle were selected for each lithological unit by combining two criteria: (i) the susceptibility class with FS ≤ 1 must include at least 50 % of landslide area of the landslide training group located on the lithological unit; and (ii) the susceptibility class with FS ≤ 1 must have the highest effective ratio.The effective ratio of a susceptibility class

Area
Under the Curve (AUC) was calculated.Additionally, the landslide susceptibility maps were classified and the effective ratio of each class was estimated.Both empirically-based and physically-based susceptibility maps were classified considering the same fraction of study area in each equivalent landslide susceptibility class.First, the IS map was classified into 5 classes based on the Factor of Safety values (≤ 1, 1 to 1.25; 1.25 to 1.5, 1.5 to 2, and > 2), which have correspondence, respectively, on the following descriptive classification of susceptibility (Very high; High; Moderate, Low; and Very low).20 In a second step, the IV map was classified into 5 classes (Very high; High; Moderate, Low; Very low), ensuring that equivalent susceptibility classes cover the same fraction of the study area in both maps.The evaluation of the spatial agreement between landslide susceptibility maps based on empirically-based and physically-based approaches was made using the Rank Difference Tool included in ArcSDM (Sawatzky et al., 2008).Lastly, empirically-based and physically-based susceptibility maps were combined into a final shallow slides susceptibility 25 map based on the intersection of the susceptibility classes in a contingency table using the Map Comparison Kit tool (e.g., Visser and Nijs, 2006), on a cell by cell comparison and Kappa statistics.scores calculated for each class of predisposing factors based on the landslide training group are 30 summarized in Table 4, and the corresponding shallow slides susceptibility map is shown in Fig. 6.The spatial distribution Nat.Hazards Earth Syst.Sci.Discuss., doi:10.5194/nhess-2016-381,2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.

Figure 9 : 5 Figure 10 :
Figure 9: Spatial agreement between IV and IS shallow slides susceptibility maps.0 means full agreement; 4 and -4 means maximum disagreement.

Table 6 .
Contingence table extracted from the overlay of IV and IS shallow slide susceptibility maps in % of the study area.Colours represent the susceptibility classes of the final map: Red -Very high; Orange -High; Yellow -Moderate; Light green -Low; Green -Very low; Grey -Uncertain, but with potential for high/very high susceptibility.5 Nat.Hazards Earth Syst.Sci.Discuss., doi:10.5194/nhess-2016-381,2016 Manuscript under review for journal Nat.Hazards Earth Syst.Sci.Published: 9 December 2016 c Author(s) 2016.CC-BY 3.0 License.