B Aim2: BOUNCEBACK RESILIENCE
rm(list=ls())
library(raster)
library(sf)
library(rasterVis)
library(tidyverse)
library(haven)
library(stringr)
library(scales)
library(colorspace)
library(MASS)
B.1 Visuals
a <- dat_year %>%
#filter(!is.na(nl_viirs_2)) %>%
filter(year>2011, year<2019) %>%
ungroup() %>%
group_by(id_loc) %>%
mutate(nl_viirs_2.baseline = first(nl_viirs_2),
VIIRS = cut(nl_viirs_2.baseline, breaks=c(-Inf,5,Inf)),
MEI_ike_udel.baseline = first(MEI_ike_udel),
MEI = cut(MEI_ike_udel.baseline, breaks=c(-Inf,0.6,Inf)),
cat = paste0(VIIRS,MEI))
B.1.1 vivax malaria
P. vivax malaria trends relative to MEI and VIIRS categories at baseline (2012). Each line represents a location (village).
a %>%
ggplot(aes(x = year, y = viv_rate, group=id_loc, col = factor(cat))) +
geom_point() +
geom_line() +
theme_bw() +
scale_color_discrete_sequential(palette = "Red-Blue") +
guides(col=F) +
labs(y="Malaria rate") +
facet_grid(VIIRS~MEI, labeller=label_both)
B.1.2 falciparum malaria
P. falciparum malaria trends relative to MEI and VIIRS categories at baseline (2012). Each line represents a location (village).
a %>%
ggplot(aes(x = year, y = fal_rate, group=id_loc, col = factor(cat))) +
geom_point() +
geom_line() +
theme_bw() +
scale_color_discrete_sequential(palette = "Red-Blue") +
guides(col=F) +
labs(y="Malaria rate") +
facet_grid(VIIRS~MEI, labeller=label_both)
B.2 Slope calcualtions
A negative binomial model per location.
c<-glm.nb(viv ~ year, data = a, na.action=na.exclude)
c$coefficients[2] <-NA
glm.nb.safe <- possibly(glm.nb, otherwise = c)
dt_s <- dat_year %>%
filter(year>2011, year<2019) %>%
ungroup() %>%
#filter(NOMBDIST=="IQUITOS") %>%
dplyr::select(id_loc, year, viv, fal, nrohab) %>%
gather(var, value, viv:fal) %>%
nest(-id_loc, -var) %>%
mutate(model = map(data, ~glm.nb.safe(value ~ year + offset(log(nrohab)), data = ., na.action=na.exclude)),
slope = lapply(model, function(x) { x$coefficients[2] })) %>%
dplyr::select(id_loc, var, slope) %>%
unnest(cols = c(slope)) %>%
spread(var,slope) %>%
rename(slope_fal = fal, slope_viv = viv)
dt_s.c <- a %>%
filter(year==2012) %>%
inner_join(dt_s, by="id_loc") %>%
inner_join(st_set_geometry(dat, NULL) %>%
filter(year==2012) %>%
distinct(id_loc, .keep_all = T), by="id_loc")
B.2.1 vivax malaria
Plotting the slope (linear predictor) vs. VIIRS at baseline
dt_s.c %>%
ggplot(aes(y = slope_viv, x = nl_viirs_2.baseline, col = MEI_ike_udel.baseline)) +
geom_point(size = 2, alpha=.4) +
scale_x_log10(label = number_format(accuracy = .1)) +
labs(y= "Negative Binomial linear predictor", x = "VIIRS at baseline (Log-scale)", col = "Malaria\nEcology\nIndex") +
scale_color_continuous_sequential(palette = "OrRd") +
theme_bw() +
facet_wrap(province.x~., ncol = 4)
B.2.2 falciparum malaria
Plotting the slope (linear predictor) vs. VIIRS at baseline
dt_s.c %>%
ggplot(aes(y = slope_fal, x = nl_viirs_2.baseline, col = MEI_ike_udel.baseline)) +
geom_point(size = 2, alpha=.4) +
scale_x_log10(label = number_format(accuracy = .1)) +
labs(y= "Negative Binomial linear predictor", x = "VIIRS at baseline (Log-scale)", col = "Malaria\nEcology\nIndex") +
scale_color_continuous_sequential(palette = "OrRd") +
theme_bw() +
facet_wrap(province.x~., ncol = 4)