H
Hicham AMAR
Guest
Time to code on R!
The theoretical part for the sensitivity analysis applied on aquifer vulnerability to pollution is described in the part I.
This article will explain the step to evaluate the sensitivity analysis applied on aquifer vulnerability to pollution using four parameters such as RASTIC (R: recharge, A: Aquifer media, S: soil media , T: topography, I: impact of vadose zone, and C:conductivity). The Sensitivity analysis by removing one layer at time OAT and the calculation of the effective weight are used in that context.
To solve that problem related to the vulnerability index, we going to flowing this plan:
1. Installation of package and import of libraries: sp for spatial data, raster/rgdal/rasterlayer for raster data management, and ggplot for graph production
2. Import of the tif file of the RSTIC parameters, we are limited on the fourth parameters
3. Explore the ARSTIC parameters, the function crs (projection system), res (resolution), …
4. Calculate the groundwater vulnerability index (GWI) using raster calculation
[1] “GVI min 58 -moy 84.1 -max 129 -sd 12.6”
[1] “R min 3 -moy 4.3 -max 5 -sd 0.9”
[1] “A min 4 -moy 6.5 -max 7 -sd 1”
[1] “S min 1 -moy 3.3 -max 10 -sd 2.8”
[1] “T min 10 -moy 10 -max 10 -sd 0”
[1] “I min 2 -moy 2.6 -max 8 -sd 1.8”
[1] “C min 4 -moy 5.8 -max 10 -sd 1.5”
5. Application of the removal map sensitivity analysis
[1] “V_Remove_R min 0 max 2.22 mean 1.1”
[1] “V_Remove_A min 0 max 2.94 mean 1.61”
[1] “V_Remove_S min 0 max 2.97 mean 1.91”
[1] “V_Remove_T min 0.11 max 1.78 mean 0.91”
[1] “V_Remove_I min 0.11 max 4.51 mean 1.12”
[1] “V_Remove_C min 0 max 3.52 mean 1”
[1] “V_Remove_T min 2.27 max 4.16 mean 3.29”
[1] “V_Remove_TC min 0.05 max 5.29 mean 2.5”
[1] “V_Remove_TCR min 0.1 max 8.03 mean 1.52”
[1] “V_Remove_TCRI min 0 max 8.24 mean 2.03”
[1] “V_Remove_TCRIA min 0.17 max 12.48 mean 7.63”
6. Application of the single parameter sensitivity analysis
[1] “R min 9.9 max 28 mean 20.9”
[1] “A min 10 max 31 mean 23.7”
[1] “S min 1.8 max 24 mean 7.5”
[1] “T min 7.8 max 17 mean 12.1”
[1] “I min 9.7 max 39 mean 14.8”
[1] “C min 14 max 34 mean 20.9”
Calculation of the effective weight of the parameters
Finally, i would like to inform you that you can apply this code for any problem with geospatial concept (flood, landslide, drought, energy panel location, …).
Hicham AMAR, Ing in Geomatic Sciences, Co-founder of Geoinfo4all.com
References
R. Azmi, H. Amar and I. Kacimi, “Photovoltaic Site Suitability Analysis using Analytical Hierarchy Process and Sensitivity Analysis Methods with GIS and Remote Sensing in Southern Morocco: Case of Draa-Tafilatet Region,” 2017 International Renewable and Sustainable Energy Conference (IRSEC), 2017, pp. 1-5, doi: 10.1109/IRSEC.2017.8477325.
Aït Sliman, A., Fekri, A., Laftouhi, N., Taj-Eddine, K., 2009. Utilisation des systèmes d’information géographique et du modele drastic pour l’évaluation de la vulnérabilité des eaux souterraines dans la plaine de berrechid, maroc. Geographia Technica, 8(1).
Neshat, A., Pradhan, B., Javadi, S., 2015f. Risk assessment of groundwater pollution using monte carlo approach in an agricultural region: An example from kerman plain, iran. Computers, Environment and Urban Systems, 50, 66-73. doi:10.1016/j.compenvurbsys.2014.11.004
Napolitano, P., 1995. Gis for aquifer vulnerability assessment in the piana campana, southern italy, using the drastic and sintacs methods.
Napolitano, P., Fabbri, A., 1996. Single-parameter sensitivity analysis for aquifer vulnerability assessment using drastic and sintacs. IAHS Publications-Series of Proceedings and Reports-Intern Assoc Hydrological Sciences, 235, 559-566.
Hicham AMAR, Ing in Geomatic Sciences, Co-founder of Geoinfo4all.com
متابعة القراءة...
The theoretical part for the sensitivity analysis applied on aquifer vulnerability to pollution is described in the part I.
This article will explain the step to evaluate the sensitivity analysis applied on aquifer vulnerability to pollution using four parameters such as RASTIC (R: recharge, A: Aquifer media, S: soil media , T: topography, I: impact of vadose zone, and C:conductivity). The Sensitivity analysis by removing one layer at time OAT and the calculation of the effective weight are used in that context.
To solve that problem related to the vulnerability index, we going to flowing this plan:
1. Installation of package and import of libraries: sp for spatial data, raster/rgdal/rasterlayer for raster data management, and ggplot for graph production
كود:
################Installation of package and import of libraries############
rm(list = ls())
##install.packages
# load the raster, sp, and rgdal packages
library(sp)
library(raster)
library(rgdal)
library(ggplot2)
library(RasterLayer)2. Import of the tif file of the RSTIC parameters, we are limited on the fourth parameters
كود:
#################Loading of the Tif files ################3
# set working directory to data folder
setwd("E:/Blog_Geoinfo4all/Janvier 2021/Article-DRASTIC/Data-Geoinfo4all")
# load raster in an R object called 'DEM'
DEM_R <- raster("RW.tif")
DEM_A <- raster("A.tif")
DEM_S <- raster("S.tif")
DEM_T <- raster("T.tif")
DEM_I <- raster("I.tif")
DEM_C <- raster("C.tif")3. Explore the ARSTIC parameters, the function crs (projection system), res (resolution), …
كود:
#############1. Data exploration ##############
#explore the metadata of the RASTIC files
#Summary of the tif file
summary(DEM_R)
#projection system
DEM_R@crs
#number of row and column
DEM_R@ncols; DEM_R@nrows
# resolution
res(DEM_R)
#extent
extent(DEM_R)4. Calculate the groundwater vulnerability index (GWI) using raster calculation
كود:
############Calculation of the GVI index using the six parameters ############
V <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
plot(V)
Maps <- list(V,DEM_R,DEM_A,DEM_S,DEM_T,DEM_I,DEM_C)
names_mps <- c("GVI","R","A","S","T","I","C")
# get the summary of each parameters
for ( x in 1:7 ) ## from R to C parameters
{
print(paste(names_mps[x],"min",cellStats(Maps[[x]], min) ,"-moy",round(cellStats(Maps[[x]], mean),1),"-max",cellStats(Maps[[x]], max) ,"-sd",round(cellStats(Maps[[x]], sd),1)))
}[1] “GVI min 58 -moy 84.1 -max 129 -sd 12.6”
[1] “R min 3 -moy 4.3 -max 5 -sd 0.9”
[1] “A min 4 -moy 6.5 -max 7 -sd 1”
[1] “S min 1 -moy 3.3 -max 10 -sd 2.8”
[1] “T min 10 -moy 10 -max 10 -sd 0”
[1] “I min 2 -moy 2.6 -max 8 -sd 1.8”
[1] “C min 4 -moy 5.8 -max 10 -sd 1.5”
كود:
#Mapping of the GVI and RASTIC parameters
par(mfrow=c(2,4))
print("GVI and RASTIC maps ")
dev.new()
par(mfrow=c(2,4))
for ( x in 1:7 )
{
eq <- names_mps[x]
plot( Maps[[x]] , main= eq, xlab= paste("min-mean-max",cellStats(Maps[[x]], min) ,"-",round(cellStats(Maps[[x]], mean),0),"-",cellStats(Maps[[x]], max)) )
}5. Application of the removal map sensitivity analysis
كود:
######################Sensitivity analysis################################
######################removal map : SA_R1################################
#unperturbed vulnerability index
V <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
########################STEEP I : removing OAT ############
print("Sensitivity analysis by removing one layer at time OAT")
################
# perturbed vulnerability index by removing R parameters : V_R
V_R <- 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
################
# perturbed vulnerability index by removing A parameters V_A
V_A <- 4*DEM_R + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
################
V_S <- 4*DEM_R + 3*DEM_A + 1*DEM_T + 5*DEM_I + 3*DEM_C
################
V_T <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 5*DEM_I + 3*DEM_C
################
V_I <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 3*DEM_C
################
V_C <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I
################
SA_R1 <- list(V,V_R,V_A,V_S,V_T,V_I,V_C)
names <- c("GVI","V_Remove_R","V_Remove_A","V_Remove_S","V_Remove_T","V_Remove_I","V_Remove_C","V")
print("mapping of the maps by removing one layer at time OAT")
dev.new()
par(mfrow=c(2,4))
for ( x in 1:8 ) ## 1:20 = > 20 first combinations
{
eq <- names[x]
plot( SA_R1[[x]] , main= eq, xlab= paste("min-mean-max",cellStats(SA_R1[[x]], min) ,"-",round(cellStats(SA_R1[[x]], mean),0),"-",cellStats(SA_R1[[x]], max)) )
}
كود:
############# S = (|V/N-V'/n|/N)*100
#Where the V and v are the unperturbed and perturbed vulnerability index
#using N and n parameters respectively.
print("calculating the sensitivity map (S) ")
#V/6, unperturbed GVI with 6 parameters
#SA_R1[[x]]/5 perturbed GVI with 5 parameters (one parameter removed)
for (x in 2:7)
{
S <- abs((V/6) - (SA_R1[[x]]/5))/V
S
print(paste(names[x],"min",round((cellStats(S, min)*100), 2)," max",round((cellStats(S, max)*100), 2)," mean",round((cellStats(S, mean)*100), 2)))
}[1] “V_Remove_R min 0 max 2.22 mean 1.1”
[1] “V_Remove_A min 0 max 2.94 mean 1.61”
[1] “V_Remove_S min 0 max 2.97 mean 1.91”
[1] “V_Remove_T min 0.11 max 1.78 mean 0.91”
[1] “V_Remove_I min 0.11 max 4.51 mean 1.12”
[1] “V_Remove_C min 0 max 3.52 mean 1”
كود:
########################STEEP II : removing of the parameter with less variation ############
#############################################
##V <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
## PARAMETERs WITH LESS VARIATION (mean)
## T(0.91), C(1), R(1.1),I(1.12), A(1.61), S(1.91)
V <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
V_T <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 5*DEM_I + 3*DEM_C
################
V_TC <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 5*DEM_I
################
V_TCR <- 3*DEM_A + 2*DEM_S + 5*DEM_I
################
V_TCRI <- 3*DEM_A + 2*DEM_S
################
V_TCRIA <- 2*DEM_S
SA_RM <- list( V_T,V_TC,V_TCR,V_TCRI,V_TCRIA, V )
names <- c("V_Remove_T","V_Remove_TC","V_Remove_TCR","V_Remove_TCRI","V_Remove_TCRIA", "V")
print("mapping of the maps by removing layers with less variation")
dev.new()
par(mfrow=c(2,4))
for ( x in 1:6 ) ## 1:20 = > 20 first combinations
{
eq <- names[x]
plot( SA_RM[[x]] , main= eq, xlab= paste("min-mean-max",cellStats(SA_RM[[x]], min) ,"-",round(cellStats(SA_RM[[x]], mean),0),"-",cellStats(SA_RM[[x]], max)) )
}
كود:
print("Sensitivity analyis by removing multiple layers at one time : with less variation ")
for (x in 1:5)
{
z <- c(5:1)
z
S_Rmv <- abs((V/7) - (SA_RM[[x]]/ z[x]))/V
S_Rmv
print(paste(names[x],"min",round((cellStats(S_Rmv, min)*100), 2)," max",round((cellStats(S_Rmv, max)*100), 2)," mean",round((cellStats(S_Rmv, mean)*100), 2)))
}[1] “V_Remove_T min 2.27 max 4.16 mean 3.29”
[1] “V_Remove_TC min 0.05 max 5.29 mean 2.5”
[1] “V_Remove_TCR min 0.1 max 8.03 mean 1.52”
[1] “V_Remove_TCRI min 0 max 8.24 mean 2.03”
[1] “V_Remove_TCRIA min 0.17 max 12.48 mean 7.63”
6. Application of the single parameter sensitivity analysis
كود:
###################################
#####W = (Pr*Pw/V)*100
#Pr and Pw are respectively the ranting, the weight of the parameters P
###################################
##V <- 4*DEM_R + 3*DEM_A + 2*DEM_S + 1*DEM_T + 5*DEM_I + 3*DEM_C
poids <- c(4,3,2,1,5,3)
name_layrs <- c("R", "A","S", "T", "I" ,"C")
layrs <- list( DEM_R, DEM_A,DEM_S, DEM_T, DEM_I ,DEM_C)
print("Calculation of the effective wiegth of the parameters ")
for (x in 1:6)
{
W_eff <- ((layrs[[x]]* poids[x])/V)*100
print(paste(name_layrs[x],"min",round((cellStats(W_eff, min)), 1)," max",
round((cellStats(W_eff, max)), 0)," mean",round((cellStats(W_eff, mean)),1)))
}[1] “R min 9.9 max 28 mean 20.9”
[1] “A min 10 max 31 mean 23.7”
[1] “S min 1.8 max 24 mean 7.5”
[1] “T min 7.8 max 17 mean 12.1”
[1] “I min 9.7 max 39 mean 14.8”
[1] “C min 14 max 34 mean 20.9”
Calculation of the effective weight of the parameters
| Parameters | Theoretical weight | Theoretical weight (%) | Effective weight (%) | Effective weight |
| R | 4 | 22,2 | 20,9 | 4 |
| A | 3 | 16,7 | 23,7 | 4 |
| S | 2 | 11,1 | 7,5 | 1 |
| T | 1 | 5,6 | 12,1 | 2 |
| I | 5 | 27,8 | 14,8 | 3 |
| C | 3 | 16,7 | 20,9 | 4 |
| Total | 18 | 100,0 | 555,6 | 18,0 |
Finally, i would like to inform you that you can apply this code for any problem with geospatial concept (flood, landslide, drought, energy panel location, …).
Hicham AMAR, Ing in Geomatic Sciences, Co-founder of Geoinfo4all.com
References
R. Azmi, H. Amar and I. Kacimi, “Photovoltaic Site Suitability Analysis using Analytical Hierarchy Process and Sensitivity Analysis Methods with GIS and Remote Sensing in Southern Morocco: Case of Draa-Tafilatet Region,” 2017 International Renewable and Sustainable Energy Conference (IRSEC), 2017, pp. 1-5, doi: 10.1109/IRSEC.2017.8477325.
Aït Sliman, A., Fekri, A., Laftouhi, N., Taj-Eddine, K., 2009. Utilisation des systèmes d’information géographique et du modele drastic pour l’évaluation de la vulnérabilité des eaux souterraines dans la plaine de berrechid, maroc. Geographia Technica, 8(1).
Neshat, A., Pradhan, B., Javadi, S., 2015f. Risk assessment of groundwater pollution using monte carlo approach in an agricultural region: An example from kerman plain, iran. Computers, Environment and Urban Systems, 50, 66-73. doi:10.1016/j.compenvurbsys.2014.11.004
Napolitano, P., 1995. Gis for aquifer vulnerability assessment in the piana campana, southern italy, using the drastic and sintacs methods.
Napolitano, P., Fabbri, A., 1996. Single-parameter sensitivity analysis for aquifer vulnerability assessment using drastic and sintacs. IAHS Publications-Series of Proceedings and Reports-Intern Assoc Hydrological Sciences, 235, 559-566.
Hicham AMAR, Ing in Geomatic Sciences, Co-founder of Geoinfo4all.com
متابعة القراءة...