Chapter 8 Spatial Analysis (U5)
8.1 Processing
rm(list = ls())
library(tidyverse)
library(spdep)
library(tmap)
library(sf)
library(leaflet)
dat_ci_map <- readRDS("./_dat/dat_ci_map_v2_u5.rds") %>%
filter(sample.size>=10)
dat_edu_map <- readRDS("./_dat/dat_edu_map_v2_u5.rds") %>%
filter(sample.size>=10) %>%
filter(!is.infinite(rii))
source("./mal_ineq_fun.R")
breaks <- c(-Inf, -2.58, -1.96, -1.65, 1.65, 1.96, 2.58, Inf)
labels <- c("Cold spot: 99% confidence", "Cold spot: 95% confidence", "Cold spot: 90% confidence", "PSU Not significant","Hot spot: 90% confidence", "Hot spot: 95% confidence", "Hot spot: 99% confidence")8.2 Wealth
dat_ci_sp <- as(dat_ci_map, "Spatial")
coords<-coordinates(dat_ci_sp)
IDs<-row.names(as.data.frame(coords))
Neigh_nb<-knn2nb(knearneigh(coords, k=1, longlat = TRUE), row.names=IDs)
dsts<-unlist(nbdists(Neigh_nb,coords))
max_1nn<-max(dsts)
Neigh_kd1<-dnearneigh(coords,d1=0, d2=max_1nn, row.names=IDs)
self<-nb2listw(Neigh_kd1, style="W", zero.policy = T)
w.getis <- dat_ci_map %>%
st_set_geometry(NULL) %>%
nest() %>%
crossing(var1 = c("ci_val", "sii", "rii")) %>%
mutate(data_sub = map2(.x = data, .y = var1, .f = ~dplyr::select(.x, .y, country, cluster_number,
prev, lat, long) %>%
rename(output = 1)),
LISA = map(.x = data_sub, .f = ~localG(.x$output, self)),
clust_LISA = map(.x = LISA, .f = ~cut(.x, include.lowest = TRUE,
breaks = breaks, labels = labels))) %>%
dplyr::select(data_sub, var1, LISA, clust_LISA) %>%
unnest()8.2.1 W-CI
dat_w_map_ci <- w.getis %>%
filter(var1 =="ci_val") %>%
st_as_sf(coords = c("long", "lat"), crs = 4326) %>%
filter(!is.na(LISA)) %>%
dplyr::mutate(lat = sf::st_coordinates(.)[,2],
long = sf::st_coordinates(.)[,1])
pal <- colorFactor("RdBu", dat_w_map_ci$clust_LISA, reverse = T)
dat_w_map_ci %>%
leaflet() %>%
addTiles() %>%
addCircles(lng = ~long, lat = ~lat, color = ~ pal(clust_LISA),
popup = paste("Country", dat_w_map_ci$country, "<br>",
"Cluster:", dat_w_map_ci$cluster_number, "<br>",
"Malaria Prevalence:", dat_w_map_ci$prev, "<br>",
"CI:", dat_w_map_ci$output, "<br>",
"LISA:", dat_w_map_ci$LISA)) %>%
addLegend("bottomright",
pal = pal, values = ~clust_LISA,
title = "LISA - Wealth CI") %>%
addScaleBar(position = c("bottomleft"))8.2.2 W-SII
dat_w_map_sii <- w.getis %>%
filter(var1 =="sii") %>%
st_as_sf(coords = c("long", "lat"), crs = 4326) %>%
filter(!is.na(LISA)) %>%
dplyr::mutate(lat = sf::st_coordinates(.)[,2],
long = sf::st_coordinates(.)[,1])
pal <- colorFactor("RdBu", dat_w_map_sii$clust_LISA, reverse = T)
dat_w_map_sii %>%
leaflet() %>%
addTiles() %>%
addCircles(lng = ~long, lat = ~lat, color = ~ pal(clust_LISA),
popup = paste("Country", dat_w_map_sii$country, "<br>",
"Cluster:", dat_w_map_sii$cluster_number, "<br>",
"Malaria Prevalence:", dat_w_map_sii$prev, "<br>",
"CI:", dat_w_map_sii$output, "<br>",
"LISA:", dat_w_map_sii$LISA)) %>%
addLegend("bottomright",
pal = pal, values = ~clust_LISA,
title = "LISA - Wealth SII") %>%
addScaleBar(position = c("bottomleft"))8.2.3 W-RII
dat_w_map_rii <- w.getis %>%
filter(var1 =="rii") %>%
st_as_sf(coords = c("long", "lat"), crs = 4326) %>%
filter(!is.na(LISA)) %>%
dplyr::mutate(lat = sf::st_coordinates(.)[,2],
long = sf::st_coordinates(.)[,1])
pal <- colorFactor("RdBu", dat_w_map_rii$clust_LISA, reverse = T)
dat_w_map_rii %>%
leaflet() %>%
addTiles() %>%
addCircles(lng = ~long, lat = ~lat, color = ~ pal(clust_LISA),
popup = paste("Country", dat_w_map_rii$country, "<br>",
"Cluster:", dat_w_map_rii$cluster_number, "<br>",
"Malaria Prevalence:", dat_w_map_rii$prev, "<br>",
"CI:", dat_w_map_rii$output, "<br>",
"LISA:", dat_w_map_rii$LISA)) %>%
addLegend("bottomright",
pal = pal, values = ~clust_LISA,
title = "LISA - Wealth RII") %>%
addScaleBar(position = c("bottomleft"))8.3 Education
dat_edu_sp <- as(dat_edu_map, "Spatial")
coords<-coordinates(dat_edu_sp)
IDs<-row.names(as.data.frame(coords))
Neigh_nb<-knn2nb(knearneigh(coords, k=1, longlat = TRUE), row.names=IDs)
dsts<-unlist(nbdists(Neigh_nb,coords))
max_1nn<-max(dsts)
Neigh_kd1<-dnearneigh(coords,d1=0, d2=max_1nn, row.names=IDs)
self<-nb2listw(Neigh_kd1, style="W", zero.policy = T)
w.getis_e <- dat_edu_map %>%
st_set_geometry(NULL) %>%
nest() %>%
crossing(var1 = c("ci_val", "sii", "rii")) %>%
mutate(data_sub = map2(.x = data, .y = var1, .f = ~dplyr::select(.x, .y, country, cluster_number,
prev, lat, long) %>%
rename(output = 1)),
LISA = map(.x = data_sub, .f = ~localG(.x$output, self)),
clust_LISA = map(.x = LISA, .f = ~cut(.x, include.lowest = TRUE,
breaks = breaks, labels = labels))) %>%
dplyr::select(data_sub, var1, LISA, clust_LISA) %>%
unnest()8.3.1 E-CI
dat_e_map_ci <- w.getis_e %>%
filter(var1 =="ci_val") %>%
st_as_sf(coords = c("long", "lat"), crs = 4326) %>%
filter(!is.na(LISA)) %>%
dplyr::mutate(lat = sf::st_coordinates(.)[,2],
long = sf::st_coordinates(.)[,1])
pal <- colorFactor("RdBu", dat_e_map_ci$clust_LISA, reverse = T)
dat_e_map_ci %>%
leaflet() %>%
addTiles() %>%
addCircles(lng = ~long, lat = ~lat, color = ~ pal(clust_LISA),
popup = paste("Country", dat_e_map_ci$country, "<br>",
"Cluster:", dat_e_map_ci$cluster_number, "<br>",
"Malaria Prevalence:", dat_e_map_ci$prev, "<br>",
"CI:", dat_e_map_ci$output, "<br>",
"LISA:", dat_e_map_ci$LISA)) %>%
addLegend("bottomright",
pal = pal, values = ~clust_LISA,
title = "LISA - Education CI") %>%
addScaleBar(position = c("bottomleft"))8.3.2 E-SII
dat_e_map_sii <- w.getis_e %>%
filter(var1 =="sii") %>%
st_as_sf(coords = c("long", "lat"), crs = 4326) %>%
filter(!is.na(LISA)) %>%
dplyr::mutate(lat = sf::st_coordinates(.)[,2],
long = sf::st_coordinates(.)[,1])
pal <- colorFactor("RdBu", dat_e_map_sii$clust_LISA, reverse = T)
dat_e_map_sii %>%
leaflet() %>%
addTiles() %>%
addCircles(lng = ~long, lat = ~lat, color = ~ pal(clust_LISA),
popup = paste("Country", dat_e_map_sii$country, "<br>",
"Cluster:", dat_e_map_sii$cluster_number, "<br>",
"Malaria Prevalence:", dat_e_map_sii$prev, "<br>",
"CI:", dat_e_map_sii$output, "<br>",
"LISA:", dat_e_map_sii$LISA)) %>%
addLegend("bottomright",
pal = pal, values = ~clust_LISA,
title = "LISA - Education SII") %>%
addScaleBar(position = c("bottomleft"))8.3.3 E-RII
dat_e_map_rii <- w.getis_e %>%
filter(var1 =="rii") %>%
st_as_sf(coords = c("long", "lat"), crs = 4326) %>%
filter(!is.na(LISA)) %>%
dplyr::mutate(lat = sf::st_coordinates(.)[,2],
long = sf::st_coordinates(.)[,1])
pal <- colorFactor("RdBu", dat_e_map_rii$clust_LISA, reverse = T)
dat_e_map_rii %>%
leaflet() %>%
addTiles() %>%
addCircles(lng = ~long, lat = ~lat, color = ~ pal(clust_LISA),
popup = paste("Country", dat_e_map_rii$country, "<br>",
"Cluster:", dat_e_map_rii$cluster_number, "<br>",
"Malaria Prevalence:", dat_e_map_rii$prev, "<br>",
"CI:", dat_e_map_rii$output, "<br>",
"LISA:", dat_e_map_rii$LISA)) %>%
addLegend("bottomright",
pal = pal, values = ~clust_LISA,
title = "LISA - Education RII") %>%
addScaleBar(position = c("bottomleft"))8.4 Summary Table
dat.ci <- dat_ci_map %>%
st_set_geometry(NULL) %>%
group_by(country) %>%
summarise(#wealth = median(wealth, na.rm = T),
sample_w = sum(sample.size, na.rm = T),
prev_w = median(prev, na.rm = T),
sii_w = median(sii, na.rm = T),
rii_w = median(rii, na.rm = T)
)## `summarise()` ungrouping output (override with `.groups` argument)dat.edu <- dat_edu_map %>%
st_set_geometry(NULL) %>%
group_by(country) %>%
summarise(#wealth = median(wealth, na.rm = T),
sample_e = sum(sample.size, na.rm = T),
prev_e = median(prev, na.rm = T),
sii_e = median(sii, na.rm = T),
rii_e = median(rii, na.rm = T)
)## `summarise()` ungrouping output (override with `.groups` argument)