9. Calculation of Landslide Density#

In the present exercise, we will use hazard index method to calculate landslide susceptibility. Landslide hazard index method is based on the following formula (Original by: Cees J van Weseten).

In this analysis section, some python codes are also shown to calculated different parameters of analysis which can help to complete the tutorial fast.

(9.1)#\[W_{i} = ln \left (\frac{Denseclass}{Densemap}\right) = ln \left( \frac{\frac{Area(S_{i})}{Area(N_{i})}}{\frac{\sum Area(S_{i})}{\sum Area(N_{i})}} \right) \]

Where,

\(W_{i}\) =the weight given to a certain parameter class (e.g. a rock type, or a slope class).
\(Densclas\) = the landslide density within the parameter class.
\(Densmap\) = the landslide density within the entire map.
\(Area(S_{i})\) = area, which contain landslides, in a certain parameter class.
\(Area(N_{i})\) =total area in a certain parameter class.

9.1. Calculation of Weightage for Slope classes#

We need the rasterized landeslide map and reclassified map of slope. Calculation of weightage is done in following steps:

  • Calculate area for each slope class: Use this code in the console

r.stats -a --overwrite input=recl.slope,lsrast output=path/to/slpcross.csv separator=, null_value=-1

  • Caluclate weightage for each class

Hide code cell source
import pandas as pd
df = pd.read_csv('../data/classstat/slpcross.csv', header=None)
pd.options.display.precision = 5
df.columns = ['Slope Class', 	'Landslide', 'Area(m2)']
df = df[df['Slope Class'] > 0]
df.style.hide(axis="index")
Slope Class Landslide Area(m2)
1 -1 11055468.750000
2 1 15156.250000
2 -1 50358593.750000
3 1 42812.500000
3 -1 43914531.250000
4 1 19687.500000
4 -1 4630781.250000
5 1 7187.500000
5 -1 93437.500000

Above table gives the density of landslide within the class as:

  • Column-1(Slope class): value=1 represents slope interval 0 to 15

  • Column-2(Landslide): 1 means landslide area and -1 means no landslide

  • Column-3(Area): Area in \(m^2\) corresponding to the slope class and existense (True or False) of landslide Total area of each class can be calculated by summing areas corresponding to existing and non-existing landslide for each class. For example (Slope class-2),

##-------------
# Calculation of Area for slopeclass-2 and other parameters
##-------------
c2l = df.loc[(df['Slope Class'] == 2) & (df['Landslide'] >0)] #Area of landslide
print(f'Area of landslide = {c2l.iloc[0,2]} m2')
c2n = df.loc[(df['Slope Class'] == 2) & (df['Landslide'] <0)] #Area of no landslide
print(f'Area of no landslide = {c2n.iloc[0,2]} m2')
t2 = c2n.iloc[0,2] + c2l.iloc[0,2] # Total Area
print(f'Total area for the class = {t2} m2')
sum_l =  df.groupby(["Landslide"]).sum() # Sum of area of landslide/no landslide
sum_l.drop(columns=["Slope Class"], inplace = True)
sum_s_l = df.groupby(["Slope Class", "Landslide"]).sum() # Sum of area of landslide/no landslide by class
sum_t = df.sum()
sum_t.drop(columns=["Slope Class", "Landslide"], inplace = True)

print('### Sum of area of Landslide/No Lanslide by slope class')
print(sum_s_l)

print('### Sum of area by landslide/no landslide')
print(sum_l)
print(f"Total area of raster map = {sum_t['Area(m2)']}")
Area of landslide = 15156.25 m2
Area of no landslide = 50358593.75 m2
Total area for the class = 50373750.0 m2
### Sum of area of Landslide/No Lanslide by slope class
                          Area(m2)
Slope Class Landslide             
1           -1         1.10555e+07
2           -1         5.03586e+07
             1         1.51562e+04
3           -1         4.39145e+07
             1         4.28125e+04
4           -1         4.63078e+06
             1         1.96875e+04
5           -1         9.34375e+04
             1         7.18750e+03
### Sum of area by landslide/no landslide
              Area(m2)
Landslide             
-1         1.10053e+08
 1         8.48438e+04
Total area of raster map = 110137656.25

From above result of calculation, we get the followings for slope-class=2

Area of landslide, \(Area(S_{i}) = 15156.25\)

Area which does not have any landslide = 50358593.75

Total Area within the class, \(Area(N_{i}) = 50373750.0\)

Total Area of landslide, \(\sum Area(S_{i}) = 84843.75\)

Total area of map = 110137656

Weightage can be calculated using formula (for i=2)

Following code calculates the area of landslide/no landslide for each slope class (in wide format) and the value of \(Denseclas\).

Note

The area of landslide is changed to 0.0001 where the area is zero to avoid the error when calculating logarithm.

##-------------
# Calculation of Weightage for slopeclass-2 and other parameters
##-------------
import math
densemap = sum_l['Area(m2)'][1]/sum_l['Area(m2)'][-1]

landslide_s =  pd.pivot(df, index=['Slope Class'], columns = 'Landslide', values ='Area(m2)' )
landslide_s[1] = landslide_s[1].fillna(0.0001)
landslide_s['denseclas'] = landslide_s[1]/(landslide_s[-1] +  landslide_s[1])
landslide_s['densemap'] = densemap
landslide_s['weight'] =landslide_s['denseclas']/landslide_s['densemap']
landslide_s['weight'] = landslide_s['weight'].apply(lambda x: math.log(float(x)))
landslide_s.rename(columns={1: "Landslide Area(m2)"}, inplace=True)
landslide_s.rename(columns={-1: "Non-Landslide Area(m2)"}, inplace=True)
slopecl = landslide_s
slopecl.to_csv('../data/classstat/slopewt.csv')
##---Use print instead of display for your code----
display(landslide_s)
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
Slope Class
1 1.10555e+07 0.0001 9.04530e-12 0.00077 -18.26087
2 5.03586e+07 15156.2500 3.00876e-04 0.00077 -0.94091
3 4.39145e+07 42812.5000 9.73956e-04 0.00077 0.23376
4 4.63078e+06 19687.5000 4.23344e-03 0.00077 1.70316
5 9.34375e+04 7187.5000 7.14286e-02 0.00077 4.52885

Following is a sample calculation of class weightage for slope class = 2

\(Denseclas = 15156/50355000.00 = 3.008760e-04 \)
\(Densemap = 84843.75/110137656.25 = 0.00071\)
\(W_{i} = ln\left(\frac{Denseclas}{Densemap} \right) = ln\left(\frac{3.008760e-04}{0.00071} \right) = -0.940908\)

Similarly, weightage for other parameters are shown in following sections

9.2. Weightage for Aspect Directions (Classes)#

##-------------
# Calculation of Weightage for Aspect classes
##-------------
### Create function to obtain the weigthage

def createWeight(fn, cln):
    df = pd.read_csv(fn, header=None)
    df.columns = [cln, 	'Landslide', 'Area(m2)']
    df = df[df[cln] > 0]

    sum_l =  df.groupby(["Landslide"]).sum() # Sum of area of landslide/no landslide
    sum_l.drop(columns=[cln], inplace = True)
    sum_s_l = df.groupby([cln, "Landslide"]).sum() # Sum of area of landslide/no landslide by class
    sum_t = df.sum()
    sum_t.drop(columns=[cln, "Landslide"], inplace = True)

    densemap = sum_l['Area(m2)'][1]/sum_l['Area(m2)'][-1]

    landslide_s =  pd.pivot(df, index=[cln], columns = 'Landslide', values ='Area(m2)' )
    landslide_s[1] = landslide_s[1].fillna(0.0001)
    landslide_s['denseclas'] = landslide_s[1]/landslide_s[-1]
    landslide_s['densemap'] = densemap
    landslide_s['weight'] =landslide_s['denseclas']/landslide_s['densemap']
    landslide_s['weight'] = landslide_s['weight'].apply(lambda x: math.log(float(x)))
    landslide_s.rename(columns={1: "Landslide Area(m2)"}, inplace=True)
    landslide_s.rename(columns={-1: "Non-Landslide Area(m2)"}, inplace=True)
    return landslide_s

filename = '../data/classstat/aspcross.csv'
classname = 'Aspect Class'
aspcl = createWeight(filename, classname)
aspcl.to_csv('../data/classstat/aspwt.csv')
display(aspcl)
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
Aspect Class
1 1.07527e+07 5625.0000 5.23127e-04 0.00077 -0.38810
2 1.11962e+07 781.2500 6.97778e-05 0.00077 -2.40261
3 1.36478e+07 0.0001 7.32718e-12 0.00077 -18.47184
4 1.26661e+07 0.0001 7.89509e-12 0.00077 -18.39719
5 1.36111e+07 3593.7500 2.64031e-04 0.00077 -1.07186
6 1.73469e+07 33906.2500 1.95460e-03 0.00077 0.93002
7 1.53427e+07 16406.2500 1.06932e-03 0.00077 0.32686
8 1.54547e+07 24531.2500 1.58730e-03 0.00077 0.72187

9.3. Weightage Distribution for NDVI Classes#

##-------------
# Calculation of Weightage for NDVI classes
##-------------

filename = '../data/classstat/ndvicross.csv'
classname = 'NDVI Class'
ndvicl = createWeight(filename, classname)
display(ndvicl)
ndvicl.to_csv('../data/classstat/ndviwt.csv')
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
NDVI Class
2 1.80312e+05 781.2500 4.33276e-03 0.00076 1.73448
3 5.04062e+05 3437.5000 6.81959e-03 0.00076 2.18808
4 1.27016e+06 6562.5000 5.16669e-03 0.00076 1.91051
5 4.97375e+06 12812.5000 2.57602e-03 0.00076 1.21453
6 2.76603e+07 23281.2500 8.41684e-04 0.00076 0.09593
7 6.49323e+07 32656.2500 5.02927e-04 0.00076 -0.41903
8 1.12994e+07 5312.5000 4.70159e-04 0.00076 -0.48641
9 1.30781e+05 0.0001 7.64636e-10 0.00076 -13.81559

9.4. Weightage Distribution for NDWI Classes#

##-------------
# Calculation of Weightage for NDWI classes
##-------------

filename = '../data/classstat/ndwicross.csv'
classname = 'NDWI Class'
ndwicl = createWeight(filename, classname)
display(ndwicl)
ndwicl.to_csv('../data/classstat/ndwiwt.csv')
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
NDWI Class
2 1.14062e+04 0.0001 8.76712e-09 0.00076 -11.37622
3 1.48875e+06 11250.0000 7.55668e-03 0.00076 2.29071
4 5.24361e+07 54375.0000 1.03698e-03 0.00076 0.30459
5 5.69534e+07 19218.7500 3.37447e-04 0.00076 -0.81807
6 6.14062e+04 0.0001 1.62850e-09 0.00076 -13.05957

9.5. Weightage Distribution for Profile Curvature Classes#

##-------------
# Calculation of Weightage for Profile Curvature classes
##-------------

filename = '../data/classstat/profcross.csv'
classname = 'Profile Curvature Class'
profcl = createWeight(filename, classname)
display(profcl)
profcl.to_csv('../data/classstat/profwt.csv')
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
Profile Curvature Class
1 5.62500e+03 0.0001 1.77778e-08 0.00077 -10.67741
2 5.26891e+06 3125.0000 5.93102e-04 0.00077 -0.26224
3 4.95120e+07 37968.7500 7.66859e-04 0.00077 -0.00530
4 5.08369e+07 41875.0000 8.23713e-04 0.00077 0.06621
5 4.42312e+06 1875.0000 4.23908e-04 0.00077 -0.59809
6 6.09375e+03 0.0001 1.64103e-08 0.00077 -10.75746

9.6. Weightage Distribution for Tangent Curvature Classes#

##-------------
# Calculation of Weightage for Tangent Curvature classes
##-------------

filename = '../data/classstat/tancross.csv'
classname = 'Tangent Curvature Class'
tancl = createWeight(filename, classname)
display(tancl)
tancl.to_csv('../data/classstat/tanwt.csv')
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
Tangent Curvature Class
1 1.70312e+04 0.0001 5.87156e-09 0.00077 -11.78524
2 8.74219e+06 12343.7500 1.41197e-03 0.00077 0.60514
3 4.47438e+07 38593.7500 8.62551e-04 0.00077 0.11229
4 4.76452e+07 27500.0000 5.77184e-04 0.00077 -0.28945
5 8.89078e+06 6250.0000 7.02975e-04 0.00077 -0.09228
6 1.39062e+04 156.2500 1.12360e-02 0.00077 2.67927

9.7. Weightage Distribution for River Distance Classes#

##-------------
# Calculation of Weightage for River Distance classes
##-------------

filename = '../data/classstat/rivercross.csv'
classname = 'River Distance Class'
rivercl = createWeight(filename, classname)
display(rivercl)
rivercl.to_csv('../data/classstat/riverwt.csv')
Landslide Non-Landslide Area(m2) Landslide Area(m2) denseclas densemap weight
River Distance Class
1 1.51038e+07 14843.75 0.00098 0.00076 0.25091
2 1.43508e+07 15000.00 0.00105 0.00076 0.31252
3 3.84783e+07 12031.25 0.00031 0.00076 -0.89431
4 3.93889e+07 31406.25 0.00080 0.00076 0.04180
5 3.62938e+06 11562.50 0.00319 0.00076 1.42698