Evaluation of Landslide Susceptibility in Cau River Basin Using a Physical-Based Model under Impact of Climate Change

This paper evaluated the probability of landslide susceptibilities through the application of the Transient Rainfall Infiltration and Grid-Based Region Slope-Stability model in Cau river basin (Vietnam) using the scenarios-based approach under the influence of the warming climate. The tested cases were developed based on various options including rainfall amount and distribution, soil depth determination, and land-cover conditions. Input data for extreme rain events included historical rainstorm in 2013, the Probable Maximum Precipitation (PMP) with the durations of 24 hours and 48 hours. The results illustrated the reduction of slope stability when the land cover changed from land-use data in 2007 (Ha12) to land-use data in 2015 (Ha22). When the whole region was assumed to be replaced by soil (Ha02), the factor of safety (Fs) decreased to lower magnitude when compared to Fs value regarding to changes in land cover condition (Ha12 & Ha22) and changes in soil-depth (Ha33). The model simulations demonstrated the agreement with the slope-failure hazard association with the destabilizing factor such as slope-cutting activities at historical landslide events. Under the same land-cover and soil depth condition, the average value of factor of safety regarding to the historical rainstorm in 2013 (Ha32) declined by 0.069 and 0.189 when compared to Fs of the 24-hour PMP with the storm distribution type 3 (1332) and Fs of the 48-hour PMP with the storm distribution type 3 (2332), respectively. The results reveal that in a warming climate, changes in extreme precipitation in terms of rain-total, rain-duration, and rain-distribution would result in the expansion of slope instability in the hilly region. This application is considered as a prevailing method for landslide susceptibility analysis and would provide important information for authorities in developing adequate land-management in the river basin. How to cite this paper: Le, T.T.T. and Kawagoe, S. (2019) Evaluation of Landslide Susceptibility in Cau River Basin Using a Physical-Based Model under Impact of Climate Change. Open Journal of Modern Hydrology, 9, 1-19. https://doi.org/10.4236/ojmh.2019.91001 Received: October 29, 2018 Accepted: November 30, 2018 Published: December 3, 2018 Copyright © 2019 by authors and Scientific Research Publishing Inc. This work is licensed under the Creative Commons Attribution International License (CC BY 4.0). http://creativecommons.org/licenses/by/4.0/ Open Access T. T. T. Le, S. Kawagoe DOI: 10.4236/ojmh.2019.91001 2 Open Journal of Modern Hydrology


Introduction
Vietnam, a tropical country locating within the Southeast Asian typhoon belt, often suffered from tropical cyclones and depressions annually.Recently, abnormal weather phenomena such as torrential rain have been occurring with more frequencies and higher intensities; water-related disasters including floods, debris flows and landslides have also borne catastrophic consequences to the society.Being parallel with the requirements of economic development and population pressure in Northeast region, a large number of arteries connecting residential zones have been constructing and rehabilitating in such hilly region.Exploiting activities from mining, deforestation and infrastructure constructions in mountainous areas have led to the increase of bared land on hill-sides.High rainfall intensities in addition to human interferences in the natural slope condition (slope-cut and tree-cutting on hilly sites), and changes in land cover condition are major triggers to the landfall occurrences in this area.The expansion of the development activities could lead to the rise in slope failures especially during triggering events such as prolonged extreme precipitation, which has been influenced greatly by climate change.
According to [1] [2] and [3], there are two major techniques for landslide susceptibility; they are qualitative method and quantitative method.The qualitative method, including inventory and heuristic approaches, refers to expert-based method; it based on landslide inventory and historical information to evaluate and determine main parameters among them.It also identifies sites with similar conditions of geology and geomorphology to evaluate the susceptibility of failure.The inventory map is implemented using satellite images, ground survey and historical database of landslide occurrence.This method highlights the location and scales of historical landslides, and is seen as one of the most susceptibility mapping techniques.It would give spatial distribution of past landslides and might be utilized for evaluating reducing risks of slope failures on a regional scale; however it does not provide prediction of landslide susceptible areas.The heuristic, which requires long-term landside data and causative parameters, estimates potential landslide from preparatory variables.Such an expert-driven approach depends on the knowledge and experiences of scientists in deciding degree and type of landslide risks.The quantitative methods, including statistical and deterministic methods, utilized mathematical model to assess the probability of slope-failure occurrence so as to identify risk regions.The statistical/stochastic method analyses the relationship between existing landslides and instability factors; it would provide prediction of landslide in region that might not have landslide in the past based on statistical analysis of landslide ground characteristics.However it Open Journal of Modern Hydrology needs a large amount of data collection and is suitable for medium scale prediction.The deterministic approach requires physical processes to analyse the equilibrium mechanical of potential block, and to calculate a factor of safety [4].It can be applied for detailed studies at large scale [5].However, difficulties in obtaining, checking and processing data especially large spatial dataset in regional slope stability studies reveal the limitation of this approach.The application of deterministic is often in case of simple landslide types in fairy homogeneous geology condition [2].
Recent studies analysed the slope instabilities in North area using the statistical method for a regional scale; the deterministic method is often applied for a local scale.In addition, the Vietnamese Government conducted a national project 2012-2020 to produce the landslide inventory maps for mountainous provinces.In fact, it is difficult to evaluate and obtain records of slope failures locating in remote areas; this project thus implemented investigations along main arteries to map landfall records that are related to the man-made slope failures in populated regions.It is noted that the Northern region is relied strongly on agriculture activities [6]; it thus is more prone to changes in weather when it has experienced with the growing number of extreme weather events during the past few years.This paper therefore would like to address the probability of landslide susceptibilities through the development of physical model in a regional scale focusing the roles of rainfall in relation to conditioning factors including changes land-use under the influence of the warming weather.

Study Area and Data
The Cau river basin belonging to Bac Kan and Thai Nguyen provinces is the fo-   was released in 2011 by the Ministry of Economy, Trade and Industry, Japan (METI) and the National Aeronautics and Space Administration (NASA) [9] with the spatial resolution of 30 × 30 m (Figure 3(a)).The slope map was derived from DEM map as described in Figure 3(b).The land use map in North Vietnam was adopted from the Japan Aerospace Exploration Agency (JAXA) [10].This high resolution land cover map was produced in 2007 and 2015 with nine categories: 1) Water; 2) Urban and built-up; 3) Rice; 4) Other crops; 5) Grasslands; 6) Orchards; 7) Barren; 8) Forest; and 9) Mangrove (Figure 4).

Methodology
In this paper, the landslide susceptibility is analysed based on the deterministic method using a Transient Rainfall Infiltration and Grid-Based Region Slope-Stability (TRIGRS) model [11].In this program, the pore-water pressure and the value of factor of safety (Fs) are calculated on a cell-by-cell basis and can be operated in the geographic information system [11].The model evaluates changes in the transient pore pressure and the factor of safety due to rainfall infiltration; it is thus often applied for analysing rain-induced landslides.The infiltration model for initially unsaturated conditions of finite depth of the bedrock was selected for this river basin; as described Figure 5(a) the soil profile was divided into two main layers including a saturated zone beneath the water table, and an unsaturated zone expanding to the ground surface [12].The factor of safety Fs at depth Z below the ground surface was derived from the balance between the downslope components of the gravitational driving stress and the resisting stress for saturated soils; Fs was defined by separating the time-variant term and the steady terms as in Equation ( 1) [11]: Open Journal of Modern Hydrology where φ is the soil friction angle (˚); δ is the slope angle (˚); c is the cohesion (kN/m 2 ); γ s is the unit weight of soil (kN/m 3 ); γ w is the unit weight of water (kN/m 3 ); Z is the depth below the ground surface (m); t is the time (sec); ψ refers to the ground water pressure head as a function of depth Z and time t (m).
In this study, a scenario-based approach is applied for landslide risk assessment under extreme weather condition; main input layers for TRIGRS model include rainfall intensities and durations in relationship with soil physics, hydro-geology and slope conditions (Figure 5(b)).In addition to the Digital Elevation Model (DEM), information of soil and precipitation was discretized and assigned in each cell in order to simulate rain infiltration as well as vertical flow through the unsaturated zone.Regarding to the thickness of topsoil-a major parameter for the slope instability evaluation, it is difficult to obtain this information and it is not feasible for field investigation upon a large scale region.The spatial distribution of this parameter is usually predicted based on empirical relationship between topography attributes and soil thickness.In this paper, we applied the Saulnier method (1997) [13] to determine soil thickness in each pixel; this approach is based on the assumption of decreasing linear function between topographic slope and soil depth as described in Equation ( 2): In which Z min and Z max refer to minimum and maximum values of effective soil depths, respectively; δ min and δ max are the minimal and maximal values of slope angle, respectively; Z i and δ i represent the soil depth and slope angle at the calculated cell i.
Moreover, in order to evaluate the role of vegetation cover on the hilly slope, Open Journal of Modern Hydrology in this article, we considered the contribution of root systems to the shear strength (s) by additional cohesion component (Δc) following the Mohr-Coulomb formulation [14] (Equation (3)): ( ) where the cohesion c in Equation ( 1) in vegetation cover area is the sum of soil cohesion c s (kN/m 2 ) and apparent cohesion provided by roots (root cohesion) Δc (kN/m 2 ); δ N refers to the normal stress on the shear plane (kN/m 2 ).For simplicity, in this study, we ignored the influence of tree-surcharge; the mechanical impact of tree cover on the shear strength of soil with vegetation cover is therefore represented only by the root cohesion (Δc).
In order to simulate shallow landslide under different extreme precipitation, it is necessary to evaluate spatial and temporal distribution of storm events.As described in Figure 1, there are 8 meteorological stations at 3 provinces in the Northeast region covering the area of the Cau basin including Bac Kan (BK), Cho Ra (CR), and Ngan Son (NS) stations in Bac Kan province, DinhHoa (DH) and Thai Nguyen (TN) stations in Thai Nguyen province, as well as ChiemHoa (CH), Ham Yen (HY), and TuyenQuang (TQ) stations in TuyenQuang province.In this paper, we analysed historical records from 1960-2016 at the 8 aforementioned meteorological stations in addition to 5 rain-gauge stations in Bac Kanprovince (Cho Don (CD) and Phu Thong (PT)) and Thai Nguyen province (Dai Tu (DT), Pho Yen (PY), and Vo Nhai (VN)) covering the entire area of this catchment.Extreme precipitation was selected from storm events that resulted in landslides.Additionally, under the impacts of climate change, the rises in temperature are believed to drive stronger downpours in many regions across the world.The magnitudes of extreme rainstorm under warmer climate were computed using the Probable maximum precipitation (PMP) relation to the moisture maximization method [15] [16].The temporal rainfall distributions were also analysed regarding to the most heavy rainstorm events in the study site.

Scenarios Setup
To evaluate the slope stability condition in the Cau river basin, we applied the scenario-based approach, in which we considered options of rainfall intensity and distribution, characteristics of surface condition, and the soil depth estimation.The actual rainstorm event was selected based on rain records from 2010 to 2016.During this period, the annual rainfall in the study area varied from 1391 mm (2011) to 2029 mm (2013); the highest values of the maximum-24-hour precipitation at most rain-stations occurred in the year 2013.Heavy rain in 2013 also resulted in a great number of landslide occurrences across the region; for example, according to the report of Bac Kan province, the road 258 from Bac Kan to Ba Be with the length of 40 km was closed in 1 month and was rehabilitated in about 3 months as a consequence of landslide occurrences triggering af-ter the downpour in the late of May 2013.In addition, the national project also conducted the landslide inventory map in Bac Kan province in 2013.Despite the fact that the investigation was implemented through the transect walks method across traffic routes in this province, and the landfall often associated with slope-cut activities; this database of landslides still provided important references for the study of sliding issues in the study area.Therefore we utilized heavy rainstorm in 2013 (Figure 6(a)) to analyse the probability of landfall in this catchment under the impact of the historical rain event.The Thiessen polygons were constructed around stations within and near the basin to determine the areas that were represented by the stations as illustrated in Figure 1(a) [17]; the Inverse Distance Weight (IDW) interpolation was also applied to estimate the spatial distribution of rainfall in the catchment as described in Figure 6(a) [16] [17] [18].
Analysing 68 heavy storm events from 1960 to 2016 in the Northeast provinces, we found that the heavy storm events often lasted in less than 48 hours; the rain durations fluctuated from 5 hours (TuyenQuang station on 1984/6/22) to 53 hours (Bac Kan station on 1990/9/21) as illustrated in Figure 7(a).In this article, proposed extreme rainfall events were hence selected based on the probable maximum precipitation regarding to the duration of 24 hours (Figure 6(b)) and 48 hours [16].The analysis of hourly rainfall regard to 68 heavy rainstorm events illustrated that the hourly distributions of precipitation varied significantly and could be grouped into three main distribution types: 1) Type 1: peak rainfall locates at the beginning of the storm; 2) Type 2: peak rain locates at the middle of the storm time; and 3) Type 3: peak precipitation appears at the end of the storm (Figure 7(b)).The values of heavy rainstorm in 2013, as well as 24-hour PMP and 48-hour PMP at 8 stations in the study site are described in Table 1.
Information of soil and hydro-geology in the catchment was gathered from recent studies [19]- [25] as well as the field investigations in 2017 and 2018, in which soil samples were tested in the geotechnics laboratory in Thuy Loi university.The average input parameters of soil properties regarding to main geology regions are summarized in Table 2.
In addition to the worst case of bared land across the region (the cohesion c in Equation ( 1) refers to the soil cohesion c s ), the land-use data obtaining from JAXAin two time periods 2007 and 2015 was put in the model to set up scenarios; the cohesion c in Equation ( 1) is the sum of the soil cohesion c s and the root cohesion Δc.The root cohesion was estimated using the trial-and-error method [26] [27] [28]; the average values varied from 0.7 KN/m 2 in orchards, savannas and shrub lands to 1.0 KN/m 2 in forest cover.
Applying the Saulnier method [13] (Equation ( 2)), information of average depth of top soil was gathered from recent studies from the government project of landslide inventory map, and was set up in three options including 1) case 1: upper value of soil depth (0.5 m -4 m) based on lithology map across the whole region; 2) case 2: average value of soil depth (0.5 m -3 m) based on lithology  map and geology distribution across the whole region; 3) case 3: average value (0.5 m -3 m) of soil depth based on geology distribution and rainfall regions (Table 3 & Figure 8).Because the information of groundwater table was unavailable, it was assumed that the initial groundwater table was at the bottom of the weather soil layer in summer time without antecedent rain.Based on options of rainstorm event amount (R), rainfall distribution (D), land cover condition (L), and soil depth assumption (S), we set up scenarios (RDLS) to compute the factor of safety for the whole region of Cau river basin as described in Table 4.

Evaluation of Landslide Susceptibility
In case of actual rainstorm distribution in 2013, we calculated the factor of safety about 60 percent of the sliding area had the factor of safety that was less than 2.0.
The case 1 of soil-depth (S = 1) still obtained the highest proportion of area having low values of Fs; while the option using the land-cover data from 2007 (Ha12) indicated only 51% of the area that had values of Fs below 2.0.In addition, throughout two field trips in September 2017 and September 2018, we investigated 56 sliding points in the catchment.Since the shape and area of each landfall was not mapped, we tried to evaluate the responses of slope stability regarding to scenario-based simulations within a region of 30 meter × 30 meter at each slide point.About 80 percent of the total area obtained Fs from 0.82 to 2.
Table 3. Summary of options in soil depth range in the study area. of cells that could be prone to the instability was higher than that of option Ha32 (Figure 12).
In addition, in order to estimate the changes in Fs, we considered the option Ha32 as a base case and computed the difference between Fs_Ha32 and factor of safety for other scenarios throughout the basin (Figure 13).In general, the values of average ΔFs in the whole basin were 0.

Conclusions
According to [29], the stability condition was classified into four groups: 1) Unstable when Fs < 1, the stability requires stabilizing factors; 2) Quasi-stable when 1 ≤ Fs < 1.25, the instability would be resulted from minor destabilizing factors; 3) Moderately stable when 1.25 ≤ Fs < 1.5, the instability would be resulted from In fact, most of historical landslides were investigated along main arteries and often associated with road-cut slopes; the instability condition of such man-made slopes were more prone to changes in the environmental condition  In addition, under the changes in climatic condition, extreme precipitation would result in higher probability of unstable regions.It is noted that under the same amount of the storm-rain and storm duration, for example 24-hour-PMP, the landslide prone areas would expanse in case of changing the storm distribution from type 1 (scenario 1132) to type 3 (scenario 1332).Comparing the rise of ΔFs in case of changing rain-event from actual storm in 2013 to 24-hour-PMP (126.9 mm -332.4 mm rise in rain amount), the factor of safety still decreased significantly when the downpour total of 24-hour-PMP increased by 21.7 mm -169.1 mm in a further 24 hours to reach the magnitude of 48-hour-PMP (Figure 13).This reveals the importance of three major components of a storm events Open Journal of Modern Hydrology including rain-amount, rain-duration, and rain-distribution to the stability condition of hilly sites.
In conclusion, the model simulations obtained uncertainties raising from the estimations of soil-depth and groundwater based on geologic and hydrologic background for a regional scale, as well as the application of average values of soil properties and the same thresholds of root cohesion Δc.Subjective errors in drawing polygons from the Google Earth soft-ware, especially in case of small landslides also hindered the evaluation of the historical landfalls.However, the findings would provide important information for authorities in developing adequate land-management in the river basin considering the relationship between landslide-prone area and changes in land-cover situations.
cus of this study.The catchment locates in the Northeast region, lying from 105.48˚ to 106.13˚ East longitude, and from 21.35˚ to 22.32˚ North latitude.According to the historical record of daily rainfall during 1962-2016 which was collected from the Vietnam Meteorological and Hydrological Administration, the average annual precipitation ranges from 1360 mm (Cho Ra) to 2572 mm (Diem Mac); the highest 24-hour rainfall varies from 208 to 475 mm (Figure 1(b)).Regarding to ecological stratification, there are two major sub forest-ecological regions: Low mountains (TV10) and midland (TV12) [7] [8].Major forest types are mixed closed evergreen humid forests (TV12), and mixed closed forests on limestone mountains/valleys (TV10, TV12).According to the national statistic data, within the period 2011-2016, this region illustrated significant rises of 1.24 -1.54 percent per year in the population size.The number of residents in Bac Kan province grew from 300.4 thousand in 2011 to 319 thousand with the density of 65.64 people per square kilometre in 2016.With total people of 1227.4 thousand in 2016, Thai Nguyen province experienced a rapid growth rate of 88 thousand people when compared to population size in 2011; it also brought about the population density of 348.03 people per square kilometre in 2016.Such significant rise in population has burdened high pressure on Open Journal of Modern Hydrology

Figure 1 .
Figure 1.(a) The meteorological and rain-gauge stations in the Thiessen network within three provinces Tuyen Quang, Bac Kan & Thai Nguyen; (b) The highest 24-hour precipitation (Rmax) in the Cau basin acorrding to the statistical records from 1962-2016.

Fs regarding to three
options of soil-depth as well as the changes in land-use option i.e. from land use data in 2007 (JAXA) (Ha12) to land-use data in 2015 (JAXA) (Ha21, Ha22, Ha31, Ha32 & Ha33) in addition to the worst option of land cover with bare soil (Ha02).Illustrations of factor of safety (Fs) relation to the historical rainfall event are described in Figure 9; Figure 10 shows the determination of Fs regarding to changes in the extreme-precipitation event including 24-hour PMP (scenarios 1132 and 1332), and 48-hour PMP (scenarios 2332 and 2333).With regard to historical rainfall event, case 1 in soil-depth determination (Ha21 and Ha31) provided lower Fs values than the results from case 2 and case 3.In case of bared soil across the whole region (Ha02), more unstable areas were seen in northern and south-eastern parts of the catchment when compared to results from the three options Ha12, Ha22, and Ha32.With regard to historical landslide data in 2013 from the Government project, we obtained shapes of 44 landslide points using the application of Google Earth.The summaries of the Fs values within those polygons were illustrated in Figure 11;

*.
For example, scenario Ha32 refers to the option applying the historical rainfall event (R = H) having an actual distribution (D = a) using JAXA land-use in 2015 with the resolution of 50 meter (L = 3), and the soil depth case 2 (S = 2).−0.046 (Ha32_Ha22).Within the range of Fs_Ha32 from 0 to 1.5, the option Ha12 and Ha22 obtained smaller threshold of ΔFs (negative values); using the same resolution of 100 meter in the land-use data, changes in the land cover such as from forest cover (Ha12 in 2007) to crop land and construction sites (Ha22 in 2015) revealed the reduction of the slope stability regarding to the range of Fs below 1.0.Interestingly when the whole region is assumed to be replaced by soil (Ha02), a significant reduction of Fs was shown when compared to the option of changes in land-use database (Ha12 and Ha22) and changes in soil-depth option (Ha33).The worst cases in assessing the slope stability still came from the case 1 of soil-depth option (Ha21 and Ha31).The differences between two aforementioned lines within the range of Fs_Ha32 from 1 -1.5 illustrated the great impact of grid-cell resolution to the model simulation.When the rainstorm amount changed from actual event (Ha32) to 24-hour-PMP (1132 and 1332) and 48-hour-PMP (2332 and 2333), the factor of safety Fs also declined in the whole basin.

Figure 11 .Figure 12 .
Figure 11.(a) Distribution of the factor of safety (Fs) in accordance with percent of accumulation area of the landslide polygons in 2013; (b) Detail of Figure 11(a) regarding to the Fs below 2.0.

Figure 13 .
Figure 13.(a) Relationship between factor of safety for option Ha32 (Fs_Ha32) and Delta Fs (ΔFs) [difference between Fs_Ha32 and factor of safety for other options (Fs_RDLS)]; (b) Detail of Figure 13(a) regarding to the Fs_Ha32 below 1.5.

Table 1 .
Total rain amount regarding to the historical rainfall event as well as the 24-hour PMP and the 48-hour PMP (mm).

Table 2 .
Average values of soil properties in Cau basin.