Chapter 1 Data
- Cases: Maria Eugenia
- Shapefile: HDX
- Nightime lights: VIIRS
- Malaria Ecology Index
rm(list=ls())
library(tidyverse)
library(sf)
library(readxl)
library(stringr)
area <- st_read("./_dat/ven_admbnda_adm3_20180502/ven_admbnda_adm3_20180502.shp") %>%
filter(ADM1_ES=="Bolivar") %>%
mutate(ADM3_ES = toupper(ADM3_ES),
ADM2_ES = toupper(ADM2_ES))
dat <- read_excel("./_dat/Malaria_Data_Southern_Venezuela_1996-2016.xlsx") %>%
rename(ADM3_ES = Parroquia, ADM2_ES = Municipio) %>%
mutate(ADM3_ES = str_replace(ADM3_ES, "S\\.C\\.", "SECCION CAPITAL"),
ADM3_ES = str_replace(ADM3_ES, "A\\.E\\. ", "ANDRES ELOY "),
ADM3_ES = str_replace(ADM3_ES, "C\\. P\\. P\\.", "CAPITAL PADRE PEDRO"),
ADM3_ES = str_replace(ADM3_ES, "CEDEñO","CEDEÑO"),
ADM3_ES = str_replace(ADM3_ES, "SABANA","SABAN"),
ADM3_ES = str_replace(ADM3_ES, "C\\.","CAPITAL"),
ADM2_ES = str_replace(ADM2_ES, "CEDEñO","CEDEÑO"),
ADM2_ES = str_replace(ADM2_ES, "Heres","HERES")) %>%
filter(ADM3_ES!="CINCO DE JULIO") %>%
filter(ADM3_ES!="SALOM")1.1 GEE pre-processing
#Validate
list.dat <- dat %>%
distinct(ADM2_ES,ADM3_ES) %>%
mutate(dat = 1)
list.map <- area %>%
st_set_geometry(NULL) %>%
distinct(ADM2_ES,ADM3_ES) %>%
mutate(map = 1)
mm <- list.dat %>%
full_join(list.map, by = c("ADM2_ES","ADM3_ES"))
#Export
sf <- area %>%
dplyr::select(ADM2_ES,ADM3_ES)
st_write(sf, "./_dat/ven_shape/ven_shape.shp")1.2 GEE
#remotes::install_github("r-spatial/rgee")
#ee_install()
library(tidyverse)
library(rgee)
library(sf)
ee_Initialize()
#rgee::ee_reattach()
# Change your FeatureCollection -------------------------------------------
years <- 2012:2019
loreto_sf <- st_read("./_dat/ven_shape/ven_shape.shp")
# Run ee_extract using a batch size ---------------------------------------
rain_data <- list()
for (year in as.character(years)) {
cat('processing data for the year:', year,"\n")
init_date <- sprintf('%s-01-01',year)
last_date <- sprintf('%s-12-31',year)
terraclimate <- ee$ImageCollection('NOAA/VIIRS/DNB/MONTHLY_V1/VCMCFG')$ # SELECT IMAGE COLLECTION
filterDate(init_date, last_date)$
map(function(x) x$select("avg_rad")) # SELECT BANDS
pp_data <- ee_extract(x = terraclimate,y = loreto_sf,scale = 5000) # fun = ee$Reducer$mean()
rain_data[[year]] <- pp_data
}
# Wrangling data ----------------------------------------------------------
## Flatten list of years
flat<- rain_data %>%
reduce(bind_cols)
# data validated by subset loreto_sf to single district. index from rain_data in the same order as loreto_sf. Correct to use bind_cols
loreto_with_rain <- loreto_sf %>%
bind_cols(flat)
plot(loreto_with_rain['X199910'])
gdata::keep(loreto_with_rain, sure = T)
saveRDS(loreto_with_rain, "./_dat/nl_ven.rds")1.3 Malaria Ecology Index (MEI) data
library(haven)
library(raster)
library(sf)
mei <- read_dta("./_dat/MEI_annualgrid_VEN_2012-2017.dta") %>%
zap_formats() %>%
filter(year>2011, year<2017) %>%
st_as_sf(coords = c("lon", "lat")) %>%
nest(-year) %>%
mutate(sp = map(data, ~as(., "Spatial")),
#change extent or resolution
raster = map(sp, ~rasterize(., raster(extent(area), res = .5), field = "MEI_udel",
fun="sum")),
extracted = map(raster, ~raster::extract(., area, fun=mean, df=TRUE, sp=T)),
data.extracted = map(extracted, ~.@data)) %>%
dplyr::select(year, data.extracted) %>%
unnest(cols = c(data.extracted)) %>%
dplyr::select(year, ADM2_ES, ADM3_ES, layer) %>%
rename(MEI_udel=layer)
#a<- mei$data.extracted[[1]]1.4 Climat covars
# remotes::install_github("healthinnovation/lis")
library(rgee)
library(lis)
ee_Initialize()
region <- st_read("./_dat/ven_shape/ven_shape.shp")
region_ee <- pol_as_ee(region , id = 'distr' , simplify = 1000)
data_tmmx <- region_ee %>% get_climate(
to = "2012-01-01", from = "2016-12-31",
by = "year", band = "tmmx", fun = "mean")
data_tmmx2 <- data_tmmx %>%
rename(tmmx2013 = tmmx.1,
tmmx2014 = tmmx.2,
tmmx2015 = tmmx.3,
tmmx2016 = tmmx.4) %>%
cbind(region %>% st_set_geometry(NULL)) %>%
gather(Year, tmmx, tmmx2012:tmmx2016) %>%
mutate(Year = as.numeric(str_remove(Year, "tmmx")),
tmmx = tmmx/10)
data_pr <- region_ee %>% get_climate(
to = "2012-01-01", from = "2016-12-31",
by = "year", band = "pr", fun = "sum")
data_pr2 <- data_pr %>%
rename(pr2013 = pr.1,
pr2014 = pr.2,
pr2015 = pr.3,
pr2016 = pr.4) %>%
cbind(region %>% st_set_geometry(NULL)) %>%
gather(Year, pr, pr2012:pr2016) %>%
mutate(Year = as.numeric(str_remove(Year, "pr")))
dd <- data_tmmx2 %>%
inner_join(data_pr2, by = c("id", "ADM2_ES", "ADM3_ES", "Year")) %>%
select(-id)
write.csv(dd, "./_dat/env_covar.csv", na = "", row.names = F)