Abstract

In this study, a simulation tool, based on the combination of Geographic Information System (GIS) and Grid-Based Regional Slope-Stability Model (TRIGRS), is developed to assess and predict the snowmelt-induced landslides in areas of sensitive marine clays in the Ottawa region (Canada). Topographic, geologic, hydrologic, and geotechnical information of the study area, in addition to snowmelt intensity data for different periods, was required to conduct this modeling study. Snowmelt intensity records for periods of 6–48 h, 3–15 days, 25 days, and 30 days, as well as the information on historical landslides in the study area, were used to examine both the timing and location of shallow landslides due to snowmelt across the Ottawa region. The developed tool is validated by comparing the predicted landslide-susceptible areas with the historical landslide maps in the study area. A good agreement between the predicted and recorded historical landslides was obtained, which suggests that the developed GIS-TRIGRS based tool can predict relatively well the snowmelt-induced landslide susceptibility in the sensitive marine clays. The modeling results show that high slope areas of sensitive marine clays are more prone to snowmelt-induced landslides. As a result, this simulation tool could be used to assess and/or predict snowmelt-induced landslides in different areas of the Ottawa region, which is an important means for decision-making processes with respect to protecting the population and infrastructure against snow-melt induced marine clay landslides.

Introduction

There are many types of soils that cause geotechnical problems and geoenvironmental hazards or disasters. Sensitive marine clay is one of these soils; they will significantly lose their shear strength if their structure is disturbed. During the last glacial era, glaciers covered many areas of the world, especially in the northern regions of the earth (e.g, Alaska, Canada, Sweden, and Norway). These areas are found to contain these problematic soils (Theenathayarl 2015). The provinces of Quebec and Ontario (particularly in the Ottawa region) in Canada have sensitive marine clays (locally called Leda clay or Champlain sea clay), which cover relatively large areas. These Canadian marine clays, younger than 12,000 years, are considered to be young glacial deposits (Haché et al. 2015; Taha, 2010).

Several landslides have previously occurred in these soils and in many countries as reported in many studies (e.g., L'Heureux 2013; Aylsworth et al. 1997). For instance, numerous landslides have occurred in the Canadian sensitive marine clay formations, particularly in the Ottawa region. A key trigger of these landslides has been identified as snow melting in the spring (Quinn 2009). During the snowmelt in the spring, huge quantities of snowmelt water is released, which has a negative effect on the stability of the slopes.

As the population of the Ottawa region is steadily rising (currently at approximately 900,000) (City of Ottawa 2015) there is a continuous growth in infrastructure facilities, including residential and other buildings, transportation infrastructure, and utilities located on the problematic marine clay. However, the presence of sensitive clay slopes in the Ottawa region and presence of snowmelt in the spring, poses a landslide hazard for both the population and infrastructure (Nader et al. 2015; Taha, 2010; Quinn 2009). A good engineering practice to reduce this landslide hazard and potential impact on people and infrastructure is to identify or assess the sensitive marine slopes that are susceptible to snowmelt-induced landslides as well as develop a tool that can assess and predict these snowmelt-induced landslides in Ottawa sensitive marine clays. Indeed, landslide susceptibility mapping is an important tool for decision-making processes with respect to protecting the population and infrastructure against marine clay landslides (Nader et al. 2015). However, until now, no tool has been developed to assess and map the sensitive clays slopes in Ottawa that are susceptible to landslides induced by snowmelt. This knowledge and technology gap is addressed in the present study.

It should be emphasized that studies on snowmelt-induced landslides are rare compared to rainfall-induced landslides. For instance, in the Clearwater River Basin in Central Idaho, Gorsevski et al. (2000) used a Digital Elevation Model (DEM) to derive analyses of snowmelt and rainfall triggered landslides on a small watershed named Rocky Point. Kimura et al. (2014) used a probabilistic approach to estimate landslide hazard. According to the topography for the study area, GIS analysis was used to classify these landslides into three groups. Khezri et al. (2013) used a GIS system to predict landslide susceptibility in the central Zabbasin of the West-Azerbaijan province in Iran. This study, claims that the key factor of mass movement is the sudden melting of snow as a result of the presence of water in the soil, increasing the pore water pressure and hydrostatic level. Kawagoe et al. (2009) analyzed snowmelt’s effect on landslide occurrence by using a probabilistic model based on multiple logistic regression analysis. The analysis was executed by using a logistic regression model to produce landslide susceptibility maps across Japan. The results revealed that across Western Japan (throughout the mountain ranges near the Japan Sea), there exists a 95% landslide probability. Ayalew et al. (2004) created a landslide susceptibility map near the Agano River in Tsugawa using a GIS based analysis of the spatial database collected for 791 landslides. Most often, landslides were triggered by heavy amounts of snow and rainfall. However, no studies on snowmelt induced landslides in the sensitive marine clays in the Ottawa have been conducted, and there is no tool to assess or predict the snowmelt-induced landslide susceptibility for the Ottawa region.

Several approaches have been established to evaluate landslide susceptibility and identify the potential of unstable slopes/regions. Such approaches can be divided into three key groups: (i) expert evaluation (i.e. inventory, heuristic approaches), (ii) statistical approaches, and (iii) deterministic approaches (Pathak 2016; Quinn 2009; Fall et al. 2006; Van Westen et al. 1997). The advantages and disadvantages of each technique have been outlined in numerous publications (e.g. Fall et al. 2006; Dai and Lee 2001). Moreover, landslide susceptibility mapping in large areas usually involves the handling, processing, and interpretation of a large amount of geospatial data. Therefore, during the past decades, the Geographical Information System (GIS) has been successfully used in several studies (Zhang and Wang 2019; Chen et al. 2015; Mukhlisin et al. 2010; Quinn 2009; Thapa and Esaki 2007; Fall et al. 2006; Chang

et al., 1999) to conduct landslide susceptibility assessment and mapping. Additionally, the infinite slope stability model has also been merged with GIS to assess shallow landslide potential and the distribution of the factor of safety in a specified region (e.g., Sorooshian et al. 2015; Quinn 2009; Fall, 2006; Van Westen et al. 1997). Furthermore, Transient Rainfall Infiltration and Grid-based Regional Slope-stability (TRIGRS) (Baum et al. 2008), a deterministic, spatially distributed model that combines a hydrologic 1-D transient infiltration model (Iverson 2000) and an infinite slope stability model to calculate the factor of safety with respect to time and location, have been successfully applied in several studies for modeling the timing and spatial distributions of shallow landslides and to compute their factor of safety with respect to time and location (e.g. Raia et al., 2014; Park et al. 2013; Salciarini et al. 2008). The TRIGRS model, when used for either saturated or unsaturated soils, is appropriate for shallow landslide susceptibility analysis (Raia et al., 2014). The application of models like TRIGRS in a GIS environment for landslide susceptibility evaluation requires digital spatial topographic, geologic, and hydrologic data in addition to any previous snowmelt and rainfall-induced landslide accounts, which are fundamental to test the model results (Baum et al. 2008; Baum et al. 2002).

Based on the facts mentioned above, in this work, the TRIGRS model combined with GIS will be applied to the study area to assess snowmelt induced landslides, and to develop snowmelt induced susceptibility landslide maps. This study aims to develop a GIS-TRIGRS model or tool and map to assess the shallow and snowmelt-induced landslide susceptibility in sensitive marine clays in Ottawa.

Study area

The Ottawa region (Fig. 1) is located between the latitudes of 45.00–45.50 N and the longitudes 75.50–76.00 W. It is bounded to the north by the Ottawa River. It is also a region consisting of gentle slopes ranging from 0 to 33 degrees; although, Western and Northern Ottawa maintain higher slopes. Along the southern bank of the Ottawa River, the city of Ottawa can be found. The oldest portion of the city lies between the Rideau Canal and the Ottawa River.

Geological setting

The stratigraphy of the study area is characterized by deltaic and estuarine soil consisting of fine to medium grained silty sands of 1.8 to 6 m thickness underlain by the deposits of Ottawa sensitive marine clays, which were deposited within the Champlain Sea basin following previous glaciations. These deposits consist of glacial till overlain by relatively thick layers of sensitive marine silty clay deposits called Champlain Sea clay (also known as Leda Clay), which were deposited as a result of the recession of the Champlain Sea. This physiographic region is underlain by different types of sedimentary rocks, such as sandstones, dolostones, limestones and shales. In addition, the final deposits of silty sand overlie Champlain Sea clay. Organic soils (like peat) are also found in certain badly drained regions (Schut and Wilson 1987). Although some may not exist in some areas, layers of peat, sand, clay, and glacial till are common types of soil material that are found in the entire area. Within the marine clays, two major units can be distinguished: the upper and lower unit. The upper unit mainly consists of sand thinly interbedded with silty sand, and locally overlay layers of clay, silty clay and silt, generally the higher parts are mottled or laminated reddish-brown and bluish-gray. The blue-gray clays dominate, and the laminations are less frequent (Schut and Wilson 1987). The lower unit is blue-gray clay, sometimes mottled. Below the upper layer of clay are coarse galciofluvial and glacial lacustrine sediments that vary in thickness (up to one meter). These deposits lie atop relatively flat limestone (Schut and Wilson 1987). The surface soils are deposits that occur as a result of the recession of the Champlain Sea.

Climatic conditions

Ottawa has a continental humid climate with an average annual relative humidity of 63%. The average monthly relative humidity ranges from 53% in May, to 75% in January. High temperatures are present in July, where the average daily maximum temperature is 23.5 °C; comparatively, in winter, the mean daily minimum temperature is − 16.4 °C (Auld et al. 2009). The annual rainfall is 714 to 869 mm (Canadian Council of Professional Engineers 2008). In Ottawa, snow and ice are typically present throughout the year except in July and August (see Fig. 2).

Generally, the most snowfall in the Ottawa region occurs from November to March with a snowfall range of 208 mm annually (Auld et al. 2009). The highest daily snowfall of 383 mm was recorded in Ottawa region in the years 1960 and 1971. The duration of the winter and snow are not steady in a typical winter in Ottawa; long-lasting snow cover starts from mid-December until early April. Freeze-thaw cycles take place during the winter with a few days lying well above the freezing point, followed by nights well below it (Energy East Pipeline Project (EEPL) 2014). Freezing rain and high wind chills are also very common. Summers are relatively warm and humid in Ottawa, even though they are typically short.

Methodology

This study takes a multi-stage approach to evaluate landslide susceptibility induced by snowmelt. The developed methodology is summarized by the six main stages presented in the flowchart, Fig. 3.

The first stage consists of collecting the input data needed for this study. This data includes information about the geotechnical characteristics as well as the climate, hydrogeological, and spatial data of the study area (Table 1). The geotechnical data was obtained from geotechnical studies or reports completed by accredited companies in different sites within the city of Ottawa as well as by government agencies with extensive work experience throughout Ottawa (e.g., Inspec-Sol Inc. Engineering Solution 2014, Houlechevrier Engineering Geotechnical Ltd 2013, Kollaard Associates Engineers 2013, Trow Associates Inc. 2010, Stantec Consulting Ltd. 2010, Golder Associates Ltd. 2008). The data retrieved from previous geotechnical investigations include standard penetration test (SPT) results, consistency limits (liquid limit, plastic limit, etc.), shear strength parameters, moisture content, depth, thickness, and the relevant types of subsurface layers (Table 1). Data for Ottawa snowmelt was gathered from Environment Canada. Snowmelt data was recorded at the Ottawa CAD station located in Ottawa from 1890 to 2007. The bar chart, shown in Fig. 4, represents the relationship between the time of a recorded landslide and the corresponding snowmelt intensity. Greatest snowmelt intensity was recorded in 1960 and 1971. Snowmelt intensities were selected for periods of 6–48 h, 3–15, 25, and 30 days.

By monitoring several hundred wells, hydrogeological data was collected by the Provincial Groundwater Monitoring Network (PGMN) Program. The PGMN program collected initial groundwater levels of the 474 wells by observing and monitoring groundwater levels on an hourly basis. The wells monitored were located at the interface between the soil and bedrock, allowing for easy retrieval of the pore water of each aquifer. Using each well’s latitude and longitude, the location and elevation of each well were identified using the Provincial Digital Elevation Model (DEM) (Rogojin 2014). For this GIS-TRIGRS model, groundwater level recordings were used for March since snowmelt begins in spring in Ottawa.

Moreover, a suitable resolution and raster scale of 1:50,000 was used in the Digital Elevation Model (DEMs) published by the Natural Resources Canada. The DEM was applied to the GIS model to evaluate shallow landslide susceptibility modeling. The original DEM raster resolution of ten by ten meters was used, as applied to the high relief terrain evident in the study area. It should be underlined that the finest DTM resolution does not necessary provide the most accurate results because too localized topography does not represent the processes governing landslide initiation (Catani et al. 2010). The terrain’s slope was calculated using the Spatial Analyst in ArcGIS (see Fig. 5a). The lower boundary raster was calculated by using the following exponential eq. (1) by Godt et al. (2008):

Where dIb represents the depth of failure (see Fig. 5b), δ is the slope angle. Additionally, raster images of the geological distribution of the sensitive clays in the study area were also generated (see Fig. 5c). The resolution of all created raster images was then transformed to 25 by 25 m to fit the TRIGRS requirement. Later, raster images for the acquired geotechnical information, and for snowmelt intensity for all periods, was generated for the entire study area. Furthermore, the digital data of the spatial distributions of the historical landslides identified in the study area was also acquired. This data is needed to test the landslide prediction ability of the developed model.

In the second stage, all data obtained from the previous stage was converted to the ASCII file format in order to run the TRIGRS model including the raster images of geotechnical parameters (e.g., ϕ′, c′, θs, θr, α), rainfall intensity, slope, depth of failure, hydraulic conductivity (Ks), and soil water diffusivity (DѰ) (see Table 1). This format is convenient for the TRIGRS model. Data conversions were performed using the ESRI ArcGIS program. The above data was entered into the TRIGRS model using the Baum et al. (2008) method.

In the third stage, the TRIGRS model was run to assess the stability of the slope in the study area, for various snowmelt durations and intensities by using the infinite slope model to compute the safety factor, as described below.

In the fourth stage, the obtained result files from the TRIGRS model (in ASCII file format) were converted back to raster format to work in GIS. Each raster was reclassified; including the number of pixels recorded using the ArcMap program. This allows for the mapping of the distribution of unstable slopes in the region. This distribution data is based on the safety factors calculated for the snowmelt durations and intensities studied in the area (unstable slopes maintain a FS < 1).

In the last stage, the modeled landslide susceptibility map was validated by comparing the predicted landslide areas, with the locations of historical landslides. An analysis of the sensitivity of the TRIGRS-GIS model to the variations in the input data was also performed.

Infinite slope stability model

The infinite slope stability model is used alongside TRIGRS to analyze the stability of the slopes. It is assumed that the landslides in the Ottawa sensitive clays have a shallow depth relative to their length and width. Figure 6 depicts a schematic representation of the conceptual infinite slope model of the study area to compute the factor of safety (and thus to assess slope stability) by using a TRIGRS model. Upon landslide occurrence, the upper part of the slope is typically left unsaturated, as observed in the previous geotechnical studies in Ottawa region (e.g., Aylsworth et al. 1977). The initial ground water table in the study area was reconstructed based on the observations of groundwater levels in wells within PGMN described previously and is interpolated in ArcView GIS analysis. Average water levels or water depths from the ground surface were extracted from the wells within 40 km boundary of the City of Ottawa limit. The average water levels were between 40 mASL and 200 mASL in this study. The water levels were then compared to ground surface to extract depth of groundwater level for the area. The equation to compute the safety factor by snowmelt water infiltration at a depth, Z, for the infinite slope model (see Fig. 6) is given by (Salciarini et al. 2008; Terzaghi 1943):

Where c′ is soil cohesion, ϕ′ is soil friction angle, Ѱ is groundwater pressure head as a function of depth Z and time t, δ is slope angle, and γw, and γs are the unit weights of water and soil, respectively. The infinite slope is stable when FS > 1, in an equilibrium state when FS =1, and FS < 1 indicates unstable conditions. Thus, the depth of landslide initiation is the depth, Z, where FS reaches 1.

Hydraulic factors

The infiltration of snowmelt water was used as the main parameter to represent the hydraulic boundary condition in this study, because the infiltration of snowmelt water results in a change in the hydraulic gradient. The hydraulic gradient (defined by head loss, h, over horizontal length, L) is derived by analyzing the soil during unsaturated infiltration and by using Richard’s equation (Cobin 2013; Kawagoe et al. 2010). The unsaturated infiltration analysis uses soil data, slope angle, and snowmelt (represents the snow water quantity) as the main input data.

Snowmelt estimation

Several methodologies can be used to quantify snowmelt. The degree-day method was used in the present study to estimate snowmelt via the degree-day equation below (Kawagoe et al. 2010; Kazama et al. 2008; Hultstrand et al. 2006).

$$ \mathrm{SM}=\mathrm{K}\times \mathrm{T} $$

(3)

Where K is the degree-day parameter (mm/°C/day) and T is the mean daily temperature (°C).

If T is less than zero degrees Celsius, snowmelt is not evaluated (because water remains solid) and the negative snowmelt value is discarded. The equation contains only one degree-day parameter, making it easy to determine the optimized degree-day parameter when correlating the estimated snow area with satellite images. With respect to input data, this method requires the following parameters to estimate snowmelt: precipitation, temperature, and elevation data. Precipitation and temperature directly affect snowmelt, because snow is generated from precipitation and temperature, and elevation directly affects the temperature of the snow (Clow et al. 2012; Kawagoe et al. 2010; Kazama et al. 2008).

Infiltration analysis

Infiltrability is known as the amount of water that can infiltrate soil as excess precipitation (known as runoff) continues to flow downhill (Godt et al. 2008; Baum et al. 2008; Raia et al. 2014). Routing surface runoff is necessary because the cells filled with excess surface water remain adjacent to downslope cells where the water can either infiltrate or flow farther downslope (known as infiltrability) (Godt et al., 2008, Baum et al. 2008). Routing surface runoff effectively can minimize the loss of excess water that cannot infiltrate the soil, and to improve performance in existing urbanized or impervious areas. When storm drains are used to divert stormwater, runoff routing is unnecessary. Runoff is assumed to occur when the precipitation and runoff supplied to a cell exceed its infiltrability (Raia et al. 2014; Godt et al. 2008; Baum et al. 2008).

The saturated hydraulic conductivity is represented by Ks. The infiltration, I, at each cell is computed as the sum of the precipitation, P, in addition to any runoff from upslope cells, Ru. The runoff is limited by infiltration and it cannot exceed the saturated hydraulic conductivity, Ks.

However, it should be underlined that the eq. 4 is a simplification, since infiltration exceeds Ks for short periods of time at the beginning of storm in the case of unsaturated soils (Cuomo and Della Sala 2013).

Provided the infiltration is less than the saturated hydraulic conductivity, Baum et al. (2008) argue that TRIGRS uses a simple method for routing surface runoff from cells with excess surface water, to adjacent cells where it can either infiltrate or flow further downslope (Baum et al. 2002). At each cell where P + Ru is greater than the saturated hydraulic conductivity, Ks, the excess flow is considered runoff, Ru, and is diverted to adjacent cells.

During analysis, overland flow is assumed to be instant; therefore, this study does not model the rate of overland flow. Storm duration during analysis should be long enough for excess surface runoff to flow downslope. The routing method used ensures mass balance throughout the storm. Ensuring mass balance confirms that the total precipitation in the cells is equal to the water infiltrating the soil cells and the excess water flowing outside the region’s border. TRIGRS models do not consider runoff betweem time steps nor do they track water exiting the system (via storm drainas). TRIGRS models for runoff analysis assumes that water infiltrates a different cell or reaches the model’s edge (Park et al. 2013; Baum et al. 2008). In addition to excess runoff from cells, water is also assumed to runoff from cells where the water table is initially at the ground surface, and the infiltration rate is initially steady and negative. In this case, water is leaving such cells and runs to adjacent downslope cells.

Alternatively, the TRIGRS model tracks water leaving cells in the mass balance calculations (Raia et al. 2014; Park et al. 2013; Baum et al. 2008). Instead of iterating to reach a solution satisfying mass balance, this model analyzes the infiltration and runoff of cells starting from the ground’s surface down to further depths, removing iteration from the analysis. In order to achieve mass balance successfully, the topographic data must be indexed correctly. Digital Elevation Models (DEM) do correctly index topographic data if GIS has already been adjusted so that the flow can be directed concurrent with the surrounding topography. Mass-balance calculations for any storm verify that water applied as an input, is accounted for as an output via infiltration or runoff (Raia et al. 2014; Park et al. 2013; Godt et al. 2008; Baum et al. 2008).

Hydrological model

Infiltration models for saturated initial conditions

Groundwater flow can be modeled in two states when initially in the saturated condition: transient or steady flow. Steady groundwater flow is dependent on the water table level, infiltration rate, saturated hydraulic conductivity, and slope angle. Steady water flow occurs in both the horizontal and vertical planes (Baum et al. 2008; Iverson 2000). The results obtained from the TRIGRS model are very sensitive to initial steady seepage conditions, so maintaining accurate initial conditions is essential to this analysis. In order to generate accurate recordings of the initial seepage conditions, field observations and steady flow models are imperative. If accurate initial conditions are lacking, TRIGRS is best applied to the modeling of hypothetical contexts (Godt et al. 2008; Baum et al. 2002).

Where t is time, Z = z/ cos δ Z is the vertical coordinate direction, z is the slope-normal coordinate direction, β = cos2δ − (IZLT/Ks), where Ks is the hydraulic conductivity in the Z direction, IZLT is the steady surface flux, and Inz is the surface flux of a given intensity of the nth time interval. The subscript LT denotes long term, D1 = D0/cos2δ, where D0 is the saturated hydraulic diffusivity, N is the total number of intervals, and H(t − tn) is the Heavy side step function tn is the time at nth time interval in the rainfall infiltration sequence (Godt et al. 2008). The function ierfc (η) is of form (Godt et al. 2008).

The first term in Eq. (8) represents steady flow, whereas all other terms represent the transient flow. When using the TRIGRS model, limits can be applied to make the maximum pressure head under gravity’s pressure not exceed the pressure head at the ground surface if the original flow direction and hydraulic gradient are maintained (Godt et al. 2008; Baum et al. 2002). This model is presented in Eq. 10:

$$ \Psi \left(\mathrm{Z},\mathrm{t}\right)\le \mathrm{Z}\beta $$

(10)

Infiltration models for unsaturated initial conditions

Unsaturated flows required an analytical solution to use the TRIGRS model to scenarios following initial conditions. This method analyzes soil as a two-layered system containing a saturated and unsaturated zone with a capillary fringe (Baum et al. 2002). The unsaturated zone is assumed to soak up the water penetrating the surface until the water table rises beyond its initial state, while also reducing surface infiltration to the deeper ground due to low permeability. This model is dependent on several input parameters to produce the soil water characteristic curve (SWCC) for the wetting of unsaturated soil: residual water content, volumetric water content, saturated hydraulic conductivity, and the alpha parameter (Srivastava and Yeh 1991). As water infiltrates the unsaturated soil, it will distribute downwards towards the bottom of the region, causing the water table to rise and pore pressure to increase. These pressure changes and waves spread quickly in very thin saturated zones. Models outside of the TRIGRS model can be used to model the transmission of the pressure waves across the saturated zones. For example, a one-dimensional form of Richard’s equation can be used to depict the infiltration of the ground water through the unsaturated zone (Baum et al. 2008). To explain the effect of a sloping ground surface on infiltration, Eq. (11) by Iverson (2000) can be used:

In this model, the dependence of hydraulic conductivity, K (Ѱ), and the fluid content, θ, on the pressure head in Richard’s equation is specified with the following formulas, see Eqs. (12) and (13) (Baum et al. 2008).

θr = is the residual water content, θs is the water content at saturation.

Results

The GIS-TRIGRS model developed was tested and validated to predict snowmelt induced landslides in sensitive marine clays in the study region before the model was applied in the study area. Uncertainty in any given input parameter (e.g., data quality or spatial variability of the data), can significantly affect the accuracy of a landslide susceptibility map; therefore, it is necessary to consider the uncertainty for each parameter and its application to the study region’s model. The model initially was used to simulate landslide susceptibility in three scenarios:

(i)

A “normal” scenario: The average values of the study area’s geotechnical parameters were used as input data (e.g., cohesion, internal friction angle). Results from scenario one, are shown in Figs. 7a and 8a, including a depiction of the landslide susceptibility maps obtained in a normal scenario for a constant snowmelt duration of 10 days.

(ii)

The “worst case” scenario: The most negative, or pessimistic, values for each geotechnical parameter within the study area was used as input data (chosen with respect to their impact on the factor of safety, FS). For example, regarding the shear strength parameters, the lowest values for both cohesion and the internal friction angles were obtained from previous geotechnical studies and used as input data for the given parameters. This is a conservative approach, but using the landslide susceptibility map obtained from this scenario will greatly reduce the risk of slope failure. Figures 7b and 8b depict the landslide susceptibility maps obtained in a worst case scenario with constant snowmelt duration of 10 days.

(iii)

An “optimistic” or “ideal” scenario: Here, the most positive, or optimistic, values for each geotechnical parameter in the study area were used as input data (with respect to their impact on the FS). For example, for the shear strength parameter the highest values of cohesion and internal friction angles were obtained from previous geotechnical studies in the study area, and used as input data. Figures 7c and 8c depict the results for this scenario via landslide susceptibility maps obtained with constant snowmelt duration of 10 days.

As expected, a comparison of Fig. 7a, b, and c with Fig. 8a, b, and c respectively, shows that there are more vulnerable areas prone to landslides in the worst case scenarios (Figs. 7b and 8b), than in a normal (Figs. 7a and 8a), or an optimistic (Figs. 7c and 8c) scenario. Likewise, there are many more vulnerable areas prone to landslides in a normal analysis scenario than in an optimistic scenario. These findings therefore confirm that the landslide susceptibility model is sensitive to any changes in the input geotechnical parameter data.

In order to confirm the results of the GIS-TRIGRS model, a comparison analysis was completed by comparing the locations of the predicted unstable sites (areas with FS < 1) with the locations of historical landslides in the study area (see Fig. 9). Note that areas with a low factor of safety (less than one), represent locations with a high possibility of future landslides induced by snowmelt. It should be emphasized that the available map of historical landslides that occurred in Ottawa, include not only landslides triggered by snowmelt, but although by other factors (e.g., rainfall, earthquake, human activities). An analysis of Fig. 9 shows that many areas predicted to have a high susceptibility of snowmelt-induced landslides, follow similar trends as historical landslides in the Ottawa region. For example, the predicted landslides are approximately on the same path as compared to the historical landslides. Moreover, an analysis of Fig. 9 indicates that much of the predicted areas (over 92% of the predicted areas) with a FS < 1, which represent potential locations of future snowmelt-induced landslides, fall within the same general trends (often elongated trends) as the areas that were affected by historical landslides in the Ottawa region (i.e. the predicted landslides are close to or on the same path or trend of the previous landslides). The results show that a high correlation exits between the predicted and historical landslide spatial distributions in the region. This validation confirms that the developed GIS-based model can predict snowmelt-induced landslide susceptibility in the marine clays of Ottawa with relatively good accuracy. However, some historical landslides in Ottawa were not triggered by snowmelt. The model’s results remain consistent with this statement and confirm that (alongside previous geotechnical investigations on Ottawa landslides) besides snowmelt, rainfall, and earthquakes are other key triggers of landslides in sensitive marine clays (e.g., Quinn 2009; Aylsworth et al., 1997; Eden and Mitchell 1970).

Fig. 9

a Details of Zones A and B, location of historical landslides and map of predicted FS (normal case scenario) with respect geotechnical data. b Details of Zones A and B, location of historical landslides and map of predicted FS (worse case scenario) with respect geotechnical data

Discussion

The model developed was also applied to the study area to create a landslide danger map identifying the regions prone to snowmelt-induced landslides. The input data used for snowmelt (snowmelt duration and intensity) used in the GIS-TRIGRS based model is shown in Fig. 4. The effect of snowmelt intensity on slope stability was analyzed for the snowmelt durations of 6 to 48 h, 3 to 15 days, 25 days, and 30 days. For both the normal and worst case scenarios, landslide predictions were conducted with snowmelt intensities of 1.53 mm per hour (mm/hr) for 6 h (see Fig. 10a, b). Figures 11a, b, 12 and 13a, b depict the results for snowmelt intensities of 1.48 mm/hr., 1 mm/hr., and 0.54 mm/hr. for durations of 2 days, 10 days, and 30 days, respectively. As illustrated in both Figs. 10 and 13, the results confirm that the majority of the slopes in the study area are not susceptible to landslides (typically maintaining a FS ≥ 1). On the other hand, Figs. 10 and 13 demonstrate that many slopes in the predicted landslide susceptible areas are potentially unstable (FS < 1). These landslide susceptibility maps depict the areas which are most susceptible to landslides. In the Ottawa region, eastern and western parts of the city are more prone to landslide due to their sensitive clay and steep slopes of Leda clay. Furthermore, areas prone to landslides (typical FS < 1) are few, in the case of high snowmelt intensity and short snowmelt duration, compared to low snowmelt intensity and long snowmelt duration. These results would suggest that, within the Ottawa area, landslides are often triggered by low intensity, long duration snowmelt events compared to high intensity, short duration snowmelt events. This result aligns with the results from previous studies on snowmelt induced landslides, which have confirmed that long duration and low intensity snowmelt are major landslide triggers, such as the landslide that occurred in Shiidomari in September 1999 and the landslide in Katanooin January 2003 (Kawagoe et al. 2009). Moreover, numerous previous studies on landslides that took place in Canadian sensitive marine clays have concluded that a steady supply of surface waters from precipitation or snowmelt is a necessary precondition for landsliding. For example, a large sensitive marine landslide happened on 16 May 1971, along the South Nation River (SN) near Casselman, Ontario (Canada), within the Champlain Sea basin. The landslide affected approximately 28 ha of flat farmland, with a total volume of about 7 × 106 m3. The crater was approximately 640 m wide, and extended 490 m inland (Eden et al. 1971). On 20 June 1993, another large sensitive marine clays landslide took place along the South Nation River near Lemieux, Ontario (Canada), (LMX), only 4.5 km downstream from the 1971 South Nation. Approximately 2.8 × 106 m3 of soil was involved, with a crater area of 17 ha, extending 680 m back from the riverbank (Gauthier and Hutchinson 2012). Each of these two examples of sensitive marine clay landslides mentioned above occurred during the spring months (i.e. period of snowmelt). This consistent with other research findings that concluded that over 70% of all the large sensitive marine clay landslides in Quebec (Eastern Canada) have occurred during the springs (Lebuis et al. 1983). Furthermore, in each landslide case discussed above heavy precipitation and/or snowmelt have been identified as the trigger of the landslide (Gauthier and Hutchinson 2012).

Fig. 10

a FS map and areas prone to landslides for snowmelt duration of 6 h and snowmelt intensity of 1.53 mm/hrs based on (normal case scenario) of the geotechnical parameters. b F. S map and areas prone to landslides for snowmelt duration of 6 h and snowmelt intensity of 1.53 mm/hrs based on (worse case scenario) of the geotechnical parameters

a FS map and areas prone to landslides for snowmelt duration of 48 h and snowmelt intensity of 1.48 mm/hrs based on (normal case scenario) of the geotechnical parameters. b FS map and areas prone to landslides for snowmelt duration of 48 h and snowmelt intensity of 1.48 mm/hrs based on (worse case scenario) scenario of the geotechnical parameters

a FS map and areas prone to landslides for snowmelt duration of 10 days and snowmelt intensity of 1 mm/hrs based on (normal case scenario) of the geotechnical parameters. b FS map and areas prone to landslides for snowmelt duration of 10 days and snowmelt intensity of 1 mm/hrs based on (worse case scenario) of the geotechnical parameters

a FS map and areas prone to landslides for snowmelt duration of 30 days and snowmelt intensity of 0.54 mm/hrs based on (normal case scenario) of the geotechnical parameters. b FS map and areas prone to landslides for snowmelt duration of 30 days and snowmelt intensity of 0.54 mm/hrs based on (worse case scenario) of the geotechnical parameters

The characteristic of low intensity, long duration snowmelt that can lead to a landslide event, is due to a factor of safety that changes over time as a result of the water’s infiltration of the soil. The factor of safety keeps changing as the infiltration process continues to increase for up to 10 to 11 days of continuous snowmelt (see Figs. 10 and 11). The higher sensitivity of the safety factor to low intensity and long duration snowmelt is due to the fact that long snowmelt durations at a low intensity encourage the infiltration of water into the soil. This infiltration results in the reduction of shear resistance in the Leda clay and reduction in the matric suction, and thus reduces the safety factor and soil stability. This argument is consistent with previous geotechnical investigations that have shown that the upper parts of Leda clay formations remain unsaturated (Haché et al. 2015; Nader et al. 2015; Taha and Fall 2014; Quinn 2009; Fall et al. 2006; Dai and Lee 2001). Accordingly, this unsaturated state was also considered in the model developed. Moreover, this finding is agreement with the conclusions of the studies on snowmelt-induced landslides in Norway by Krogli et al. (2018). These authors concluded that intense or long-duration water supply, caused by rain and/or snowmelt, increases the water content in the soil or snow. The cohesiveness of soil or snow particles decreases with higher water content, increasing the landslide susceptibility. Similar conclusions were also made by Wayllace et al. (2019), who investigated the hydrological behavior of an infiltration-induced landslide in Colorado, USA. The authors stated that due to the infiltration rainwater and/or snowmelt water into hillslopes, there is variation of the water content in the hillslope and the water table level. Consequently, matric suction, suction stress, total unit weight, and effective stress change throughout the hillslope, and thus, the stability of the slope is affected.

Conclusion

A GIS-TRIGRS modelling approach was adopted in this study to assess and predict snowmelt-induced landslides in the Ottawa sensitive marine clays. Then, the model was successfully validated by comparing the predicted snowmelt induced landslide areas with the historical landslides maps of the study area. After validation, the model was applied to determine areas prone to snowmelt-induced landslides in the Ottawa region using geotechnical and hydrological information, and a DEM of the study area. Some new areas subjected to landslides were detected by comparing the model results to historically recorded landslides. Coupling GIS-TRIGRS together is a useful predictor for snowmelt-induced landslides. It enables users to create maps showing the distribution of the factor of safety. Thus, the developed approach could be considered a good method for assessing snowmelt induced landslide susceptibility in the Ottawa region. This approach can be improved, by considering more factors that may influence slope stability in different conditions, such as the effect of earthquake activity in seismically active regions. Further studies are recommended to be conducted to continue development or improvement of this tool.

Availability of data and materials

The datasets used and/or analysed during the current study are available from the corresponding author on reasonable request.

Quinn, P., 2009. Large landslides in sensitive clay in eastern Canada and the associated hazard and risk to linear infrastructure. A thesis submitted to the Department of Geological Sciences and Geological Engineering, Queen’s University Kingston, Ontario, Canada

Acknowledgments

The authors would like to thank the Natural Sciences and Engineering Research Council of Canada (NSERC), the University of Ottawa, and the Ministry of Higher Education and Scientific Research of Iraq for their financial support.

Funding

The authors gratefully acknowledge funding from Natural Sciences and Engineering Research Council of Canada (NSERC), University of Ottawa, and the Ministry of Higher Education and Scientific Research of Iraq.

Contributions

Each author has made substantial contribution to the conception and design of the work, the development of the simulation tool, the analysis and interpretation of the results or data as explained below. MA has made the most substantial contribution to the development of the simulation tool, the analysis and interpretation of the results or data. He has written the first draft of the manuscript. MF has made the most substantial contribution to the conception and design of the work. He also significantly contributed to the development of the simulation tool, the analysis and interpretation of the results or data. He was also scientific leader and supervisor of the present project. He also provided substantial inputs in the writing of the manuscript and reviewed the draft manuscript. BD has made substantial contribution to the conception and design of the work, the development of the GIS-TRIGRS simulation tool, the analysis and interpretation of the GIS and TRIGRS data and results. He also provided substantial inputs in the writing of the manuscript and review of the draft manuscript. All authors read and approved the final manuscript.

Corresponding author

Ethics declarations

Competing interests

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.