validating_fldrisk
CDL_Jeeno
2022-10-29
Background
This documents explains how the differnet layers of socio-economic variables are brought together to evaluate the overall vulnerability. The overall layers is multiplied with the monthly flood maps to get the monthly flood risk. These risk maps are useful for (1) identifying the hot spots where more attention are required , and (2) for identifying the critical months of the selected administrative unit.
For the work we have to load the following libraries:
packages <- c("terra","tidyverse","dplyr","sp","sf","rgdal","tidyverse","reshape")First, we calculate the overall vulnerability for all the years of study. Here we read the damages from the
frims (currently these are crop area affected, population affected and inmates in relief camps but will be soon replaced all all the data collected in FRIMS);
the demographic variables (Population , Elederly population, Below 5 years populaiton, Sex Ratio, AAY population, deprived population (Deprived households - Only one room with kucha walls and kucha roof; No adult member between age 16 to 59; Female headed households with no adult male member between age 16 to 59; Disabled member and no able bodied adult member; SC/ST households; No literate adult above 25 years; Landless households deriving major part of their income from manual casual labour), households with neither land line or mobile phone, population with no drinking water or sanitation facilities;
the physical infrastructure variables (Arterial Roads, Local Roads, Railway network, built density, vegetation density, distribution of hospitals))
Here we calculate vulnerability by adding up all the layers. However, there is a plan to : (1) use the damages to weigh the socio-economic layers , or (2 add all variables as part of structural equation model)
Plans ahead - Add animal stats
#run the loop for all years
for (yr in 2018:2018) {
#DAMAGE
cropDam <- rast(sprintf("outputs_vulnerability/State_crop_affected_%s_RS.tif",yr))
plot(cropDam)
popAff <- rast(sprintf("outputs_vulnerability/State_pop_affected_%s_RS.tif",yr))
inmatesRef <- rast(sprintf("outputs_vulnerability/State_inmated_Damage_%s_RS.tif",yr))
#overall damage
damagelayer <- sum(cropDam ,popAff,popAff, na.rm = T)
plot(damagelayer, main = sprintf("Overall Damages - %s",yr))
plot(rast(sprintf("outputs_vulnerability/damagelayer_%s_Aug.tif",yr)))
#DEMOGRAPHY
aged_2020 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_AgedPop_%s.tif",yr))
pop_2022 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_pop_%s.tif",yr))
popden <- rast(sprintf("conditioningfactors/RESAMPL_Aug_popDen_2022.tif"))
young_2020 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_YoungPop_%s.tif",yr))
sexRatio2020 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_sexratioPop_%s.tif",yr))
aay_2020 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_aayPop_%s.tif",yr))
noph_2021 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_nophPop_%s.tif",yr))
DEPRV_2021 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_deprPop_%s.tif",yr))
nodrnkng2021 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_nodrnkPop_%s.tif",yr))
nosanit2021 <- rast(sprintf("conditioningfactors/RESAMPL_Aug_nosanitPop_%s.tif",yr))
#overall demographic vulnerability
demolayer <- sum(pop_2022, popden, aged_2020 ,young_2020, sexRatio2020, aay_2020,
noph_2021, DEPRV_2021,nodrnkng2021,nosanit2021, na.rm = T)
plot(demolayer, main = sprintf("Overall demographic characteristics"))
#INFRA
hospital <- rast(sprintf("outputs_vulnerability/%s_State_005_VulRS.tif",'hosptial'))
embnk <- rast(sprintf("outputs_vulnerability/%s_State_005_VulRS.tif",'embankment'))
arterial <- rast(sprintf("outputs_vulnerability/%s_State_005_VulRS.tif",'arterial'))
local <- rast(sprintf("outputs_vulnerability/%s_State_005_VulRS.tif",'local'))
rail <- rast(sprintf("outputs_vulnerability/%s_State_005_VulRS.tif",'rail'))
rd_dist <- rast("conditioningfactors/Method2RoadDist_RESAMPL_Aug_RL_wt_NORM.tif")
rd_intrsc <- rast("conditioningfactors/Method2Roadintersection_RESAMPL_Aug_RL_wt_NORM.tif")
natresrv_dist <- rast("conditioningfactors/Method2NatReserv_RESAMPL_Aug_RL_wt_NORM.tif")
ndvi <- rast(sprintf("conditioningfactors/Method2veg_RESAMPL_Aug_%s_RL_wt.tif",yr))
ndbi <- rast(sprintf("conditioningfactors/Method2built_RESAMPL_Aug_%s_RL_wt.tif",yr))
infrlayer <- sum(hospital, ndvi, ndbi, embnk, arterial, local, rail, rd_dist, rd_intrsc, natresrv_dist, na.rm = T)
plot(infrlayer, col = rev(topo.colors(5)))
#overall vulnerability - add all the layers
overall_vul <- sum(damagelayer, demolayer, infrlayer)
plot(overall_vul, col = rev(rainbow(7)), main = sprintf("Composite infrastructre map - %s",yr))
}
Crop to the district required district.
districts <- vect("assam_district_35.shp")
#view(districts)
#Codes of pilot districts
#Morigoan - 3
#Kamrup (R) - 30
#Jorhat - 22
#Darrang - 26
###****change to district code###
d <- 22
districts[d,2]$district## [1] "JORHAT" #change to district number
crop_district <- (districts[d,2])
plot(crop_district)
s_files <- dir("~/bhuvanscraper/FPM")Prepare a user-defined monthly inundation months from daily flood images scraped from Bhuvan. The function selects images of a particular year and month. The functions with no images are assigned a dummy raster with ‘0’ values. The months with flood images are added together into a single image and saved.
montlyfldprob <- function(s_files, mnth){
sample_img <- rast("conditioningfactors/slope.tif")
sample_img <- terra::crop(sample_img , crop_district)
sample_img <- terra::mask(sample_img , crop_district , overwrite = T)
dummy <- rast(ext=ext(sample_img),crs=crs(sample_img, proj=T, describe=FALSE, parse=FALSE),
nrow=dim(sample_img)[1],ncol=dim(sample_img)[2])
values(dummy) <- 0
dummy <- terra::crop(dummy , crop_district)
dummy <- terra::mask(dummy , crop_district , overwrite = T)
names(dummy) <- 'Inundation'
plot(dummy)
#mnth <- '04'
split <- strsplit(s_files,"_")
pos <- c(1:1)
for (i in 1:length(split)) {
#i <- 1
pos[i] <- ifelse(split[[i]][1] == yr,i,NA)
}
pos <- which(!is.na(pos[]))
s_mnth_files <- s_files[pos]
s_mnth_files
split <- strsplit(s_mnth_files,"_")
for (i in 1:length(split)) {
#i <- 1
pos[i] <- ifelse(split[[i]][2] == mnth,i,NA)
}
pos <- which(!is.na(pos[]))
class(pos)
s_mnth_files <- s_mnth_files[pos]
s_mnth_files
while (length(pos) == 0) {
crop_flood_mnth <- dummy
return(0)
}
while (length(pos) != 0) {
extension <- "/home/jeeno/bhuvanscraper/FPM/"
#combining flood images of month "06"
final_files <- paste(extension,s_mnth_files,sep ="")
flooded_2022_mnth <- rast(final_files)
flooded_2022_mnth <- terra::crop(flooded_2022_mnth , crop_district)
flooded_2022_mnth <- terra::mask(flooded_2022_mnth , crop_district , overwrite = T)
img <- flooded_2022_mnth
plot(img)
img[img == 0] <- NA
crop_flood <- img#/ max(img[], na.rm = T)
crop_flood_mnth <- resample(crop_flood, sample_img)
plot(crop_flood_mnth)
return(0)
plot(rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,mnth,districts[d,2]$district)), main = sprintf("%s",yr))
}
} The user-defined function is called for each months from April to October for different years. The overall vulnerability map is cropped to the selected district. The flood map is multiplied with the overall vulnerability to get the flood risk. Extract the mean value of flood frequency and flood risk for each sub-district. Also, the month in which a subdistrict had faced the highest flood risk in each year is saved in a table.
s_files <- dir("~/bhuvanscraper/FPM")
#run for all years
#for (yr in 2018:2018) {
yr <- 2018
montlyfldprob(s_files, '04')## [1] 0 fld_04 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"04",districts[d,2]$district))
plot(fld_04)
montlyfldprob(s_files, '05')## [1] 0 fld_05 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"05",districts[d,2]$district))
plot(fld_05)
montlyfldprob(s_files, '06')
## [1] 0 fld_06 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"06",districts[d,2]$district))
plot(fld_06)
montlyfldprob(s_files, '07')
## [1] 0 fld_07 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"07",districts[d,2]$district))
plot(fld_07)
montlyfldprob(s_files, '08')
## [1] 0 fld_08 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"08",districts[d,2]$district))
plot(fld_08)
montlyfldprob(s_files, '09')
## [1] 0 fld_09 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"09",districts[d,2]$district))
plot(fld_09)
montlyfldprob(s_files, '10')
## [1] 0 fld_10 <- rast(sprintf("conditioningfactors/floodedareas_grtr1in5_%s_%sAug_RS_%s.tif",yr,"10",districts[d,2]$district))
plot(fld_10)
vul <- rast(sprintf("outputs_vulnerability/overall_vulnerability_%s_Aug.tif",yr))
vul <- terra::crop(vul , crop_district)
vul <- terra::mask(vul , crop_district , overwrite = T)
plot(vul)
risk_04 <- fld_04 * vul
plot(risk_04)
risk_05 <- fld_05 * vul
plot(risk_05)
risk_06 <- fld_06 * vul
plot(risk_06)
risk_07 <- fld_07 * vul
plot(risk_07)
risk_08 <- fld_08 * vul
plot(risk_08)
risk_09 <- fld_09 * vul
plot(risk_09)
risk_10 <- fld_10 * vul
plot(risk_10)
#extract fld and risk sub-district wise
subdistricts <- vect("ASSAM_SUBDISTRICT_4326.shp")
subdistricts$dtname <- toupper(subdistricts$dtname)
#enter the district
k <- which(subdistricts$dtname == districts[d,2]$district)
sel_subdist <- subdistricts[k,]
plot(sel_subdist, main = "Sub-districts")
#model - subdistrict_flooded areas
subd_m <- terra::extract(fld_04, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_04 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_04, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_04 <- subd_a$Inundation
subd_m <- terra::extract(fld_05, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_05 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_05, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_05 <- subd_a$Inundation
subd_m <- terra::extract(fld_06, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_06 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_06, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_06 <- subd_a$Inundation
subd_m <- terra::extract(fld_07, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_07 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_07, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_07 <- subd_a$Inundation
subd_m <- terra::extract(fld_08, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_08 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_08, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_08 <- subd_a$Inundation
subd_m <- terra::extract(fld_09, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_09 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_09, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_09 <- subd_a$Inundation
subd_m <- terra::extract(fld_10, sel_subdist, 'mean', na.rm=T)
sel_subdist$fld_10 <- subd_m$Inundation
#actual - subdistrict_flooded areas
subd_a <- terra::extract(risk_10, sel_subdist, 'mean', na.rm=T)
sel_subdist$risk_10 <- subd_a$Inundation
sel_subdist$year <- rep(yr, length(sel_subdist))
#prepared table
sel_subdist## class : SpatVector
## geometry : polygons
## dimensions : 5, 28 (geometries, attributes)
## extent : 93.95672, 94.60326, 26.34746, 27.18374 (xmin, xmax, ymin, ymax)
## coord. ref. : lon/lat WGS 84 (EPSG:4326)
## names : OBJECTID stcode11 dtcode11 sdtcode11 stname dtname sdtname
## type : <num> <chr> <chr> <chr> <chr> <chr> <chr>
## values : 2056 18 312 02077 ASSAM JORHAT Jorhat West
## 2057 18 312 02078 ASSAM JORHAT Jorhat East
## 2058 18 312 02079 ASSAM JORHAT Teok
## Subdt_LGD Dist_LGD State_LGD (and 18 more)
## <int> <int> <int>
## 2077 290 18
## 2078 290 18
## 2079 290 18 critical_mnths <- sel_subdist[,c("sdtname", "risk_10", "risk_09", "risk_08", "risk_07", "risk_06", "risk_05", "risk_04", "year")]
critical_mnths <- data.frame(critical_mnths)
critical_mnths[is.na(critical_mnths)] <- 0
critical_mnths_db <- melt(critical_mnths, id.vars = c("sdtname","year"), variable_name = "fld")
critical_mnths_names <- critical_mnths_db %>% group_by(sdtname) %>%
filter(value == max(value)) %>%
dplyr::select(sdtname, fld, value)
critical_mnths_names$year = rep(yr, nrow(critical_mnths_names))
#critical months
critical_mnths_names## # A tibble: 11 × 4
## # Groups: sdtname [5]
## sdtname fld value year
## <chr> <fct> <dbl> <dbl>
## 1 Majuli risk_10 0 2018
## 2 Jorhat West risk_09 4.61 2018
## 3 Jorhat East risk_09 4.48 2018
## 4 Majuli risk_09 0 2018
## 5 Majuli risk_08 0 2018
## 6 Teok risk_07 4.12 2018
## 7 Titabor risk_07 3.21 2018
## 8 Majuli risk_07 0 2018
## 9 Majuli risk_06 0 2018
## 10 Majuli risk_05 0 2018
## 11 Majuli risk_04 0 2018 #}