Assessment of landslide susceptibility for Meghalaya in North Eastern Region of India using bivariate and multi-criteria decision analysis models


 The state of Meghalaya of the North Eastern Region (NER) of India, situated in the India Himalayan Region (IHR), is the rainiest place in the country and falls under seismic zone V. The Himalayan ranges account for 80% of total landslide hazards in India. The main goal of the present study is to generate the GIS-based landslide susceptibility map (LSM) of Meghalaya by using frequency ratio (FR), Shannon entropy (SE), analytical hierarchy process (AHP), and fuzzy-AHP (FAHP) models and compare these models for the study area. Fifteen landslide conditioning factors are used for susceptibility mapping includes a slope, aspect, elevation, plan curvature, stream power index (SPI), topographic wetness index (TWI), land use land cover (LULC), normalized difference vegetation index (NDVI), distance from the river, road and faults, rainfall (30 years mean annual rainfall), soil texture, geomorphology, and lithology. Landslide inventory of 1330 landslide events is prepared and mapped from various sources. The inventory dataset is randomly split in a 70/30 ratio to make the training dataset (70%) used in the model and testing dataset (remaining 30%) for validation purposes. The southern escarpment, the southeast region of the study area, and hillslope along the roadside show high susceptibility for landslide occurrence in all four models. The LSMs produced in the present study are validated using the area under curve (AUC) value. The presented LSMs can help concerned authorities and planners to make sustainable development plans and formulate risk mitigation strategies keeping in mind the critical areas for landslide hazards.


Introduction
Landslide is a natural disaster, defined as the movement of a mass of rock, debris, or soil mass Himalayas, is highly prone to seismic hazards (seismic zone V), and experiences heavy rainfall. 39 The region has numerous faults, shear zones, and other tectonic features. Together rainfall, 40 high seismicity, and numerous tectonic features make the region highly susceptible to hazard 41 like a landslide. 42 To reduce the adverse impact of landslides, prepare risk mitigation strategies and plan the 43 infrastructural development accordingly, the landslide susceptibility studies are proven to be  Kanungo 2004). The conditioning factors are the factors associated with topography, 50 geomorphology, geology, land use land cover (LULC), anthropogenic activity, rainfall, 51 seismicity, etc. (Shano et al. 2020) and are responsible for the slope failure. The relation of 52 these factors with the past landslides forms the basis for estimating the future susceptibility of 53 landslide occurrence (Chimidi et al. 2017). 54 In recent times, with the use of GIS and remote sensing, several landslide susceptibility studies 55 have been carried out worldwide using various methods/models (Sarkar and Kanungo 2004;    Being in the IHR, it is one of the most tectonic-active regions and rainiest places globally 95 (Prokop 2014). The area received an average yearly rainfall of 1234.31 to 7467.48 mm between 6 1991 and 2020 (30-year period) (Fig. 1). The southern escarpment received the highest rainfall, 97 as high as 12000 mm annual rainfall (recorded in Cherrapunji). The elevation of the southern 98 escarpment of the study area is about 1200-1500 m and is related to the Dauki fault (along the 99 southern boundary), which is much steeper than the northern slope. Due to this sudden rise in 100 elevation over a short distance, the southern escarpment controls rainfall distribution over the 101 region. In the study area, the slope ranges from 0° to 76°.

102
The study area is covered by various lithologic formations, including Proterozoic Pleistocene) (Cn) types of formations (Fig. 2), the details of which are given in Table 1. The 107 region also consists of many lineaments and structural discontinuities and is associated with 108 active tectonics. With respect to land use land cover, most of the study area is covered by dense  Table 1). These topological, geological, and other geoenvironmental 112 factors make the study area more prone to disastrous events like landslides.     The plan curvature is derived from DEM using ArcMap 10.8 with a resolution of 10 m.

159
Curvature influences the surface erosion processes, especially during the rainfall, by either  Where ꞵ is the local slope ( in degrees). Where a is upslope catchment area, and tan(ꞵ) is the slope angle.

178
The present study prepared the SPI and TWI map using SAGA GIS tools in QGIS and classified      Table 2). The overall accuracy is 85.33%, while the kappa coefficient (k) 211 value is 0.824. The value of k > 0.8 shows that the used map is reasonably accurate.

283
These FR values of different classes (Table 5) are then used to obtain the prediction rate (PR) 284 of each factor which depicts the weightage of individual factors, using Equations 4-6.
Where RF is relative frequency, MAX(RFi,j) is the maximum value of RF of j th factor, 289 MIN(RFi,j) is the minimum value of RF of j th factor, PRj is the prediction rate of j th factor.

290
The PRj will be the weight of the j th factor, i. Entropy is the quantitative measurement of deviation, variability, instability, and uncertainty 296 of a system and can be used to predict the future trend of a specified system (Lotfi and   Table 5 shows entropy weights obtained for all the conditioning 310 factors. These are normalized and used to get the LSM shown in Fig. 6.  these values can vary inversely (1/9 to 1) when the element on the horizontal axis is more 324 dominant than that on a vertical axis (Table 3). In the present study, for assigning the degree Where CI = consistency index, λmax = principal Eigenvalue, and n = order of the matrix. And 333 RI = random consistency index that depends upon the order of the matrix (Table 4).

334
As per Saaty (2008), CR should be less than 0.10, only then the formed comparison matrix is 335 consistent, and if not so, it represents inconsistency in the factor ratings. One must revise the 336 matrix until it becomes consistent. In the present study, for the pairwise comparison matrix of 337 conditioning factors, the CR is equal to 0.049. Also, for the comparison matrix of classes of 338 each factor, the CR value is less than 0.10 ( Where Wj,AHP = weight of j th conditioning factors and wij,AHP = weight of an i th class of the j th 342 factor using AHP. Fig. 8 shows the LSM using this model.     Step 2: Calculation of fuzzy geometric mean value (ri) for i th criteria

378
In this step, the fuzzy weights are de-fuzzified using the center of area (COA) method Where wi is non-fuzzy weights.  The landslide susceptibility maps obtained using all the methods are classified into five 387 susceptibility classes (very low, low, moderate, high, and very high) based on the natural breaks 388 classification system (Pourghasemi et al. 2012b) (Fig. 5, 6, 7 & 8).  to each factor and their classes are calculated using FR and SE method, listed in Table 5. The

404
FR value shows a spatial correlation between factors and landslide inventory. Therefore, it is 405 assumed that the higher the FR, the larger the influence of a particular factor on the landslide.

406
In the present study, pixels with slopes equal to or greater than 30° have higher FR than others.
In AHP and FAHP models, the subcategory of 30°-40° and >40° slope also show more 408 significant influence than others ( as slope, rainfall, distance from road, lithology, and LULC are found with higher weight share 416 than others, while the distance from fault is found with the least weightage (Table 7). In the 417 FAHP model, the dominant landslide factors remain the same as AHP (Table 8). According to the FR model (Fig. 5), 2.17%, 5.98%, and 13.10% areas of 425 the total study region are classified as very high, high, and moderate susceptibility categories, 426 respectively (Fig. 7). For the SE model (Fig. 6), 2.07%, 5.38%, and 10.87% areas have very 427 high, high, and moderate susceptibility classes. Similarly, using the AHP model (Fig. 8) shows moderate susceptibility to landslide, the highest among all four models (Fig. 7). Along with the southern escarpment and southeast region of the study area, these classes are 433 concentrated along highways of the study area in the case of AHP and FAHP models (Fig. 8   434 and 9). The LSM produced using adopted models is validated using the receiver operating      prediction accuracy). The produced LSMs reveals that the southern escarpment of the study 516 area, the area in the southeast, and hillslopes along the roads possess great susceptibility for 517 future landslides. If the road network gets affected due to landslide events, the intra-518 district/state, inter-district/state connectivity get hampered and impart substantial economic 519 losses to the population in the region. Therefore, the presented LSM for the considered study 520 area can help the authorities and decision-makers to plan and manage the risk mitigation 521 strategies for future landslides and plan the sustainable infrastructure development in the region 522 accordingly.

524
The authors declare that they have no known competing financial interests or non-financial 525 interests or personal relationships that are directly or indirectly related to the work submitted 526 for publication that could have appeared to influence the work reported in this paper.