Sensitivity and identifiability of rheological parameters in debris flow modeling

Over the past decades, several numerical models have been developed to understand, simulate and predict debris flow events. Typically, these models simplify the complex interactions between water and solids using a singlephase approach and different rheological models to represent flow resistance. In this study, we perform a sensitivity analysis on the parameters of a debris flow numerical model (FLO-2D) for a suite of relevant variables (i.e., maximum flood area, maximum flow velocity, maximum height and deposit volume). Our aims are to (i) examine the degree of model overparameterization and (ii) assess the effectiveness of observational constraints to improve parameter identifiability. We use the Distributed Evaluation of Local Sensitivity Analysis (DELSA) method, which is a hybrid local–global technique. Specifically, we analyze two creeks in northern Chile (∼ 29 S, 70W) that were affected by debris flows on 25 March 2015. Our results show that SD (surface detention) and β1 (a parameter related to viscosity) provide the largest sensitivities. Further, our results demonstrate that equifinality is present in FLO-2D and that the final deposited volume and maximum flood area contain considerable information to identify model parameters.


Introduction
In steep mountain environments, intense and localized storms can trigger the sudden movement of sediments, generating flash floods with solid volumetric concentrations up to 40 %-60 % (Takahashi, 1981;O'Brien and Julien, 1988;Calvo and Savi, 2009). These events, also known as debris flows, differ from water floods because -in addition to fluid stress -solid-fluid and solid-solid interactions dominate the flow motion (Takahashi, 1981;Iverson et al., 1997). In recent years, debris flows have been recognized as a major natural hazard (Calvo and Savi, 2009), affecting infrastructure, economic activities and human life. For instance, debris flow events in Switzerland produced 24 fatalities and overall losses of USD 380 million during the period 1972-2007 (Hilker et al., 2009). In Chile, estimated economic losses associated with the five biggest debris flow events recorded over 1980-2017 were at least USD 1.6 billion with nearly 1000 people dead or missing (Servicio Nacional de Geología y Minería, 2017).
Over the last decades, numerical models have emerged as a powerful tool to understand the behavior and magnitude of debris flow events, since they allow for the quantification of key variables used by engineers and decision-makers for risk management (Quan Luna et al., 2011;Frey et al., 2016;Calvo and Savi, 2009) and urban planning (Hürlimann et al., 2006;Lucà et al., 2014;Naef et al., 2006;Arattano et al., 2006). However, the application of debris flow models requires several assumptions and simplifications that make results diverge from reality at various levels (Sosio et al., 2007). For example, uncertainties in terrain elevation models (e.g., satellite product and horizontal resolution), physical parameters (e.g., rheology parameters and solid concentration) and hydrological fluxes (e.g., precipitation and streamflow) used to force debris flow simulations can substantially impact relevant variables, such as flood area, sediment volumes or maximum flow depth.
The use of debris flow models for practical problems typically requires the implementation of single-phase numerical models Naef et al., 2006) that solve one-dimensional or two-dimensional Saint Venant equations, using different rheological approaches to account for frictional stress S f . Many studies have reported good agreement between debris flow model results and post-event measurements, e.g., runout distance, flow velocity, deposit depth and flood area (D'Agostino and Tecca, 2006;Naef et al., 2006;Sosio et al., 2007;Rickenmann et al., 2006;Cesca and D'Agostino, 2008;Lin et al., 2011;Hungr, 1995). Nevertheless, it is recognized that, because complex debris flow dynamics change in time and space (Coussot and Meunier, 1996), the appropriate choice of rheological parameters is critical for a good agreement between debris flow model output and field data (Sosio et al., 2007). In this context, various approaches have been adopted to characterize the sensitivity of debris flow model results to variations in model parameters, i.e., the coefficients in the model equations. For example, D' Agostino and Tecca (2006) compared FLO-2D (O'Brien and Garcia, 2009) simulations performed with two sets of rheological parameters and three values of the laminar coefficient K (six simulations in total), concluding that K controls the flood area and that rheological parameters control the maximum depth. Boniello et al. (2010) compared FLO-2D model results from a set of 12 back-calculated rheological parameters selected from previous studies, with another set of parameter values obtained from laboratory rheological analyses, finding a better representation of debris flow behavior with back-calculated parameters. Chow et al. (2018) conducted simulations with FLO-2D using 26 different sets of rheological parameters obtained from previous studies, combined with different values of volumetric sediment concentration C v , specific gravity G s and the surface detention SD, a parameter used in FLO-2D to represent flow detention. They found that the most important parameters were C v , SD and β 1 , which characterizes fluid viscosity. All these studies used fixed sets of rheological parameters in their numerical experiments, and, therefore, the relative importance of such coefficients on relevant simulated variables -specifically, flow depth, flow velocity, deposit volume and flood area -remains unknown. Therefore, this paper addresses the following questions: 1. How sensitive are debris flow model results to uncertain -and typically fixed -rheological parameters?
2. What are the most effective post-event measurements to constrain the parameter search towards more realistic simulations?
To answer these questions, we perform a sensitivity analysis on the parameters of a numerical debris flow model and examine the effects of using post-event in situ measurements on expected parameter ranges. In particular, we analyze a debris flow event that occurred in two creeks located in the Atacama region (northern Chile; ∼ 29 • S, 70 • W) during March 2015. This event was the consequence of heavy precipitation over a 3 d period, which exceeded 60 mm at several locations (Bozkurt et al., 2016), resulting in the loss of human lives and massive infrastructure damage. Therefore, our intention is to provide guidance on the choice of uncertain rheological parameters, contributing to more reliable numerical simulations for debris flow risk assessments and land use planning. The remainder of this paper is organized as follows: Sect. 2 describes the case study creeks and data; Sect. 3 describes the numerical debris flow model, sensitivity analysis and parameter search strategies; Sect. 4 presents the results and discussion; and Sect. 5 summarizes our main findings.

Study domain and data
We choose two ephemeral creeks located nearby in the upper Huasco River basin, Acerillas and La Mesilla (Fig. 1a), where debris flows were triggered by an extreme precipitation event on 24-26 March 2015 (Bozkurt et al., 2016;Ortega et al., 2019). The Huasco Valley is a semi-arid fluvial system located at the southern edge of the Atacama region, Chile. This valley is characterized by perennial rivers that only exist in the trunk valleys, while tributaries only show ephemeral streams. In these areas, heavy rainfall events may induce catastrophic debris flows and mud floods that greatly contribute to erosion . The Acerillas creek (15 km 2 basin area; Fig. 1b) has a markedly narrow channel with almost no alluvial fan, allowing for the transportation of sediments towards the El Carmen River. Post-event measurements indicate a deposited sediment volume of 6000 m 3 ) and a maximum flood area of 37 000 m 2 . Conversely, the La Mesilla creek (2.5 km 2 basin area; Fig. 1c) is characterized by a big alluvial fan where considerable sedimentation occurs, and postevent measurements show a deposited sediment volume of 102 000 m 3  and a maximum flood area of 246 500 m 2 . These flood areas were estimated by comparing pre-and post-event satellite Google Earth imagery. Also, a post-event topography lidar scan (acquired in February-March 2017 by the Chilean Ministry of Public Works) was available for this study. This dataset has a 1 × 1 m 2 horizontal resolution and was post-processed in order to eliminate vegetation and buildings.
Flow discharge data at the outlet of each creek were obtained from a distributed hydrological model, HEC-HMS (Hydrologic Engineering Center Hydrologic Modeling System) version 4.2 (USACE, 2015). The model was configured for the entire Huasco River basin (7242 km 2 ), upstream of the Santa Juana irrigation reservoir, as part of a debris flow mitigation project for the Chilean Ministry of Public Works. Hydrologic model simulations were forced using data from point measurements at 14 meteorological stations, spatially distributed with the inverse-distance-weighting (IDW) inter- polation method (Teegavarapu et al., 2006). Total rainfall records range from 20 to 76 mm, with a maximum registered intensity of 16 mm h −1 . The HEC-HMS model parameters were calibrated against hourly streamflow observed at two gauge stations located in the upper part of the basin -Río El Carmen en El Corral and Río Conay en Las Lozas (Fig. 2) -obtaining a Nash-Sutcliffe efficiency (NSE; Nash and Sutcliffe, 1970) of 0.78 for the former and 0.64 for the latter. Although all the other gauging stations were buried or destroyed by debris flows, simulated total water volumes were similar with those captured by the Santa Juana reservoir, whose levels were low before the event. Estimated peak flow discharges for the event analyzed are 8 m 3 s −1 at La Mesilla and 12.7 m 3 s −1 at Acerillas.

Debris flow model
We use the two-dimensional FLO-2D debris flow model (O'Brien et al., 1993), configured at a 10 m horizontal resolution. FLO-2D is a finite difference model that simulates water or debris flows in channels or unconfined surfaces. The governing equations solved by FLO-2D are the depth-averaged continuity and momentum conservation (Eqs. 1 and 2), and the flood wave progression is controlled by topography and flow resistance (O'Brien and Julien, 1988). FLO-2D can also simulate debris flow rheologies using a "quadratic" rheological model that combines components associated with creep, viscous, dispersive (collisions) and turbulent stresses (O'Brien and Julien, 1988;O'Brien and Garcia, 2009;Naef et al., 2006). Based on this quadratic rheology, the friction slope S f is estimated as (Eq. 3): where h is the local flow depth, t is time, V x and V y are depthaveraged velocity components along the x and y directions, g is gravitational acceleration, S f is the friction slope, S o is the bed slope, τ y is the yield stress, K is a laminar strength parameter, η is the interstitial fluid dynamic viscosity, and n td is the conventional Manning's roughness coefficient corrected by C v (n td = 0.0538ne 6.0896C v ; O' Brien and Julien, 1988). O'Brien et al. (1993) proposed the following empirical relationships to calculate the viscosity and yield stress as a function of the volumetric sediment concentration c v (Eqs. 4 and 5): where α 1,2 and β 1,2 are experimentally defined empirical coefficients (O'Brien et al., 1993;O'Brien and Garcia, 2009).

Sensitivity analysis
Sensitivity analysis (SA) is a powerful tool to characterize the effects of variations in input factors on environmentalmodel responses (Razavi and Gupta, 2015;Gupta and Razavi, 2018). When the factors of interest are the model parameters, SA helps to identify those that are redundant for the modeling purposes, contributing to a more efficient parameter search (Mendoza et al., 2015). Different types of SA techniques have been proposed in the literature depending on specific objectives or even the meaning of sensitivity (see, for example, reviews by Pianosi et al., 2016). In this work, we apply the Distributed Evaluation of Local Sensitivity Analysis (DELSA) method (Rakovec et al., 2014), which is a frugal local-global hybrid technique, to identify the parameters that have the largest impact on simulated debris flow variables. Although our implementation only examines first-order sensitivities across the parameter space -as in Rakovec et al. (2014) -it should be noted that DELSA has considerable unexplored potential to characterize parameter interactions, which could be achieved by including additional terms in the prediction total variance, as suggested by Sobol and Kucherenko (2010).
First-order sensitivities are obtained using local gradients that quantify the sensitivity of a modeled output relative to individual variations of a parameter θ j . The local gradients ∂ ∂θ j | k are used to compute the first-order sensitivity of each parameter j at each point k of the parameter space: where V K ( ) is the total local variance at point k The first-order sensitivity measures S1 j k vary between 0 and 1, and the sum of first-order sensitivities from all parameters is equal to 1 at each sample point. In this work, parameter sampling is performed using the Latin hypercube sampling (LHS) method. LHS is a statistical method to generate an almost random sample of parameter values from a multidimensional distribution and has proven to be more efficient than other methods like Monte Carlo sampling (Olsson et al., 2003;Olsson and Sandberg, 2002).
Local sensitivities can be analyzed in a disaggregated manner -through their cumulative frequency distribution across the parameter space -or aggregated by computing a specific statistical property, e.g., the median of all local sensitivity measures for a particular pair of parameter and target variables. We use both approaches to analyze SA results.
In this paper, we focus on the effects of debris flow model parameters on four response variables: maximum average runoff speed V mean (m s −1 ), maximum average runoff height H mean (m), maximum flood area A max (m 2 ) and deposited volume Vol dep (m 3 ). These response variables are calculated using the outputs from FLO-2D as where j is the cell index; N WC is the total number of wet cells (h > 0); dx and dy indicate the cell size along the x and y axis (numerical grid); and V max (j ), H max (j ) and H final (j ) are the maximum flow speed, maximum flow depth and final runoff height of the cell j , respectively. The FLO-2D parameters considered for DELSA are those that describe the fluid rheology (Table 1): α 1,2 , β 1,2 , C v , K, n, the SD parameter and the total volume of sediments mobilized Vol T .
The detention coefficient SD is a model parameter that controls flow detention. The FLO-2D user's manual and previous studies (D'Agostino and Tecca, 2006) suggest that SD acts as the minimum physically plausible flow depth (i.e., flow stops if flow depth < SD). Further, D' Agostino and Tecca (2006) noted that this coefficient has a strong influence on the results, and it can be used as a surrogate of the rheology. The estimated total sediment volume mobilized Vol Test by each debris flow event was calculated with the equation proposed by Chang et al. (2011) (Eq. 12).
where A w is the watershed area, A L is the landslide area (zero in these cases), GI is the geological index (where a value of 2.5 is assumed based on a study zone report made for the Chilean Ministry of Public Works), D is the rainfall duration (D = 48 h for this event) and C R is the cumulative rainfall (C R = 76 mm). We obtain Vol Test values of 185 000 m 3 for Acerillas and 154 000 m 3 for La Mesilla. Debris flow concentration is assumed to vary with streamflow between a minimum concentration C v min = 0.1-0.4 and a maximum concentration C v max at the time of peak flow, which is treated as a model parameter. To this end, we propose the following function for C v : where φ is a coefficient that changes the shape of the concentration curve. φ and C v min are calculated in order to match the total volume and minimize C v min value: Rheological parameter ranges are obtained from previous debris flow studies (O'Brien and Julien, 1988;Boniello et al., 2010;D'Agostino and Tecca, 2006;Sosio et al., 2007;Rickenmann et al., 2006;Chang et al., 2011;O'Brien and Garcia, 2009) and are summarized in Table 1. However, additional restrictions are imposed for τ y and η, with maximum values of 35 000 dyn cm −2 (dyne; g cm s −2 ) and 100 000 P (poise; 0.1 kg m −1 s −1 ), respectively . Since τ y and η are functions of rheological parameters (Eqs. 5 and 4), such limits impose restrictions for α 1,2 , β 1,2 and C v max that we ensure are implemented in DELSA. To this end, we develop a Python script that allows for running FLO-2D in parallel and sequentially, reducing computational cost considerably. For Vol T , we assume a range of variation of ±30 % with respect to estimated values (Eq. 12).

Parameter selection via constrained search
We explore the effects of parameter uncertainty on simulated debris flow variables at the two case study creeks. We also examine the utility of using reference values for specific variables to constrain the search of physically plausible parameter sets. Such values are obtained from a reference simulation conducted by Zegers (2017), who reproduced flood area and sediment volume observed during the 2015 debris flow events at Acerillas and La Mesilla creeks using FLO-2D. Zegers (2017) calibrated model parameters by contrasting results against measured flood areas, deposited volumes and flow velocity estimated from a video captured with a cell phone by a local person at Acerillas. Such a validation strategy has provided reliable results for several other creeks in the area. Parameter values used for the reference simulation are provided in Table 2. Modeled flood areas and maximum flow depth are shown in Fig. 1, and discrepancies with respect to observations can be attributed to the use of postevent topography. Based on the reference simulation, we choose reference values of V mean = 1 m s −1 and H mean = 1.5 m at both creeks. Additionally, we estimate the reference maximum flood area using Google Earth imagery -obtaining values of 246 500 m 2 for La Mesilla and 37 000 m 2 for Acerillasand use deposited sediment volumes reported by Cabré et al. (2020) as reference values, which correspond to 102 000 and 6000 m 3 for La Mesilla and Acerillas, respectively.

Choice of sample size
First, we seek to identify the minimum sample size N k for which stable DELSA results can be obtained in order to minimize computational cost (Rakovec et al., 2014). Therefore, we explore the effects of the choice of N k on the cumulative frequency distributions (CDFs) of DELSA firstorder sensitivity indices. Since we include N j = 9 parameters, the total number of simulations required for each case is N t = (N j + 1)N k . Figure 3 illustrates the sensitivity of  DELSA results to variations in sample size N k for four variables simulated by FLO-2D, with respect to parameters β 1 and C v max , at the Acerillas creek. Since the curves obtained for N k = 500 and N k = 1000 show slight differences, departing from the CDF for N k = 100, we conclude that a sample size of N k = 500 is adequate for further analyses. The sensitivity of DELSA results to N k was also examined at La Mesilla creek, obtaining the same conclusion regarding sample size (not shown). Figure 3 also shows that, depending on the target variable and parameter analyzed, first-order sensitivity indices can be highly heterogeneous across the parameter space. In particular, the modeled response is highly sensitive to variations of β 1 , with first-order sensitivities larger than 0.2 for approximately 60 % of cases and sensitivities greater than 0.5 for 20 %-40 % depending on the variable analyzed. On the other hand, the modeled variables are less sensitive to C v max in most cases, with DELSA indices smaller than 0.1 for approximately 70 % of the parameter sets. Figure 4 displays the median of the full frequency distribution (obtained with N k = 500) of local first-order sensitivity indices for the two study domains: (i) Acerillas creek (ACE) and (ii) La Mesilla creek (MES). The uncertainty bands are obtained by performing bootstrapping with replacement (resampled 1000 times). In general, β 1 provides the largest sensitivities for the simulated variables analyzed, which means that the fluid rheology, in particular the viscosity coefficient (η = η(β 1 )), is a main parameter controlling flow behavior.

Sensitivity of model responses to model parameters
Moreover, the results of Acerillas seem to be more sensitive than those obtained at La Mesilla with respect to β 1 . This could be better explained when analyzed together with SD, the detention coefficient.
As expected, the simulated deposited volume Vol dep is very sensitive to SD because this parameter controls flow detention. For the remaining simulated variables, SD also rises as an important parameter at La Mesilla, but it shows secondary importance at Acerillas, which can be explained by catchment differences. While the fluid rheology explains flow behavior at Acerillas creek (mainly sensitive to β 1 ), depositional or detention processes -represented by SD -gain importance across the larger alluvial fan of La Mesilla creek.
Another parameter that provides large sensitivities -especially in simulated mean flow velocity V mean -is the total sediment volume Vol T , whose sensitivities have the same order of magnitude as those produced by β 1 . The high sensitivity of Vol T on V mean is explained because the former influences fluid rheology through C v (t) (see Eqs. 4, 5 and 13). However, Vol T does not produce large variations in the total deposited volume, which could be explained because, as deposition occurs, the flow is channelized between the deposited margins, preventing flow spreading. For example, when increasing SD values, A max decreases as the flow is forced to stop at deeper heights. This could be a structural weakness of FLO-2D, which lacks proper representations of complex depositional and dewatering processes.
The large sensitivities in model response to variations of β 1 suggest that the viscous stress (second term in Eq. 3) is the main contributor to sensitivities in simulated frictional slope. On the other hand, DELSA sensitivity indices asso-  ciated with yield stress τ y and Manning's roughness coefficient -the other components of the frictional slope -are of second-order importance. Even more, model results are practically insensitive to Manning's roughness coefficient.
Our results show some differences with related previous work. D 'Agostino and Tecca (2006) compared FLO-2D simulations performed with two sets of rheological parameters and three values of K (i.e., six simulations), concluding that the latter controls the flood area, while rheological parameters control the maximum flow depth. On the other hand, our results indicate that the laminar coefficient K provides small sensitivities in simulated flood areas and that the most influential rheological parameter on simulated maximum flow depth is β 1 . Although D' Agostino and Tecca (2006) provided the first insights into FLO-2D parameter sensitivities, the number of parameters involved and the sample size were not large enough to draw robust conclusions. Recently, Chow et al. (2018) conducted simulations with FLO-2D using 26 different sets of rheological parameters obtained from 45 previous studies. They found that the most influential parameters -in order of importance -were C v , SD and β 1 . These results are different from ours due to discrepancies in the parameter sampling method, the use of fixed sets of rheological parameters and their parameter ranking definition, which does not consider separate effects on key simulated variables. For example, Chow et al. (2018) ignored the effect of the total sediment volume when analyzing C v , whereas we examine the separate effects of maximum sediment concentration and total sediment volume, obtaining that the latter is more important than the maximum sediment concentration. Figure 5 shows the full range of variation in model responses (produced by N k = 500 parameter sets), normalized by reference values obtained from post-event measurements (Vol dep and A max ) and results obtained from the reference models (for H mean and V mean ). These are compared with the simulated ensemble that results from screening model outputs imposing five different observational constraints: (i) ±20 % reference mean flow velocity (FVEL), (ii) ±20 % reference mean flow depth (FH), (iii) ±20 % reference maximum flood area (FAREA), (iv) ±40 % reference volume deposit (FVOL), and (v) ±20 % reference maximum flood area and ±40 % reference volume deposit (FAREAVOL). To be clear, constraint (i) results from keeping all those parameter sets that provide a simulated mean flow velocity within the range 0.8-1.2 V ref . Constraints (ii)-(iv) work in a similar way for other observed variables, while constraint (v) filters all parameter sets that simultaneously provide flood areas and deposit volumes within the ranges 0.8-1.2 A ref and 0.6-1.4 Vol ref , respectively. We assume a weaker observational constraint in the case of the deposited volume because of possible uncertainties associated with different measurements techniques. Figure 5 also shows that, in general, the effects of parametric uncertainty on simulated variables are considerable and the ensemble median can be substantially different from the reference boundaries. This is somewhat expected, since the literature provides large ranges for model parameters (see Table 1 for details). Overall, uncertainties arising from the original parameter samples (top panels) are larger in Acerillas, except for mean flow velocity. This could be explained by the larger sensitivity of flow velocity to flow rheology, while the rest of simulated variables are more sensitive to deposition processes, mainly represented by SD.

Parameter uncertainty effects
Most simulations overestimate deposited volumes and flow depth, especially at Acerillas. For flow velocity, the ensemble of parameter sets provides mixed results in both creeks, with underestimation in most cases (median values lower than the reference values); however, there are still several parameter sets that produce an overestimation of flow velocity. The results obtained for maximum area reveal differences among both creeks: in Acerillas, most parameter sets tend to overestimate the flood area, whereas most simulated values are within the expected range at La Mesilla. This could happen because Vol dep Vol T at Acerillas, while Vol dep ∼ Vol T at La Mesilla; moreover, the maximum flood area at La Mesilla is approximately 6 times larger than in Acerillas. Thus, small variations in Vol dep imply important fractional changes with respect to the volume reference values at Acerillas. Figure 5 shows that the velocity constraint FVEL does not have an impact on the rest of simulated variables. Nevertheless, the application of alternative observational constraints helps to reduce the spread of the remaining variables. For example, the height filter FH improves simulations of maximum flood area, although it does not have much effect on velocity or deposit volume. The area restriction FAREA improves simulated flow depth only at Acerillas, as most of the original ensemble members were already inside the expected reference boundaries at La Mesilla. The volume constraint FVOL reduces the uncertainty in all variables, with the smallest improvement for flow velocity. This is because the volume is directly linked to flow height and flood area, but not to flow velocity. Finally, the largest reductions in ensemble spread are obtained when parameter sets are constrained by using area and volume observations (FAREAVOL).
The maximum flood area and deposited volume are relatively easy to measure and are probably the most used post-event measurements for calibrating debris flow models (Chow et al., 2018;Cesca and D'Agostino, 2008;D'Agostino and Tecca, 2006;Sosio et al., 2007;Frey et al., 2016;Quan Luna et al., 2011). Figure 6 illustrates the effects of applying observational constraints, specifically flood area and deposited volume constraints, on parameter identifiability. Results show that the re- sampled values of α 1,2 , β 2 , C v max , Vol T , K and n cover practically the entire original range. However, the application of observational restrictions provides substantial reductions in the ranges of β 1 and SD, which mainly explain the flow rheology and depositional processes in our study areas. Further, lower parameter values are obtained in comparison to the full range, especially SD at the Acerillas creek. These results indicate that viscosity η = η(α 1 , β 1 , C v max ) is the most restricted parameter when applying these constraints, discarding all medium-high values. On the other hand, resampled values of τ y = τ y (α 2 , β 2 , C v max ) cover almost the entire original range.

Parameter identifiability
The reference SD value is close to the upper range in La Mesilla after applying the FAREAVOL constraint, and also much larger than the resulting maximum values filtered at Acerillas. This result is somewhat expected, since SD does not provide large model sensitivities in that domain. Similarly, the reference β 1 value is within the upper range of filtered values (> Q 75 ). However, filtered η values around the baseline model parameter result from the compensation of low values of α 1 (near the minimum in the filtered range) and C v max (below the median of filtered values). A similar effect is observed for τ y , whose filtered range results from the compensation of α 2 , β 2 and C v max . In summary, different combinations of α 1,2 , β 1,2 and C v can generate viscosity and yield stress values that are suitable to reproduce the 2015 debris flow events in Acerillas and La Mesilla. This is a well-known problem in environmental models -referred to as the equifinality, nonuniqueness or nonidentifiability of model parameters -that has been widely discussed for more than 3 decades in the hydrology literature (e.g., Beven, 2006;Kelleher et al., 2015) but not carefully addressed in the debris flow modeling community. Figure 6 also shows different behavior in other parameters. For example, the reference values for K and n are in the lower range of the filtered ensembles, while the reference Vol T is in the upper body of the boxplot. Low K values produce low S f 2 , the second term at the right hand of Eq. (3), representing viscous stress. This could be compensated for using larger values of SD, as in the reference model. These results demonstrate that equifinality in FLO-2D does not only involve rheological parameters and that SD could be an important parameter to correct the unrealistic model representation of rheology (D'Agostino and Tecca, 2006).
The main goal of this study was to characterize the sensitivity of model responses to variations in uncertain rheo- Figure 6. Effects of applying a flood area-volume observational constraint (FAREAVOL) on parameter identifiability. Results are displayed for Acerillas (black) and La Mesilla (red) creeks. The vertical bold line in the boxplots is the median; the body of each boxplot shows the interquartile range (Q 75 -Q 25 ); and the whiskers represent the sample minima and sample maxima. The grey-black diamonds represent parameter values for the reference simulations. logical parameters, using only independent information on parameter values (i.e., the situation comparable to Sobol, as proposed by Rakovec et al., 2014). Hence, the only verification dataset available (flood area and sediment volume from the March 2015 event) was used to examine the identifiability of model parameters. However, the relative importance of additional observations could be assessed through the Observation-Prediction (OPR) statistic (Tiedeman and Hsieh, 2004), and the potential new information provided by field data (e.g., sedimentological and morphological characteristics) for a specific parameter (e.g., α 1 and β 1 ) could be quantified with the Parameter-Prediction (PPR) statistic (Tonkin et al., 2007). It should be noted that, in both cases, the equation for the total local variance (Eq. 7) would be different, as additional information should be incorporated (see Appendix A in Rakovec et al., 2014).

Summary and conclusions
We performed a sensitivity analysis on the parameters of a widely used numerical debris flow model (FLO-2D) and assessed the effects of applying observational constraints on parameter identifiability. Our study domains are two morphologically different ravines, Acerillas and La Mesilla, located in the Atacama region (∼ 29 • S, 70 • W), Chile. While Acerillas is characterized by a straight and well-defined channel with almost non-alluvial fans, La Mesilla has a big alluvial fan where deposition is prone to occur.
We found that β 1 -a parameter used to estimate the fluid mixture viscosity -provides the largest sensitivities in the variables analyzed, followed by SD -a model parameter used to represent flow detention. Interestingly, the relative importance of β 1 and SD depends on the study site, with the former being more important for Acerillas and the latter being more important for La Mesilla. These results suggest that, while rheological processes dominate flow behavior at Acerillas (straight channel with small alluvial fan), sedimentation and detention processes control flow in La Mesilla (big alluvial fan). Although the total mobilized sediment Vol T does not have an effect on Vol dep , it is important for representing flow velocity, as Vol T is used to estimate C v (t), a key parameter for fluid rheology. Finally, model results seem to be almost insensitive to Manning's roughness coefficient n, while DELSA sensitivities for the remaining parameters are of second-order importance and provide similar indices.
The comparison between the original model parameter ranges (N = 500) and the ensemble resulting from applying observational restrictions shows that SD and β 1 (i.e., η) are the parameters whose identifiability is mostly improved, while others practically preserve their original range. In addition, we obtain that different combinations of model parameters (including those that describe rheology) can provide very similar results, indicating that equifinality is present in FLO-2D. Our results also support the idea that single-phase rheological models lack a strong physical basis (Iverson, 2003), and, therefore, their determination requires expert knowledge. However, an encouraging finding is that the final deposited volume Vol dep and maximum flood area A max contain considerable information for identifying model parameters.
We obtain that SD strongly affects model results at La Mesilla, having also large effects on simulated deposited volumes at Acerillas. Moreover, this study provides evidence that SD is one of the most important parameters controlling flow behavior and could possibly surrogate rheology in the model (D'Agostino and Tecca, 2006). One-phase debris flow models still lack robust representations of complex process interactions during flow stopping that produce temporal and spatial changes in fluid rheology. Thus, these rheological changes have been replaced by simpler approaches (e.g., the incorporation of SD).
Future investigations should aim to improve the structure of debris flow models and hence achieve better simulations of deposition and erosion processes, stopping phases, and changing rheologies. Further, the development of computationally frugal methods to improve the understanding of parameter interactions in environmental models emerges as an attractive avenue for future research.