Estimating return period of landslide triggering by Monte Carlo simulation

Accepted Manuscript Estimating return period of landslide triggering by Monte Carlo simulation D.J. Peres, A. Cancelliere PII: DOI: Reference:

S0022-1694(16)30141-X http://dx.doi.org/10.1016/j.jhydrol.2016.03.036 HYDROL 21140

To appear in:

Journal of Hydrology

Please cite this article as: Peres, D.J., Cancelliere, A., Estimating return period of landslide triggering by Monte Carlo simulation, Journal of Hydrology (2016), doi: http://dx.doi.org/10.1016/j.jhydrol.2016.03.036

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Estimating return period of landslide triggering by Monte Carlo simulation D.J. Peresa,∗, A. Cancellierea a

Department of Civil Engineering and Architecture, University of Catania, Via Santa Sofia, 64 – 95123 Catania (Italy)

Abstract Assessment of landslide hazard is a crucial step for landslide mitigation planning. Estimation of the return period of slope instability represents a quantitative method to map landslide triggering hazard on a catchment. The most common approach to estimate return periods consists in coupling a triggering threshold equation, derived from an hydrological and slope stability processbased model, with a rainfall intensity-duration-frequency (IDF) curve. Such a traditional approach generally neglects the effect of rainfall intensity variability within events, as well as the variability of initial conditions, which depend on antecedent rainfall. We propose a Monte Carlo approach for estimating the return period of shallow landslide triggering which enables to account for both variabilities. Synthetic hourly rainfall-landslide data generated by Monte Carlo simulations are analysed to compute return periods as the mean interarrival time of a factor of safety less than one. Applications are first conducted to map landslide triggering hazard in the Loco catchment, located in highly landslide-prone area of the Peloritani Mountains, Sicily, ∗

Corresponding author Email address: [email protected] (D.J. Peres)

Italy. Then a set of additional simulations are performed in order to evaluate the traditional IDF-based method by comparison with the Monte Carlo one. Results show that return period is affected significantly by variability of both rainfall intensity within events and of initial conditions, and that the traditional IDF-based approach may lead to an overestimation of the return period of landslide triggering, or, in other words, a non-conservative assessment of landslide hazard. Keywords: Debris flow, Hazard mapping, TRIGRS, Neyman Scott, rainfall intensity-duration-frequency, antecedent rainfall

1

1. Introduction

2

Landslide susceptibility and hazard mapping can be effectively used as

3

an aid for urban and landslide mitigation planning, which often requires a

4

multidisciplinary approach (Carrara, 1983; Carrara et al., 1991; Van Westen

5

et al., 1997; Guzzetti et al., 1999; Lee, 2004; Ayalew and Yamagishi, 2005;

6

H¨ urlimann et al., 2006; Gorsevski et al., 2006; Conoscenti et al., 2008; Dewitte

7

et al., 2010; Oh and Lee, 2011; Conforti et al., 2012; Pradhan, 2013; Regmi

8

et al., 2013; Bregoli et al., 2015). Several authors map landslide hazard in

9

terms of return period of landslide triggering (Borga et al., 2002; D’Odorico

10

et al., 2005; Rosso et al., 2006; Salciarini et al., 2008; Tarolli et al., 2011;

11

Lanni et al., 2012; Schilir`o et al., 2015). To this end, models considering

12

at least both rainfall intensity and duration as control factors in landslide

13

triggering are suitable (Wu and Sidle, 1995; Baum et al., 2002; Iverson, 2000;

14

D’Odorico et al., 2005; Rosso et al., 2006; Simoni et al., 2008; Baum et al.,

15

2008, 2010; Sorbino et al., 2010; Greco et al., 2013; Capparelli and Versace,

2

16

2014). In these works the estimation of return period is generally carried

17

out by coupling hydrological rainfall infiltration and geomechanical slope-

18

stability physically-based models with rainfall intensity duration-frequency

19

(IDF) relationships, these providing the link between rainfall events and their

20

long-term frequency of occurrence (see Stedinger et al., 1993; Burlando and

21

Rosso, 1996). Simplistic assumptions commonly made within this approach

22

include representation of rainfall events as uniform (i.e. of constant intensity)

23

hyetographs, and the use of prefixed initial conditions. On the other hand,

24

as shown by D’Odorico et al. (2005) and by Peres and Cancelliere (2014), the

25

shape of the hyetograph or, in other words, the variability of instantaneous

26

rainfall intensity within events, may have a significant effect on the triggering

27

of landslides. Use of only rainfall duration and average intensity to character-

28

ize the rainfall events’ potential to trigger landslides, though very common in

29

literature (Guzzetti et al., 2007), may not be sufficient. On the other hand,

30

initial conditions are not properly taken into account, since in most of the

31

cited studies on hazard mapping they are fixed with no regard to their prob-

32

ability of occurrence, which generally may affect return period of landslide

33

triggering. For instance, in Rosso et al. (2006) two return-period maps are

34

presented making the assumption of an initial water table height of 0 and of

35

0.15 m, without taking into account the different probability associated to

36

these two different initial conditions, whilst in the work by D’Odorico et al.

37

(2005) the initial conditions are derived by the model of Montgomery and

38

Dietrich (1994), but no probability is assigned to the steady-state rainfall

39

required by such a model.

40

In this paper we use the Monte Carlo simulation procedure presented

3

41

in Peres and Cancelliere (2014) to show quantitatively how the two above-

42

mentioned hydrological factors may affect the estimation of the return period

43

of shallow landslide triggering. The Monte Carlo simulation approach essen-

44

tially consists in combining a stochastic rainfall model able to generate fine-

45

resolution (hourly) rainfall data for a physically-based hillslope hydrological

46

model. These latter model is suited to compute initial conditions based on an-

47

tecedent hillslope response, the built of transient pressure head due to rainfall

48

events, and finally geomechanical slope stability. Specifically, the following

49

models are used are in our framework: a seasonal Neyman-Scott rectangular

50

pulses rainfall stochastic model (Neyman and Scott, 1958; Rodriguez-Iturbe

51

et al., 1987a,b; Cowpertwait, 1991; Cowpertwait et al., 1996) and the TRI-

52

GRS v.2 unsaturated model (Baum et al., 2008, 2010), combined with a water

53

table recession model based on the linear reservoir hypothesis to compute ini-

54

tial conditions implicitly linked to antecedent rainfall. Finally, return period

55

of shallow landsliding is estimated based on the analysis of the generated

56

synthetic pressure head data. Results obtained by the Monte Carlo method

57

are compared to those obtained by the ”traditional” IDF-based approach,

58

in order to demonstrate and quantify how the simplified assumptions of this

59

latter approach can affect return period estimation.

60

An application is carried out to map shallow landslide triggering hazard

61

in the Loco catchment, located in the Peloritani Mountains, Sicily, Italy.

62

Then, sensitivity analyses are conducted in order to verify the generality of

63

considerations about the reliability of the traditional IDF-based method.

4

64

2. Methods

65

2.1. The Monte Carlo method and return period estimation

66

Generally speaking, the Monte Carlo method consists in the use of a

67

stochastic model for generating the input to a mathematical model which

68

represents the behavior of the physical system under study, and then to

69

analyse statistically the output (Salas, 1993).

70

The Monte Carlo simulation technique for synthetic rainfall-landslide

71

data generation is illustrated briefly in Fig. 1 – for a more detailed descrip-

72

tion see Peres and Cancelliere (2014), where the method has been used with

73

the aim of deriving landslide-triggering thresholds suitable for early warning.

74

First, NRE individual rainfall events are generated from a 1000-year long

75

hourly synthetic rainfall time series, obtained as a Neyman-Scott rectangu-

76

lar pulses process (see Appendix A). For isolating the events from the whole

77

series the following criterion is adopted: when two wet spells are separated

78

by a dry time interval less than ∆tmin =24 hours, these are considered to

79

belong to the same rainfall event; otherwise they are considered as separate.

80

24 hours is the minimum time interval necessary to avoid overlapping of the

81

response produced by subsequent rainfall events for the analyzed hydraulic

82

properties (Peres and Cancelliere, 2014) (see Tabs. 1 and 4) – a similar ap-

83

proach is adopted by Balistrocchi et al. (2009) and Balistrocchi and Bacchi

84

(2011). Then, the hillslope response to the sequence of generated events

85

i = 1, 2, . . ., NRE is computed, by the following steps:

86

1. The TRIGRS unsaturated model is used during each event to compute

87

the transient pressure head ψ1 (Baum et al., 2008, 2010) (see Appendix

88

B, Eqn. B.1). Since pressure head may continue rising after the end 5

89

of rainfall, the computation of transient pressure head is prolonged for

90

∆ta = ∆tmin −1 hours after the ending instant tend,i of any given rainfall

91

event.

92

2. The instant tf,i = max(tend,i , tmax,i ), where tmax being the time instant

93

at which maximum transient pressure head occurs, is computed. It

94

follows that the final response to rainfall event i, in terms of water

95

table height, is ψ(dLZ , tf,i )/ cos2 δ (slope parallel flow is assumed), dLZ

96

being soil depth.

97

computed, with ti+1 being the instant at which rainfall event i + 1

98

begins.

(in)

Moreover, the time interval ∆ti+1 = ti+1 − tf,i is (in)

99

3. The water table height at the beginning of rainfall event i + 1 is com-

100

puted by a sub-horizontal drainage model formula (see Appendix B,

101

Eqn. B.2) which uses ψ(dLZ , tf,i )/ cos2 δ and ∆ti+1 .

102

The result of the whole procedure above described, is a series of maximum

103

pressure head responses, each corresponding to a given rainfall event. An

104

event is considered to trigger a landslide if it yields a factor of safety for

105

slope stability less than one, i.e. if the critical pressure head is exceeded: ( ) ϕ′ c′ − 1 − tan γS dLZ sin δ cos δ tan δ ψCR = (1) γw tan ϕ′

106

where c′ is soil cohesion for effective stress, ϕ′ is the soil friction angle for

107

effective stress, γw is the unit weight of groundwater, γs is the soil unit weight,

108

δ is the slope angle and dLZ is soil depth. Equation 1 is derived from infinite

109

slope stability analysis (see Appendix B) (Taylor, 1948).

110

of landslide triggering is then computed by its definition, i.e. as the mean

111

inter-arrival time of a pressure head ψ > ψCR . This return period computed 6

Return period

112

as here described is denoted as TR in the ensuing text.

113

The procedure described for a single infinite hillslope includes assump-

114

tions that allow us to conveniently apply it on a spatially distributed fashion,

115

based on information derived from digital terrain models (DTMs). In partic-

116

ular, infinite slope stability analysis, and the neglect of the possible additional

117

infiltration on a point due to runoff caused by rainfall exceeding infiltration

118

capacity on upslope areas (see, .e.g. Rosso et al., 2006; Baum et al., 2010),

119

implies that the possible failure of cells within a catchment is assumed in-

120

dependent from the geomechanical and hydrological processes that occur on

121

the other catchment cells. Hence, once grids of spatially-distributed charac-

122

teristics are derived, the return period grid can be computed by applying cell

123

to cell the Monte Carlo results of the corresponding infinite slopes. To this

124

aim, in order to reduce the computational time required for the mapping,

125

interpolation can be used to derive return period using the results obtained

126

for single infinite hillslopes varying the spatially-distributed characteristics

127

– i.e., slope δ, upslope specific contributing area A/B, soil depth dLZ and

128

eventually the other controlling parameters (see Fig. 1).

129

Multi-dimensional analysis of slope stability may change totally the re-

130

turn period analysis for a real catchment in comparison to that of an isolated

131

infinite slope, since failure mechanism conceptually changes and requires the

132

iterative search of the least-stable group of cells (cf. Bellugi et al., 2015).

133

It may be worthwhile to mention that we are focusing on the objective

134

to derive the return period of landslide triggering for given catchment con-

135

ditions. These implies that no modeling of soil parameters temporal change,

136

and in particular soil depth – which generally changes after a landslide event

7

137

– has to be made, as return periods are computed to have the meaning of an

138

hazard metric for given catchment conditions. The same stands in regards

139

of land-use temporal change.

140

2.2. Return period estimation by the traditional approach: use of IDF curves

141

In order to assess the accuracy of the IDF-based approach, we compare its

142

results with the Monte Carlo simulation results. In the following we describe

143

how the IDF-based approach is here implemented. Choices given below, in-

144

cluding the simple-scaling formulation of IDFs, the use of GEV distribution,

145

and the way the deterministic thresholds are derived, are common to sev-

146

eral studies (Borga et al., 2002; D’Odorico et al., 2005; Rosso et al., 2006;

147

Salciarini et al., 2008; Tarolli et al., 2011).

148

149

IDF curves may be expressed in the following form, in the case of a rainfall process for which simple-scaling is valid (Burlando and Rosso, 1996): IT (D) = pT P¯ (1)Dn−1

(2)

150

where IT (D) is the rainfall mean intensity associated to rainfall duration D

151

and return period T , P¯ (1) is the mean annual maxima rainfall on a hourly

152

duration, n is the scaling exponent, which is 0 < n < 1, and pT is the di-

153

mensionless rainfall quantile corresponding to a non-exceedance probability

154

equal to

155

tribution is fitted to a rainfall depth dimensionless data set Pi (Dj )/P¯i (Dj ),

156

for i = 1, 2, . . . , M , where M is the number of available annual maxima data

157

of rainfall depth Pi for each considered duration j = 1, 2, . . . , D;P¯i (τj ) is

158

the mean of annual maxima precipitation for duration τj .The dimensionless

T −1 . T

To determine such dimensionless quantiles, a probability dis-

8

159

160

quantile pT is derived from the Generalized Extreme Value (GEV) distribution, whose cumulative distribution function (cdf) is as follows: { ( )− ξ1 } 0 x − µ0 FX (x) = exp − 1 + ξ0 σ0

(3)

161

where µ0 , σ0 and ξ0 ̸= 0 are parameters to be estimated from the sample,

162

for instance by the maximum likelihood method (cf. Kotz and Nadarajah,

163

2000). For consistency with the Monte Carlo procedure, we estimate the

164

IDF curve parameters from the annual maxima series, of various durations,

165

extracted from the NSRP synthetic hourly series (see Sect. 3.2).

166

Then the IDF model is combined with a physically-based threshold that

167

is derived by assuming rainfall intensity I constant over duration D. If we de-

168

note the physically-based threshold as ICR = f (DCR ) then the return period

169

of each point is given by the return period of the IDF curve that intersects

170

it in that point, so the physically-based threshold has not a constant return

171

period. Hence, the minimum return period is assumed as the return period

172

of shallow landslide triggering. Graphically (cf., e.g. Fig. 12 of Rosso et al.,

173

2006), this corresponds to the IDF curve which is tangent to the physically-

174

based threshold. In the ensuing text we denote as TR0 the return period

175

estimated by this procedure. In the case of the TRIGRS unsaturated model

176

(which is briefly described in Appendix B), these physically-based thresh-

177

olds cannot be expressed in closed form, so we derive them numerically. In

178

doing this we have assumed an initial water table height at the soil-bedrock

179

interface, in order to properly compare the results with the Monte Carlo sim-

180

ulations in the no-memory ψ0 = 0 case. Some examples of these curves are

181

shown in Fig. Appendix B (discussed on later Sect. 4.1). 9

182

Another issue that we may focus on concerns post-event analysis (back-

183

analysis) of occurred landslides. In particular, consistently to the tradi-

184

tional IDF-based approach, the return period of a landslide event, to which

185

a triggering-rainfall critical intensity ICR and duration DCR can be associ-

186

ated, may be computed as the return period the IDF curve passing for the

187

(DCR , ICR ) point, which follows from T = 1/(1 − FX (pT )), and Eqns. 2 and

188

3: T (DCR , ICR ) =

 (  1 − exp − 1 + 

1 ICR −µ0 (k−1) ¯ (1)D P CR σ0

)− ξ1  . 0 ξ0

(4)



189

The critical duration DCR is assumed as the time interval that starts at the

190

beginning of the rainfall event and finishes at the instant at which ψCR (Eqn.

191

1) is reached, and the critical intensity is ICR =

192

accumulated over duration DCR . In the case that ψCR is reached after the end

193

of the rainfall event, then DCR = Dtot and ICR =

194

ti

195

over duration Dtot . This commonly-applied procedure of back-analysis (cf.

196

Schilir`o et al., 2015) neglects the return period of initial conditions, and

197

thus may generally lead to a non correct assessment of the return period of

198

the landslide event. In order to put into evidence the potential magnitude of

199

these errors, the minimum value of back-analysis return period given by Eqn.

200

4, denoted in the following as TR2 , is also computed and compared to the

201

Monte Carlo TR , which represents the most correct way to estimate return

202

period of landsliding for the given conditions. In other words, the return

203

period associated by IDF analysis with the lowest-magnitude Monte Carlo

204

simulated triggering rainfall event (virtual landslide), TR2 , is compared to

(in)

WCR , DCR

Wtot , Dtot

being WCR rainfall

being Dtot = ti,end −

the total duration of the rainfall event and Wtot rainfall accumulated

10

205

the return period of landslide triggering derived by Monte Carlo simulation,

206

TR .

207

3. Application

208

3.1. The Loco catchment, Peloritani Mountains, Italy

209

An application of the proposed approach has been carried out with ref-

210

erence to the Loco catchment in the Peloritani Mountains, located nearby

211

the coastline of northeastern Sicily, Italy, which has an extent of 0.14 km2 .

212

Figure 2 shows location of the catchment, the 2×2 resolution digital terrain

213

model (DTM) of year 2007 used in our application, and of the slides triggered

214

on 1 October 2009, which evolved as debris flows, killing 37 persons.

215

Maps of data for model application are shown in Fig. 3, and in particular:

217

slope δ, depth to basal boundary dLZ , flow accumulation map (number of √ upslope draining cells Nd ) and ε = dLZ cos δ/ A. Depth dLZ of the permeable

218

strata, mainly composed of loamy sands with a high proportion of gravel,

219

has been measured to be ranging from about 2.5 m to 1.5 m, on points of

220

slope varying from about 35◦ to 45◦ . Following other researchers (DeRose,

221

1996; Salciarini et al., 2006, 2008; Baum et al., 2010), we assume a negative

222

exponential relationship between dLZ and δ, which reproduces the available

223

measurements. In particular, the relationship dLZ = 32 exp(−0.07δ) has been

224

adopted, which yields dLZ (35◦ ) = 2.76 m and dLZ (45◦ ) = 1.37 m. Bedrock

225

outcropping has been assumed in areas with slope δ > 50◦ , while an upper

226

bound of dLZ = 5 m has been considered. Flow accumulation map expressing

227

the number of draining cells Nd , and hence upslope contributing area A has

228

been computed using the single direction (D8) method (O’Callaghan and

216

11

229

Mark, 1984), which presents the advantages of robustness, while being not-

230

dispersive and hence providing calculation of upslope specific catchment area

231

232

that is consistent with its definition (Tarboton, 1997). √ The map of the ratio ε = dLZ cos δ/ A was calculated in order to check

233

the applicability of the vertical infiltration models, which requires that ε ≪ 1

234

(Iverson, 2000) (See also Appendix B). From the map it can be inferred

235

that for most of the potentially unstable areas ε is in the range 0.1 − 0.5,

236

which is acceptably less than one (cf.

237

than 0.1 in the drainage network. A small proportion of the catchments

238

has a ratio in the range 0.5 − 1, and few cells have a ratio greater than 1

239

(maximum value is 2.5). Overall the condition ε ≪ 1 is acceptably satisfied.

240

Tab. 1 shows the other parameters required by the adopted models, namely

241

the saturated hydraulic conductivity KS , the saturated hydraulic diffusivity

242

D0 , the residual water content θr , the saturated hydraulic content θs , the α

243

parameter of the Gardner SWRC (cf. Eqn. B.1), soil and water unit weights

244

γs and γw = 9800 N/m3 , soil friction angle for effective stress ϕ′ and cohesion

245

c′ .

246

3.2. Rainfall modelling

Baum et al., 2010). Ratio is less

247

The NSRP model has been calibrated using the rainfall series measured

248

at Fiumedinisi rain gauge installed nearby the catchment, available for the

249

period 02/21/2002 - 02/09/2011 (almost 9 years) at a 10 minutes resolution.

250

The rainfall series presents six homogeneous rainfall seasons (see Peres and

251

Cancelliere, 2014): (i) September and October, (ii) November, (iii) Decem-

252

ber, (iv) January - March, (v) April and (vi) May - August. Thus, separate

253

sets of parameters of the NSRP model have been estimated for each of the 12

254

four rainy seasons, while the last two seasons have been considered to be of

255

negligible rainfall for our scopes. In particular, model parameters λ, ν, β, η

256

and ξ have been calibrated by the method of moments (cf., e.g., Rodriguez-

257

Iturbe et al., 1987a,b; Cowpertwait et al., 1996; Calenda and Napolitano,

258

1999), while b has been determined by trials in the range 0.6 ≤ b ≤ 0.9,

259

following Cowpertwait et al. (1996). Parameters obtained from calibration

260

are shown in Tab. 2. Although calibration has been carried out taking into

261

account seasonality, by calibrating the model separately for the various ho-

262

mogeneous seasons within the year, it is noteworthy to point out that the

263

generated series are globally stationary, being the final aim to assess land-

264

slide hazard under the present climatic condition. It is therefore beyond our

265

scope to estimate the effect of climatic change on landslide hazard (see also

266

end of sect. 2.1).

267

For the derivation of the IDF curves, we extracted annual maxima of

268

various durations from the generated synthetic rainfall series. In particular,

269

simple scaling has been assumed for rain durations 3 ≤ D ≤ 96 hours, and

270

the durations of 3, 6, 12, 24, 48 and 96 hours were considered for calibration

271

of IDF curves, fitting a GEV distribution (see Eqn. 3) for the renormal-

272

ized variate. Maximum likelihood estimators yield the following values of

273

the distribution parameters: µ0 = 0.1787, σ0 = 0.2845 and ξ0 = 0.7742.

274

Furthermore, P¯ (1) = 66.74 mm and n = 0.178. Figure 4 shows in a double-

275

logarithmic I-D plane the derived IDF curves for several return periods.

13

276

4. Results and discussion

277

4.1. Landslide hazard mapping

278

Figure Appendix B shows selected outcomes of the several Monte Carlo

279

simulations performed to map return period in the Loco catchment (Fig.

280

6). In particular the points indicate the rainfall characteristics duration-

281

average intensity of triggering events. In the same figure the TRIGRS-v2

282

deterministic threshold, derived for constant hyetographs and a fixed initial

283

water table height hi = 0 is shown with the related tangent IDF curve. In

284

analysing our simulations, we assume that 1000-years of simulation enable a

285

reliable assessment of return period not greater than 100 years, which roughly

286

corresponds to at least 10 triggering events in 1000 years. Thus return periods

287

greater than 100 years are simply indicated as T > 100. Analogous results

288

are also shown in table 3, which includes the ratios between traditionally-

289

estimated return period TR0 and the one estimated via the Monte Carlo

290

approach TR .

291

Plots on the first column of Fig. Appendix B show, for increasing slope

292

values δ, the results for the no-memory case ψ0 = 0, when a null water

293

table height is assumed at the beginning of rainfall events (hi = 0), see Eqn.

294

B.2. In this case, the effect of antecedent rainfall is neglected and therefore

295

the spread of the points only reflects the effect of rainfall intensity variability

296

within events. To correctly interpret the results it may be worthwhile to note

297

that soil depths dLZ = 32 exp{−0.07δ} and critical pressure head – or critical

298

wetness ratios ζCR (Eqn. B.3) – are a function of the slope (see Tab. 3). As it

299

can be seen, to neglect rainfall intensity variability within events can yield a

300

return period that is 2 - 3 times overestimated (TR /TR 0 = 0.3÷0.5). Plots in 14

301

the second column of Fig. Appendix B show the effect of antecedent rainfall

302

(pressure head memory), which increases with increasing specific upslope

303

contributing area A/B. As A/B increases the return period TR of shallow

304

landslide triggering drastically decreases. In this case the ratio

305

lower than 1/3 (Tab. 3).

TR TR0

is even

306

Thus it can be concluded that the traditional IDF-based approach may

307

yield significantly overestimated return periods, since it does not account for

308

the variability of rainfall intensity and of initial conditions.

309

The plots indicate also that the IDF-based return periods associated with

310

the Monte Carlo simulation triggering points, computed by Eqn. 4 result

311

drastically lower than the ones computed by the traditional IDF approach.

312

In Fig. 6, the map of the return periods of shallow landslide triggering for

313

the Loco catchment is shown. From the TR map of Fig. 6, it can be inferred

314

that areas corresponding to TR < 25 years fall within the areas affected by

315

landslides on 1 October 2009. Nevertheless, there are regions within the 1

316

October 2009 slided areas of return period greater than 100 years.

317

Regarding the predictive skill of the physically-based landslide analysis

318

model the map of Figure 6 has to be considered as a preliminary assessment

319

of landslide hazard, though return mapping using the same detail of infor-

320

mation is common to many other studies (e.g., Borga et al., 2002; D’Odorico

321

et al., 2005; Rosso et al., 2006; Salciarini et al., 2008; Baum et al., 2010;

322

Tarolli et al., 2011; Schilir`o et al., 2015). Lack of real slide data for a signif-

323

icant number of independent landslide events, and the uncertainty affecting

324

the data, related also to difficulties in separating mass wasting area due to

325

landsliding from those due to bed erosion occurring in debris flow propaga-

15

326

tion (Stancanelli and Foti, 2015), hamper a throughout assessment of model’s

327

predictive skill; nevertheless, preliminary analyses (Peres, 2013) have shown

328

that predictive performance of TRIGRS v.2 for the Peloritani area is similar

329

to that assessed for other catchments (cf. Baum et al., 2010), and that the

330

use the water table recession model (Eqn. B.2) to estimate initial conditions

331

leads to a significant improvement in the replication of the 1 October 2009

332

event. One major source of inaccuracy is the use of infinite slope stabil-

333

ity analysis, that tends to underestimate the factor of safety; in fact diverse

334

studies have demonstrated predictive power improvement obtained thanks to

335

advanced three dimensional slope stability analysis (Lehmann and Or, 2012;

336

Milledge et al., 2014; Bellugi et al., 2015; Anagnostopoulos et al., 2015). An-

337

other factor that may have importance in the mapping accuracy, is the role of

338

vegetation on slope stability, whose main effect is the generation of spatially-

339

variable root cohesion (Wu and Sidle, 1995; Wilkinson et al., 2002; Hwang

340

et al., 2015); here this effect has not been taken into account for lack of data.

341

Nevertheless, while these issues may be relevant for the accuracy of return

342

period mapping, they do not affect the results related to the evaluation of

343

the traditional IDF-based procedure (our main focus), since their effect con-

344

sists in a modified mapping of the critical pressure head ψCR , which would

345

be the same for both the Monte Carlo and the traditional method. Another

346

point is that a complete assessment of debris-flow hazard should include both

347

landslide triggering and propagation aspects (Berenguer et al., 2015). How-

348

ever, the assessment of return period is not affected by the run-out process,

349

generally assumed deterministic, and thus results of the comparison between

350

TR and TR0 are still valid in the case of scenarios that consider also debris

16

351

flow propagation.

352

4.2. Further simulations

353

Following the same line of thought of Peres and Cancelliere (2014) a

354

sensitivity analysis has been conducted with respect to the following variables

355

(Table 4): the hydraulic conductivity KS , the leakage ratio cd , the soil depth

356

dLZ and the critical wetness ratio ζCR (Eqn. B.3). In particular, plots similar

357

to those of Fig. Appendix B can be derived for each set of values of such

358

variables. In Figs. 7–11 plots of return periods TR0 and TR are shown as a

359

function of ζCR , for various values of KS (cf. Figs. 7, 8 and 9) and cd (cf.

360

Figs. 7, 10 and 11).

361

Results generally confirm the fact that the traditional procedure leads

362

to significantly non-conservative estimations of landslide triggering hazard.

363

As expected, the presence of memory of antecedent rainfall A/B = 10 m

364

increases the overestimation of return period, and such an effect is stronger

365

as the saturated conductivity KS increases. With regards to the no-memory

366

ψ0 = 0 case, it can be observed that the effect of the variation of KS is

367

less dramatic. Also the leakage ratio cd affects significantly the frequency of

368

landslide triggering, that decreases with the increase of cd . Return period

369

increases with soil depth dLZ as well.

370

It can be observed how return period increases almost exponentially with

371

ζCR , till a constant return period from a certain ζCR is reached, for both

372

the considered methodologies (TR0 and TR ). This ζCR corresponds to the

373

height of the capillary fringe, which for the model assumptions is equal to

374

1/α = 1/3.5 = 29 cm (see Appendix B).

17

375

5. Conclusions

376

In this study we have used a Monte Carlo simulation framework for the

377

quantitative assessment of landslide hazard, in terms of return period of

378

shallow landslide triggering. The approach enables to take into account the

379

effect of intensity variability within rainfall events, as well as pore pressure

380

memory which determines the conditions at the beginning of rainfall events

381

in dependence of antecedent rainfall. These aspects are not neglected in

382

commonly applied (traditional) methods based on coupling rainfall Intensity-

383

Duration-Frequency (IDF) curves with physically-based landslide-triggering

384

thresholds.

385

Applications based on soil parameters and topographic data measured on

386

a landslide-prone catchment in the Peloritani Mountains (Italy) have shown

387

that to neglect rainfall intensity time-variability during events and the effect

388

of antecedent rainfall on initial conditions may lead to significant overesti-

389

mation of return period of landslide triggering.

390

The main conclusions are: a) the effect of rainfall intensity variability dur-

391

ing events may significantly affect the return period of landslide triggering.

392

Commonly applied approaches, which implicitly assume constant-intensity

393

rainfall during events generally lead to non-conservative hazard assessment;

394

based on our simulations we found that the return period may be overes-

395

timated by a factory greater than two if the variability of rainfall intensity

396

within events is neglected; b) in several studies return period of landslide

397

triggering is computed by arbitrarily fixing initial conditions, i.e. the initial

398

water table depth. Such a way to proceed has the drawback of neglecting

399

the probability of occurrence of that initial condition. As shown by our sim18

400

ulations, depending on the strength of memory, which for a given climate

401

depends on the τM parameter, a pre-fixed water table initial condition may

402

be more or less probable and thus significantly modify hazard assessment.

403

In particular, common assumption of an initial water table depth at the

404

base of the impervious layer may lead to a dramatic overestimation of the

405

landslide-triggering return period; c) due to rainfall intensity variability and

406

to antecedent rainfall memory, the return period of a rainfall event associated

407

to a landslide (used in back-analysis as an estimation of occurred landslide

408

events), may be arbitrarily underestimated if computed as the return period

409

of the IDF curve passing through the critical rainfall duration and intensity

410

point, since the probability of antecedent rainfall is not taken into account.

411

One of the main disadvantages of the Monte Carlo approach is that it is

412

way more computationally demanding than the IDF-based one. Hence such a

413

traditional approach may still be useful for preliminary assessments of land-

414

slide hazard taking into account extreme rainfall climate features. Simplified

415

procedures to include the effects of both rainfall intensity variability and

416

initial conditions are thus further directions of our research.

417

Acknowledgments This work has been partially funded by the projects

418

PON no.

01 01503 ”Integrated Systems for Hydrogeological Risk Mon-

419

itoring, Early Warning and Mitigation Along the Main Lifelines”, CUP

420

B31H11000370005 and PON02 000153 2939551 ”Development of innovative

421

technologies for energy saving and environmental sustainability of shipyards

422

and harbour areas” (SEAPORT), of the Italian Education, University and

423

Research Ministry (MIUR).

19

424

Appendix A. Neyman-Scott rectangular pulses rainfall model

425

Stochastic rainfal models such as the Neyman-Scott rectangural pulses

426

(NSRP) model, belong to the so-called class of cluster models (Neyman

427

and Scott, 1958; Kavvas and Delleur, 1975; Waymire and Gupta, 1981a,b,c;

428

Rodriguez-Iturbe et al., 1987a,b; Salas, 1993; Cowpertwait et al., 1996). The

429

NSRP process is obtained by the following steps:

430

• First, clusters are originated by a Poisson process of parameter λt

431

• For each cluster origin, rectangular pulses (rain cells) are generated.

432

The number of pulses C associated with each storm is extracted from

433

another separate Poisson distribution. In particular, to have realiza-

434

tions of C not less than one, it is assumed that C ′ = C − 1, with

435

436

437

438

c′ = 0, 1, 2, . . . (which implies c = 1, 2, 3, . . .), is Poisson distributed with mean ν − 1 • Each cell has origin at time τi,j with j = 1, 2, . . ., ci measured from ti , according to an exponential random variable of parameter β

439

• A rectangular pulse of duration di,j and intensity xi,j corresponds to

440

each rain cell. Pulses have duration exponentially distributed with

441

parameter η, while intensities X are extracted from a Weibull distribu-

442

tion (cf. Cowpertwait et al., 1996), which has cumulative distribution

443

function F (x; ξ, b) = 1 − exp(−ξxb )

444

445

• Finally, the total intensity at any point in time is given by the sum of the intensities of all active cells at that point.

20

446

Appendix B. Landslide-triggering model

447

Following various researchers (cf. Iverson, 2000; D’Odorico et al., 2005),

448

we consider pressure head ψ response to rainfall events as the superimposition

449

of an initial ψ0 and a transient part ψ1 , the former being related to antecedent

450

rainfall, while the latter being strictly related to single rainfall events. In

451

452

case of shallow soils for which the ratio between soil thickness and the square √ root of upslope draining area is small, ε = dLZ cos δ/ A ≪ 1, ψ1 can be

453

computed by 1-D vertical infiltration equations, while ψ0 is mainly related to

454

sub-horizontal drainage occurring during no-rainfall periods (Iverson, 2000).

455

To compute the transient part ψ1 , we use the TRIGRS v. 2 model (Baum

456

et al., 2008, 2010) which is based on the Richards’ (1931) vertical-infiltration

457

equation for a sloping surface, and the Gardner’s (1958) exponential soil–

458

water retention curve K(ψ) = Ks exp{α(ψ − ψcf )}: α1 (θs − θr ) ∂K ∂2K ∂K = − α , 1 KS ∂t ∂Z 2 ∂Z

(B.1)

459

where Ks is the saturated hydraulic conductivity, α is the SWRC parameter,

460

ψcf = −1/α is the pressure head at the top of the capillary fringe, θr is

461

the residual water content, θs is the water content at saturation, and α1 =

462

α cos2 δ.

463

Initial conditions ψ0 to equation B.1 are updated in terms of initial water

464

table depth hi at the beginning of each rainfall event i, using the following

465

equation: hi =

(

ψ(dLZ , tf,i−1 ) Ks sin δ exp − A ∆ti 2 cos δ (θ − θr ) B s

) =

i ψ(dLZ , tf,i−1 ) − ∆t τM e , cos2 δ

(B.2)

466

where A is the contributing area draining across the contour length B of the

467

lower boundary of the hillslope, δ is the inclination of the hillslope, Ks is the 21

468

saturated hydraulic conductivity, and θs − θr is soil porosity. The pressure

469

head at the end of preceding rainfall event ψ(dLZ , tf,i−1 ) and the inter-arrival

470

time ∆ti are defined in Sect. 2.1. The time constant τM regulates the pressure

471

head memory from one event to another. The ratio A/B, which can be

472

computed based on a digital terrain model (DTM), is the well-known specific

473

upslope contributing area (cf. Montgomery and Dietrich, 1994; Holmgren,

474

1994; Borga et al., 2004). Since we use the non-dispersive single direction

475

(D8) method (O’Callaghan and Mark, 1984), it is A/B = B Nd , where Nd

476

is the number of cells draining into the local one.

477

The factor of safety FS for slope stability is computed using a infinite slope ψ dLZ cos2 δ

478

model. Equivalently, slope failure occurs when the wetness ratio

479

exceeds the critical one, defined as the ratio between critical pressure head

480

(i.e. corresponding to a factor of safety for slope stability FS = 1) and

481

pressure head at saturation

482

as follows: ζCR

γs = γw

[(

ψCR , dLZ cos2 δ

which for the infinite slope formula is

) ] c′ tan δ −1 +1 . γs dLZ sin δ cos δ tan ϕ′

(B.3)

483

The ζCR varies from 0 to 1, respectively, for an unconditionally unstable

484

and a unconditionally stable hillslope (Montgomery and Dietrich, 1994), and

485

hence it is a metric of the natural degree of stability of the hillslope.

486

In order to understand the controlling factors of shallow landslide trig-

487

gering, it is useful to separate the analysis of the response to rainfall in terms

488

of the transient part only. This may be done by performing the simulations

489

assuming a water depth at the soil-bedrock interface as an initial condition

490

for all rainfall events hi = 0 (denoted throughout the paper as the ψ0 = 0,

491

no memory, case), while for A/B > 0, there is the presence of pressure head 22

492

493

494

memory, i.e. antecedent rainfall determines generally ψ0 > 0. Further details of the model components briefly described in these appendixes are given in Peres and Cancelliere (2014).

23

495

List of Symbols

496

A

497

A/B upslope specific contributing area

498

B

contour (stream tube) length

499

D0

soil saturated hydraulic diffusivity

500

IT

rainfall mean intensity of return period T

501

K(ψ) hydraulic conductivity

502

KS

soil saturated hydraulic conductivity

503

NRE

Number of generated synthetic rainfall events

504

TR

return period of landslide-triggering estimated via the Monte Carlo approach

505

506

upslope drainage area

TR0

return period of landslide triggering as computed by the traditional (IDF) approach

507

508

TR2

back-analysis IDF-based return period of a given landslide

509

WCR , ICR , DCR critical (corresponding to slope failure) rainfall event cumulative depth, intensity and duration

510

511

Z

vertical depth measured from ground surface

512

P¯ (1) mean annual maxima rainfall depth on a hourly duration

513

c′

soil cohesion for effective stress 24

514

cd

leakage ratio

515

dLZ

soil depth

516

hi

water table height at the beginning of rainfall event i

517

n

rainfall scaling exponent, relative to the rainfall IDF curve

518

pT

dimensionless rainfall depth quantile of return period T

519

ti

520

tend,i Ending instant of i-th synthetic rainfall event

521

tmax,i Time instant at which the maximum transient pressure head occurs

(in)

Time instant at which i-th rainfall event begins

for i-th rainfall event

522

523

α

soil water retention curve (SWRC) parameter

524

δ

terrain slope respect to an horizontal reference

525

γS

unit weight of soil

526

γw

unit weight of water

527

λ, ν, β, η, ξ,b parameters of the Neyman-Scott rectangular pulses (NSRP) stochastic rainfall model

528

529

µ0 , σ0 , ξ0 Generalized Extreme Value (GEV) distribution parameters

530

ϕ′

soil friction angle for effective stress

531

ψ

pressure head

25

532

ψ0

initial (at the beginning of rainfall events) part of pressure head

533

ψ1

transient part of pressure head

534

ψCR

critical pressure head, corresponding to slope incipient unstability

535

τM

water table recession model time constant

536

θr

soil residual water content

537

θs

soil saturated water content

538

ζCR

critical soil wetness

26

539

Anagnostopoulos, G. G., Fatichi, S., Burlando, P., 2015. An advanced

540

process-based distributed model for the investigation of rainfall-induced

541

landslides: The effect of process representation and boundary conditions.

542

Water Resources Research 51 (9), 7501–7523.

543

Ayalew, L., Yamagishi, H., 2005. The application of GIS-based logistic regres-

544

sion for landslide susceptibility mapping in the Kakuda-Yahiko Mountains,

545

Central Japan. Geomorphology 65 (1-2), 15–31.

546

Balistrocchi, M., Bacchi, B., 2011. Modelling the statistical dependence of

547

rainfall event variables through copula functions. Hydrology and Earth

548

System Sciences 15 (6), 1959–1977.

549

Balistrocchi, M., Grossi, G., Bacchi, B., 2009. An analytical probabilistic

550

model of the quality efficiency of a sewer tank. Water Resources Research

551

45 (12), n/a–n/a.

552

Baum, R. L., Godt, J. W., Savage, W. Z., 2010. Estimating the timing and

553

location of shallow rainfall-induced landslides using a model for transient,

554

unsaturated infiltration. Journal of Geophysical Research 115, F03013.

555

Baum, R. L., Savage, W. Z., Godt, J. W., 2002. TRIGRS – A FORTRAN

556

Program for Transient Rainfall Infiltration and Grid-Based Regional Slope-

557

Stability Analysis. U.S. Geological Survey Open-File Report 02-0424, Re-

558

ston, Virginia.

559

Baum, R. L., Savage, W. Z., Godt, J. W., 2008. TRIGRS – A FORTRAN

560

program for transient rainfall infiltration and grid-based regional slope-

27

561

stability analysis, version 2.0. U.S. Geological Survey Open-File Report

562

2008-1159, Reston, Virginia.

563

Bellugi, D., Milledge, D. G., Dietrich, W. E., McKean, J. A., Perron, J. T.,

564

Sudderth, E. B., Kazian, B., 2015. A spectral clustering search algorithm

565

for predicting shallow landslide size and location. Journal of Geophysical

566

Research: Earth Surface 120 (2), 300–324.

567

Berenguer, M., Sempere-Torres, D., H¨ urlimann, M., 2015. Debris-flow fore-

568

casting at regional scale by combining susceptibility mapping and radar

569

rainfall. Natural Hazards and Earth System Sciences 15 (3), 587–602.

570

Borga, M., Dalla Fontana, G., Cazorzi, F., 2002. Analysis of topographic

571

and climatic control on rainfall-triggered shallow landsliding using a quasi-

572

dynamic wetness index. Journal of Hydrology 268 (1-4), 56–71.

573

Borga, M., Tonelli, F., Selleroni, J., 2004. A physically based model of the

574

effects of forest roads on slope stability. Water Resources Research 40 (12),

575

1–9.

576

Bregoli, F., Medina, V., Chevalier, G., H¨ urlimann, M., Bateman, A., 2015.

577

Debris-flow susceptibility assessment at regional scale: Validation on an

578

alpine environment. Landslides 12 (3), 437–454.

579

Burlando, P., Rosso, R., 1996. Scaling and muitiscaling models of depth-

580

duration-frequency curves for storm precipitation. Journal of Hydrology

581

187 (1-2), 45–64.

582

Calenda, G., Napolitano, F., 1999. Parameter estimation of Neyman-Scott 28

583

processes for temporal point rainfall simulation. Journal of Hydrology

584

225 (1-2), 45–66.

585

Capparelli, G., Versace, P., 2014. Analysis of landslide triggering conditions

586

in the sarno area using a physically based model. Hydrology and Earth

587

System Sciences 18 (8), 3225–3237.

588

Carrara, A., 1983. Multivariate models for landslide hazard evaluation. Jour-

589

nal of the International Association for Mathematical Geology 15 (3), 403–

590

426.

591

Carrara, A., Cardinali, M., Detti, R., Guzzetti, F., Pasqui, V., Reichenbach,

592

P., 1991. GIS techniques and statistical models in evaluating landslide

593

hazard. Earth Surface Processes & Landforms 16 (5), 427–445.

594

Conforti, M., Robustelli, G., Muto, F., Critelli, S., 2012. Application and

595

validation of bivariate GIS-based landslide susceptibility assessment for the

596

Vitravo river catchment (Calabria, south Italy). Natural Hazards 61 (1),

597

127–141.

598

Conoscenti, C., Di Maggio, C., Rotigliano, E., 2008. GIS analysis to assess

599

landslide susceptibility in a fluvial basin of NW Sicily (Italy). Geomorphol-

600

ogy 94 (3-4), 325–339.

601

Cowpertwait, P., 1991. Further developments of the Neyman-Scott clustered

602

point process for modeling rainfall. Water Resources Research 27 (7), 1431–

603

1438.

604

Cowpertwait, P. S. P., O’Connell, P. E., Metcalfe, A. V., Mawdsley, J., 1996. 29

605

Stochastic point process modelling of rainfall. I. Single-site fitting and val-

606

idation. Journal of Hydrology 175 (1-4), 17–46.

607

DeRose, R., C., 1996. Relationships between slope morphology, regolith

608

depth, and the incidence of shallow landslides in eastern Taranaki hill

609

country. Z. Geomorphol. 105, 49–60.

610

Dewitte, O., Chung, C.-J., Cornet, Y., Daoudi, M., Demoulin, A., 2010.

611

Combining spatial data in landslide reactivation susceptibility mapping:

612

A likelihood ratio-based approach in W Belgium. Geomorphology 122 (1-

613

2), 153–166.

614

D’Odorico, P., Fagherazzi, S., Rigon, R., 2005. Potential for landsliding: De-

615

pendence on hyetograph characteristics. Journal of Geophysical Research:

616

Earth Surface 110 (F1).

617

Gardner, W. R., 1958. Some steady-state solutions of the unsaturated mois-

618

ture flow equation with evaporation from a water table. Soil Science 85,

619

228–232.

620

Gorsevski, P., Gessler, P., Foltz, R., Elliot, W., 2006. Spatial prediction of

621

landslide hazard using logistic regression and ROC analysis. Transactions

622

in GIS 10 (3), 395–415.

623

Greco, R., Comegna, L., Damiano, E., Guida, A., Olivares, L., Picarelli, L.,

624

2013. Hydrological modelling of a slope covered with shallow pyroclastic

625

deposits from field monitoring data. Hydrology and Earth System Sciences

626

17 (10), 4001–4013.

30

627

Guzzetti, F., Carrara, A., Cardinali, M., Reichenbach, P., 1999. Landslide

628

hazard evaluation: A review of current techniques and their application in

629

a multi-scale study, Central Italy. Geomorphology 31 (1-4), 181–216.

630

Guzzetti, F., Peruccacci, S., Rossi, M., Stark, C. P., 2007. Rainfall thresholds

631

for the initiation of landslides in central and southern Europe. Meteorology

632

and Atmospheric Physics 98, 239–267.

633

Holmgren, P., 1994. Multiple flow direction algorithms for runoff modelling

634

in grid based elevation models: an empirical evaluation. Hydrological Pro-

635

cesses 8 (4), 327–334.

636

H¨ urlimann, M., Copons, R., Altimir, J., 2006. Detailed debris flow hazard as-

637

sessment in Andorra: A multidisciplinary approach. Geomorphology 78 (3-

638

4), 359–372.

639

Hwang, T., Band, L. E., Hales, T. C., Miniat, C. F., Vose, J. M., Bol-

640

stad, P. V., Miles, B., Price, K., 2015. Simulating vegetation controls

641

on hurricane-induced shallow landslides with a distributed ecohydrological

642

model. Journal of Geophysical Research: Biogeosciences 120 (2), 361–378.

643

Iverson, R. M., 2000. Landslide triggering by rain infiltration. Water Re-

644

sources Research 36, 1897–1910.

645

Kavvas, M., L., Delleur, J., W., 1975. The stochastic and cronologic structure

646

of rainfall sequence: Application to Indiana. Tech. Rep. 1325, 129 p., Tech.

647

Rep. 57, Water Resour. Res. Cent., Purdue University, West Lafayette, Ind.

648

Kotz, S., Nadarajah, S., 2000. Extreme value distributions. Vol. 31. World

649

Scientific. 31

650

Lanni, C., Borga, M., Rigon, R., Tarolli, P., 2012. Modelling shallow landslide

651

susceptibility by means of a subsurface flow path connectivity index and

652

estimates of soil depth spatial distribution. Hydrology and Earth System

653

Sciences 16 (11), 3959–3971.

654

Lee, S., 2004. Application of likelihood ratio and logistic regression models

655

to landslide susceptibility mapping using GIS. Environmental Management

656

34 (2), 223–232.

657

Lehmann, P., Or, D., 2012. Hydromechanical triggering of landslides: From

658

progressive local failures to mass release. Water Resources Research 48 (3),

659

n/a–n/a.

660

URL http://dx.doi.org/10.1029/2011WR010947

661

Milledge, D. G., Bellugi, D., McKean, J. A., Densmore, A. L., Dietrich,

662

W. E., 2014. A multidimensional stability model for predicting shallow

663

landslide size and shape across landscapes. Journal of Geophysical Re-

664

search: Earth Surface 119 (11), 2481–2504.

665

Montgomery, D. R., Dietrich, W. E., 1994. A physically based model for the

666

topographic control on shallow landsliding. Water Resources Research 30,

667

1153–1171.

668

Neyman, J., Scott, E. L., 1958. Statistical approach to problems of cosmology.

669

Journal of the Royal Statistical Society. Series B (Methodological) 20 (1),

670

1–43.

671

O’Callaghan, J. F., Mark, D. M., 1984. The extraction of drainage networks

32

672

from digital elevation data. Computer Vision, Graphics, and Image Pro-

673

cessing 28 (3), 323–344.

674

Oh, H.-J., Lee, S., 2011. Landslide susceptibility mapping on Panaon Island,

675

Philippines using a geographic information system. Environmental Earth

676

Sciences 62 (5), 935–951.

677

Peres, D. J., 3 2013. The hydrologic control con shallow landslide triggering:

678

empirical and Monte Carlo physically-based approaches. Ph.D. thesis, Uni-

679

versity of Catania, Department of Civil and Environmental Engineering,

680

Catania.

681

Peres, D. J., Cancelliere, A., 2014. Derivation and evaluation of landslide-

682

triggering thresholds by a Monte Carlo approach. Hydrology and Earth

683

System Sciences 18 (12), 4913–4931.

684

Pradhan, B., 2013. A comparative study on the predictive ability of the

685

decision tree, support vector machine and neuro-fuzzy models in landslide

686

susceptibility mapping using GIS. Computers and Geosciences 51, 350–365.

687

Regmi, N., Giardino, J., McDonald, E., Vitek, J., 2013. A comparison of

688

logistic regression-based models of susceptibility to landslides in western

689

Colorado, USA. Landslides, 1–16.

690

691

Richards, L., 1931. Capillary conduction of liquids through porous mediums. Physics 1, 219–233.

692

Rodriguez-Iturbe, I., Cox, D. R., Isham, V., 1987a. Some Models for Rainfall

693

Based on Stochastic Point Processes. Proceedings of the Royal Society of 33

694

London. Series A, Mathematical and Physical Sciences 410 (1839), 269–

695

288.

696

Rodriguez-Iturbe, I., De Power, B. F., Valdes, J. B., 1987b. Rectangular

697

Pulses Point Process Models for Rainfall: Analysis of Empirical Data.

698

Journal of Geophysical Research 92 (D8), 9645–9656.

699

Rosso, R., Rulli, M. C., Vannucchi, G., 2006. A physically based model for

700

the hydrologic control on shallow landsliding. Water Resources Research

701

42 (6), 1–16.

702

703

Salas, J. D., 1993. Analysis and Modeling of Hydrologic Time Series. McGraw-Hill, Ch. 19, Handbook of Hydrology.

704

Salciarini, D., Godt, J. W., Savage, W. Z., Baum, R. L., Conversini, P., 2008.

705

Modeling landslide recurrence in Seattle, Washington, USA. Engineering

706

Geology 102 (3-4), 227–237.

707

Salciarini, D., Godt, J. W., Savage, W. Z., Conversini, P., Baum, R. L.,

708

Michael, J. A., 2006. Modeling regional initiation of rainfall-induced shal-

709

low landslides in the eastern Umbria Region of central Italy. Landslides 3,

710

181–194.

711

Schilir`o, L., Esposito, C., Scarascia Mugnozza, G., 2015. Evaluation of shal-

712

low landslide-triggering scenarios through a physically based approach: an

713

example of application in the southern messina area (northeastern sicily,

714

italy). Natural Hazards and Earth System Science 15 (9), 2091–2109.

34

715

Simoni, S., F. Zanotti, F., Bertoldi, G., Rigon, R., 2008. Modelling the prob-

716

ability of occurrence of shallow landslides and channelized debris flows

717

using GEOtop-FS. Hydrol. Processes 22, 532–545.

718

Sorbino, G., Sica, C., Cascini, L., 2010. Susceptibility analysis of shallow

719

landslides source areas using physically based models. Natural Hazards

720

53 (2), 313–332.

721

Stancanelli, L. M., Foti, E., 2015. A comparative assessment of two different

722

debris flow propagation approaches blind simulations on a real debris flow

723

event. Natural Hazards and Earth System Science 15 (4), 735–746.

724

Stedinger, J. R., Vogel, R. M., Foufoula-Georgiou, E., 1993. Frequency Anal-

725

ysis of Extreme Events. McGraw-Hill, Ch. 18, Handbook of Hydrology.

726

Tarboton, D. G., 1997. A New Method For The Determination Of Flow

727

Directions And Upslope Areas In Grid Digital Elevation Models. Water

728

Resources Research 33 (2), 309–319.

729

Tarolli, P., Borga, M., Chang, K.-T., Chiang, S.-H., 2011. Modeling shallow

730

landsliding susceptibility by incorporating heavy rainfall statistical prop-

731

erties. Geomorphology 133 (3-4), 199–211.

732

Taylor, D., 1948. Fundamentals of Soil Mechanics. John Wiley, New York.

733

Van Westen, C., Rengers, N., Terlien, M., Soeters, R., 1997. Prediction of

734

the occurrence of slope instability phenomena through GIS-based hazard

735

zonation. International Journal of Earth Sciences 86 (2), 404–414.

35

736

Waymire, E., Gupta, V. K., 1981a. The mathematical structure of rainfall

737

representations: 1. A review of the stochastic rainfall models. Water Re-

738

sources Research 17 (5), 1261–1272.

739

Waymire, E., Gupta, V. K., 1981b. The mathematical structure of rainfall

740

representations: 2. A review of the theory of point processes. Water Re-

741

sources Research 17 (5), 1273–1285.

742

Waymire, E., Gupta, V. K., 1981c. The mathematical structure of rainfall

743

representations: 3. Some applications of the point process theory to rainfall

744

processes. Water Resources Research 17 (5), 1287–1294.

745

Wilkinson, P. L., Anderson, M. G., Lloyd, D. M., 2002. An integrated hydro-

746

logical model for rain-induced landslide prediction. Earth Surface Processes

747

and Landforms 27 (12), 1285–1297.

748

URL http://dx.doi.org/10.1002/esp.409

749

750

Wu, W., Sidle, R. C., 1995. A Distributed Slope Stability Model for Steep Forested Basins. Water Resources Research 31, 2097–2110.

36

Figure captions

751

752

Figure 1

Scheme of the adopted Monte Carlo simulation frame-

753

work: a stochastic model generates a long synthetic

754

rainfall time series which is used as input to a hydro-

755

logical and slope stability model. Specifically, 1000

756

years of hourly rainfall data are obtained by a NSRP

757

model, then the hillslope response to the rainfall events

758

is computed by a combination of a transient infiltra-

759

tion model (TRIGRS v.2) in which the initial condi-

760

tions are computed by a sub-horizontal linear-reservoir

761

drainage model. Subsequently, based on a infinite

762

slope stability analysis, a rainfall event is considered

763

as triggering (FS≤1) or not. Frequency analysis of

764

triggering events enables to finally estimate return pe-

765

riod of landslide triggering.

766

Figure 2

Map showing the location and 2-meters-resolution dig-

767

ital terrain model of the study area, Loco catchment,

768

Sicily (Italy). Slides occurred on 1 October 2009 are

769

also indicated.

37

770

Figure 3

GRID maps for the Loco catchment: (a) slope com-

771

puted from the DTM, (b) soil depth, estimated as

772

dLZ = 32 × exp(−0.07δ), (c) flow accumulation ac-

773

cording to D8 method and (d) map of the ε =

774

to verify reliability of the TRIGRS model in assessing

775

triggering conditions.

776

Figure 4

dLZ √cos δ A

Rainfall intensity-duration-frequency (IDF) curves de-

777

rived from the NSRP model, considered for computa-

778

tion of return periods TR0 and TR2

779

Figure 5

Plots showing the computation of landslide-triggering return

780

period by the Monte Carlo approach TR and the IDF-based

781

traditional method TR0 . In particular for each of the plots the

782

points indicate the synthetic rainfall events which correspond

783

to landslide triggering. The blue curve is the so-called de-

784

terministic threshold, which is derived from the TRIGRS v.2.

785

model by assuming rainfall events of constant-intensity and an

786

initial water table height at the soil-bedrock interface. Return

787

period TR0 is given by the IDF curve which is tangent to such

788

a threshold (green dotted line). Minimum return period TR2

789

computed as the IDF curve intersecting the simulated trigger-

790

ing points (Eqn. 4) is also shown. Plots on the left are relative

791

to the no-memory ψ0 = 0 case, while plots on the right are

792

relative to increasing specific upslope contributing areas A/B.

793

794

Figure 6

Return period of landslide triggering for the Loco catchment estimated by the Monte Carlo approach

38

795

Figure 7

Return period of landslide triggering as a function the

796

critical wetness ζCR : comparison of the traditional

797

(IDF-based) and Monte Carlo methodologies. Case

798

KS = 72 mm/h and cd = 0.10. (a) dLZ = 1 m (b) dLZ

799

= 1.5 m (c) dLZ = 2 m.

800

Figure 8

Same as Fig. 7 but for the following data: KS = 36

801

mm/h and cd = 0.10. (a) dLZ = 1 m (b) dLZ = 1.5 m

802

(c) dLZ = 2 m.

803

Figure 9

Same as Fig. 7 but for the following data: KS = 108

804

mm/h and cd = 0.10. (a) dLZ = 1 m (b) dLZ = 1.5 m

805

(c) dLZ = 2 m.

806

Figure 10

Same as Fig. 7 but for the following data: KS = 72

807

mm/h and cd = 0.05. (a) dLZ = 1 m (b) dLZ = 1.5 m

808

(c) dLZ = 2 m.

809

Figure 11

Same as Fig. 7 but for the following data: KS = 72

810

mm/h and cd = 0.20. (a) dLZ = 1 m (b) dLZ = 1.5 m

811

(c) dLZ = 2 m.

39

Figure 1: Scheme of the adopted Monte Carlo simulation framework: a stochastic model generates a long synthetic rainfall time series which is used as input to a hydrological and slope stability model. Specifically, 1000 years of hourly rainfall data are obtained by a NSRP model, then the hillslope response to the rainfall events is computed by a combination of a transient infiltration model (TRIGRS v.2) in which the initial conditions are computed by a sub-horizontal linear-reservoir drainage model. Subsequently, based on a infinite slope stability analysis, a rainfall event is considered as triggering (FS≤1) or not. Frequency analysis of triggering events enables to finally estimate return period of landslide triggering.

40

Figure 2: Map showing the location and 2-meters-resolution digital terrain model of the study area, Loco catchment, Sicily (Italy). Slides occurred on 1 October 2009 are also indicated.

41

Figure 3: GRID maps for the Loco catchment: (a) slope computed from the DTM, (b) soil depth, estimated as dLZ = 32 × exp(−0.07δ), (c) flow accumulation according to D8 method and (d) map of the ε =

dLZ √cos δ A

to verify reliability of the TRIGRS model in

assessing triggering conditions.

42

3

10

TR0 = 2 years TR0 = 5 years TR0 = 10 years TR0 = 25 years

2

10 Intensity [mm/h]

TR0 = 50 years TR0 = 100 years

1

10

0

10 0 10

1

10 Duration [h]

Figure 4: Rainfall intensity-duration-frequency (IDF) curves derived from the NSRP model, considered for computation of return periods TR0 and TR2

43

δ = 40° ψ0 = 0

δ = 30° ψ = 0 TR > 100 TR2>100

0

10

TR0 > 100 0

Intensity [mm/h]

Intensity [mm/h]

0

TR = 5.35 TR2=4.7

0

10

TR0 = 9.97 0

2

10

10

10 Duration [h]

δ = 40° A/B = 10 m

T > 100

Intensity [mm/h]

Intensity [mm/h]

δ = 35° ψ0 = 0 R

TR2=36.0 0

10

T

R0

0

> 100

T = 3.70 R

TR2=1

0

10

TR0 = 9.97

2

10

0

10 Duration [h]

10

T = 0.40 R

TR2=1

0

TR0 = 1.22 0

10

TR = 2.53 TR2=1

0

10

TR0 = 9.97 0

2

10

10 Duration [h]

TR = 0.19 T =1.0 R2

0

TR0 = 1 0

10

2

10 Duration [h]

δ = 40° A/B = 50 m Intensity [mm/h]

Intensity [mm/h]

δ = 50° ψ0 = 0

10

2

10 Duration [h]

δ = 40° A/B = 20 m Intensity [mm/h]

Intensity [mm/h]

δ = 45° ψ0 = 0

10

2

10 Duration [h]

T = 1.06 R

TR2=1

0

10

TR0 = 9.97 0

2

10

10 Duration [h]

Triggering (FS = 1)

Det. Thr.

2

10 Duration [h]

TR0 IDF curve

Figure 5: Plots showing the computation of landslide-triggering return period by the Monte Carlo approach TR and the IDF-based traditional method TR0 . In particular for each of the plots the points indicate the synthetic rainfall events which correspond to landslide triggering. The blue curve is the so-called deterministic threshold, which is derived from the TRIGRS v.2. model by assuming rainfall events of constant-intensity and an initial water table height at the soil-bedrock interface. Return period TR0 is given by the IDF curve which is tangent to such a threshold (green dotted line). Minimum return period TR2 computed as the IDF curve intersecting the simulated triggering points (Eqn. 4) is also shown. Plots on the left are relative to the no-memory ψ0 = 0 case, while plots on the right are relative to increasing specific upslope contributing areas A/B.

44

Figure 6: Return period of landslide triggering for the Loco catchment estimated by the Monte Carlo approach

45

ψ0 = 0

2

Return period [years]

Return period [years]

10

1

10

0

10

−1

10

0.5 ζCR

2

Return period [years]

Return period [years]

0.5 ζCR

1

0

0.5 ζCR

1

0

0.5 ζCR

1

2

1

0

10

−1

1

10

0

10

−1

0

b)

0.5 ζCR

10

1

2

2

10 Return period [years]

10 Return period [years]

0

10

10

1

10

0

10

−1

c)

0

10

10

1

10

10

1

10

−1

0

a)

10

A/B = 10 m

2

10

1

10

0

10

−1

0

0.5 ζCR

10

1

Traditional (IDF)

Monte Carlo

Figure 7: Return period of landslide triggering as a function the critical wetness ζCR : comparison of the traditional (IDF-based) and Monte Carlo methodologies. Case KS = 72 mm/h and cd = 0.10. (a) dLZ = 1 m (b) dLZ = 1.5 m (c) dLZ = 2 m.

46

ψ0 = 0

2

Return period [years]

Return period [years]

10

1

10

0

10

−1

10

0.5 ζCR

2

Return period [years]

Return period [years]

0.5 ζCR

1

0

0.5 ζCR

1

0

0.5 ζCR

1

2

1

0

10

−1

1

10

0

10

−1

0

b)

0.5 ζCR

10

1

2

2

10 Return period [years]

10 Return period [years]

0

10

10

1

10

0

10

−1

c)

0

10

10

1

10

10

1

10

−1

0

a)

10

A/B = 10 m

2

10

1

10

0

10

−1

0

0.5 ζCR

10

1

Traditional (IDF)

Monte Carlo

Figure 8: Same as Fig. 7 but for the following data: KS = 36 mm/h and cd = 0.10. (a) dLZ = 1 m (b) dLZ = 1.5 m (c) dLZ = 2 m.

47

ψ0 = 0

2

Return period [years]

Return period [years]

10

1

10

0

10

−1

10

0.5 ζCR

2

Return period [years]

Return period [years]

0.5 ζCR

1

0

0.5 ζCR

1

0

0.5 ζCR

1

2

1

0

10

−1

1

10

0

10

−1

0

b)

0.5 ζCR

10

1

2

2

10 Return period [years]

10 Return period [years]

0

10

10

1

10

0

10

−1

c)

0

10

10

1

10

10

1

10

−1

0

a)

10

A/B = 10 m

2

10

1

10

0

10

−1

0

0.5 ζCR

10

1

Traditional (IDF)

Monte Carlo

Figure 9: Same as Fig. 7 but for the following data: KS = 108 mm/h and cd = 0.10. (a) dLZ = 1 m (b) dLZ = 1.5 m (c) dLZ = 2 m.

48

ψ0 = 0

2

Return period [years]

Return period [years]

10

1

10

0

10

−1

10

0.5 ζCR

2

Return period [years]

Return period [years]

0.5 ζCR

1

0

0.5 ζCR

1

0

0.5 ζCR

1

2

1

0

10

−1

1

10

0

10

−1

0

b)

0.5 ζCR

10

1

2

2

10 Return period [years]

10 Return period [years]

0

10

10

1

10

0

10

−1

c)

0

10

10

1

10

10

1

10

−1

0

a)

10

A/B = 10 m

2

10

1

10

0

10

−1

0

0.5 ζCR

10

1

Traditional (IDF)

Monte Carlo

Figure 10: Same as Fig. 7 but for the following data: KS = 72 mm/h and cd = 0.05. (a) dLZ = 1 m (b) dLZ = 1.5 m (c) dLZ = 2 m.

49

ψ0 = 0

2

Return period [years]

Return period [years]

10

1

10

0

10

−1

10

0.5 ζCR

2

Return period [years]

Return period [years]

0.5 ζCR

1

0

0.5 ζCR

1

0

0.5 ζCR

1

2

1

0

10

−1

1

10

0

10

−1

0

b)

0.5 ζCR

10

1

2

2

10 Return period [years]

10 Return period [years]

0

10

10

1

10

0

10

−1

c)

0

10

10

1

10

10

1

10

−1

0

a)

10

A/B = 10 m

2

10

1

10

0

10

−1

0

0.5 ζCR

10

1

Traditional (IDF)

Monte Carlo

Figure 11: Same as Fig. 7 but for the following data: KS = 72 mm/h and cd = 0.20. (a) dLZ = 1 m (b) dLZ = 1.5 m (c) dLZ = 2 m.

50

Table captions

812

813

Table 1

cial strata of Loco catchment

814

815

Material strength and hydraulic properties for surfi-

Table 2

Parameters of the NSRP rainfall model resulting from

816

calibration on Fiumedinisi rainfall data, for the four

817

homogeneous rainy seasons – the Weibull shape pa-

818

rameter has been fixed to b = 0.6 (after Peres and

819

Cancelliere, 2014).

820

Table 3

Comparison between return period estimated by the

821

traditional IDF-based procedure TR0 and from Monte

822

Carlo simulations TR , for simulations perfomed to

823

map return period of landslide triggering on the Loco

824

catchment (see Fig. 6). In bold are indicated cases

825

considered in Fig. 5.

826

827

Table 4

Varied soil properties considered for sensitivity analysis (see also Peres and Cancelliere, 2014).

51

828

Table 5

Comparison between return period estimated by the

829

traditional IDF-based procedure TR0 and from Monte

830

Carlo simulations TR , in the ψ0 = 0 (no-memory

831

case), i.e. isolating the sole effect of rainfall intensity

832

variability. For all these simulations slope δ = 40◦ ,

833

soil depth dLZ = 2 m, soil unit weight γS = 19 000

834

N/m3 , internal friction angle ϕ′ = 39◦ and cohesion

835

c′ = 4 kPa. Several values of the soil hydraulic con-

836

ductivity KS and leakage ratios cd are considered (see

837

Tab. 4).

838

Table 6

Same as Tab. 5 but for the case that pressure-head

839

memory is present (A/B = 10 m), i.e. by considering

840

the effect of antecedent rainfall.

52

Table 1: Material strength and hydraulic properties for surficial strata of Loco catchment

ϕ′ [◦ ] 39

c′

γs

[kPa] [N/m3 ] 4

19000

θs

KS

θr

α

D0

[-]

[m/s]

[-]

[m−1 ]

[m2 /s]

0.35

2 × 10−5

0.045

3.5

5 × 10−5

53

Table 2: Parameters of the NSRP rainfall model resulting from calibration on Fiumedinisi rainfall data, for the four homogeneous rainy seasons – the Weibull shape parameter has been fixed to b = 0.6 (after Peres and Cancelliere, 2014).

Parameter λ [h−1 ] ν β [h−1 ] η [h−1 ] ξ [hb mm−b ]

Jan, Feb, Mar

Sep, Oct

Nov

Dec

0.002295 0.021195

0.001485

0.003185

44.28

1.57

42.41

42.61

0.010161

2.1179

0.0059551

0.0098760

0.72113

0.83999

0.94053

0.67735

1.13441

0.46260

0.69261

1.03521

54

Table 3: Comparison between return period estimated by the traditional IDF-based procedure TR0 and from Monte Carlo simulations TR , for simulations perfomed to map return period of landslide triggering on the Loco catchment (see Fig. 6). In bold are indicated cases considered in Fig. 5.

dLZ

δ

TR

TR TR0

[years] [years]

[–]

KS

cd

A/B

ζCR

[m]

[◦ ] [mm/h]

[–]

[m]

[–]

3.92

30

72 0.1

0

0.930

>100

>100

–

2.76

35

72 0.1

0

0.611

>100

>100

–

1.95

40

72 0.1

0

0.371

9.97

5.35

0.54

1.37

45

72 0.1

0

0.238

1.22

0.40

0.33

0.97

50

72 0.1

0

0.245

1.01

0.19

0.18

3.92

30

72

0.1

10

0.930

>100

>100

–

2.76

35

72

0.1

10

0.611

>100

>100

–

1.95

40

72 0.1

10

0.371

9.97

3.70

0.37

1.37

45

72

0.1

10

0.238

1.22

0.32

0.27

0.97

50

72

0.1

10

0.245

1.01

0.16

0.16

3.92

30

72

0.1

20

0.930

>100

>100

–

2.76

35

72

0.1

20

0.611

>100

>100

–

1.95

40

72 0.1

20

0.371

9.97

2.53

0.25

1.37

45

72

0.1

20

0.238

1.22

0.27

0.22

0.97

50

72

0.1

20

0.245

1.01

0.16

0.15

3.92

30

72

0.1

50

0.930

>100

>100

–

2.76

35

72

0.1

50

0.611

>100

63.48

–

1.95

40

72 0.1

50

0.371

9.97

1.06

0.11

1.37

45

72

0.1

50

0.238

1.22

0.25

0.20

0.97

50

72

0.1

50

0.245

1.01

0.20

0.20

55

TR0

Table 4: Varied soil properties considered for sensitivity analysis (see also Peres and Cancelliere, 2014).

KS

D0

cd

dLZ

τM (A/B = 10 m)

[m2 s−1 ]

[−]

[m]

[days]

1 × 10−5 (36 mm h−1 )

2.5 × 10−5

0.1

1, 1.5, 2

5.5

2 × 10−5 (72 mm h−1 )

5 × 10−5

0.05, 0.1, 0.2

1, 1.5, 2

2.7

3 × 10−5 (108 mm h−1 ) 7.5 × 10−5

0.1

1, 1.5, 2

1.8

[m s−1 ]

56

Table 5: Comparison between return period estimated by the traditional IDF-based procedure TR0 and from Monte Carlo simulations TR , in the ψ0 = 0 (no-memory case), i.e. isolating the sole effect of rainfall intensity variability. For all these simulations slope δ = 40◦ , soil depth dLZ = 2 m, soil unit weight γS = 19 000 N/m3 , internal friction angle ϕ′ = 39◦ and cohesion c′ = 4 kPa. Several values of the soil hydraulic conductivity KS and leakage ratios cd are considered (see Tab. 4).

TR

TR TR0

[years] [years]

[–]

dLZ

KS

cd

ζCR

[m]

[mm/h]

[–]

[–]

1.0

36

0.1

0.789

4.04

1.78

0.44

1.5

36

0.1

0.502

13.82

4.40

0.32

2.0

36

0.1

0.359

13.00

5.33

0.41

1.0

72

0.05

0.789

4.04

1.30

0.32

1.5

72

0.05

0.502

6.19

2.40

0.39

2.0

72

0.05

0.359

5.30

2.53

0.48

1.0

72

0.1

0.789

5.84

2.45

0.42

1.5

72

0.1

0.502

9.85

4.80

0.49

2.0

72

0.1

0.359

11.29

5.25

0.46

1.0

72

0.2

0.789

8.33

4.41

0.53

1.5

72

0.2

0.502

19.02

9.54

0.50

2.0

72

0.2

0.359

23.91

12.24

0.51

1.0

108

0.1

0.789

7.61

3.25

0.43

1.5

108

0.1

0.502

13.82

5.73

0.41

2.0

108

0.1

0.359

13.00

6.37

0.49

57

TR0

Table 6: Same as Tab. 5 but for the case that pressure-head memory is present (A/B = 10 m), i.e. by considering the effect of antecedent rainfall.

TR

TR TR0

[years] [years]

[–]

dLZ

KS

cd

ζCR

[m]

[mm/h]

[–]

[–]

1.0

36

0.1

0.789

4.04

0.43

0.11

1.5

36

0.1

0.502

13.82

1.14

0.08

2.0

36

0.1

0.359

13.00

1.82

0.14

1.0

72

0.05

0.789

4.04

0.69

0.17

1.5

72

0.05

0.502

6.19

1.26

0.20

2.0

72

0.05

0.359

5.30

1.48

0.28

1.0

72

0.1

0.789

5.84

1.39

0.24

1.5

72

0.1

0.502

9.85

3.18

0.32

2.0

72

0.1

0.359

11.29

3.62

0.32

1.0

72

0.2

0.789

8.33

3.04

0.37

1.5

72

0.2

0.502

19.02

7.44

0.39

2.0

72

0.2

0.359

23.91

10.54

0.44

1.0

108

0.1

0.789

7.61

2.48

0.33

1.5

108

0.1

0.502

13.82

4.73

0.34

2.0

108

0.1

0.359

13.00

5.26

0.40

58

TR0

Estimating return period of landslide triggering by Monte Carlo simulation D. J Peresa,* and A. Cancellierea   a

Department of Civil Engineering and Architecture, University of Catania, Catania, Italy

*corresponding author - email: [email protected]

Highlights 

Monte Carlo simulations are performed to compute return period of landsliding



Commonly-used methods for landslide return period estimation are evaluated



Common assumption of constant hyetographs is non-conservative



Initial conditions have a probability to occur which affects return period



Use of rainfall IDF curves may lead to significant overestimations of return period