Making maps in R using ggplot and sf

# Activating packages suppressWarnings(suppressMessages({ library(tidyverse) # Modern data science workflow library(stringr) # Functions designed to work with strings library(sf) # Encodes spatial vector data library(RColorBrewer) # ColorBrewer Palettes library(leaflet) # Interactive web maps library(ggplot2) # Modern visualization tools library(ggthemes) # For extra themes, scales and geoms for 'ggplot2' library(patchwork) # make it simple to combine separate ggplots library(plotly) # Create interactive web graphics from 'ggplot2' library(ggrepel) ## For displaying labels on ggplot2 object })) # Adjustments knitr::opts_chunk$set(fig.align = 'center') options(warn = -1) options(scipen = 10000) options(repr.plot.width = 6, repr.plot.height = 4)

df <- read.csv("https://raw.githubusercontent.com/haginta/regional-convergence-indonesia-satellite-gdp/main/java118_.csv") glimpse(df)

mymap <- read_sf("https://github.com/haginta/regional-convergence-indonesia-satellite-gdp/raw/main/java_118.gpkg") glimpse(mymap)

st_crs(mymap)

df[,4:7] <- log(df[4:7]) head(df)

df <- df %>% mutate(sat_g = sat2020 - sat2019, gdp_g = gdp2020 - gdp2019) head(df)

df_map <- merge(mymap, df, by=c("districtID", "district", "province")) df_map <- st_make_valid(df_map) glimpse(df_map)

df_map %>% ggplot() + geom_sf()

#Create dataframe of selected cities city <- c("Kota Jakarta Pusat", "Kota Yogyakarta", "Kota Surabaya") h_city <- df_map %>% filter(district %in% city) #Add coordinates coordinates <- as.data.frame(st_coordinates(st_centroid(h_city))) h_city$X <- coordinates$X h_city$Y <- coordinates$Y

glimpse(h_city)

ggplot() + geom_sf(data = h_city, aes(fill = district), lwd = 0) + geom_sf(data = df_map, aes(fill = NA, alpha=0), show.legend = F) + geom_sf(data = h_city, aes(fill = district, alpha=0), show.legend = F, color="red") + geom_label_repel(data = h_city, aes(x = X, y = Y, label = district), size=2, color="blue", point.padding = NA) + theme(legend.position="None") + labs(x = "Longitude", y="Latitude")

ggplot(df_map) + geom_sf(aes(fill=sat2019), color = "white", lwd = 0.05) + scale_fill_viridis_c(option = "magma", name = "Radiance (log)") + theme_map() + theme(legend.direction="horizontal") + labs(title = "Districts of Java seen from space at night", subtitle = "Nighttime light (NTL) intensity across 118 districts in Java, Indonesia, 2019", caption = "Geometries: GADM \n Data: Earth Observation Group")

ggplot() + geom_sf(data = h_city, aes(fill = sat2019), color = "white", lwd = 0.05) + geom_sf(data = df_map, aes(fill = sat2019), color = "white", lwd = 0.05, show.legend = T) + geom_sf(data = h_city, aes(fill = sat2019), show.legend = F, color="red") + scale_fill_viridis_c(option = "magma", name = "Radiance (log)") + theme_map() + geom_label_repel(data = h_city, aes(x = X, y = Y, label = district), size=2, color="blue", point.padding = NA) + theme(legend.direction="horizontal", legend.position=c(0.05, -0.12))+ labs(x = "Longitude", y="Latitude")

p1 <- ggplot(df_map) + geom_sf(aes(fill=sat2020), color = "white", lwd = 0.05) + scale_fill_viridis_c(option = "magma", name = "Radiance (log)") + theme_map() + theme(legend.direction="horizontal", legend.position=c(-0.0, -0.15))+ labs(subtitle = "(a) Spatial distribution of NTL per capita") p2 <- ggplot(data = df_map) + geom_sf(aes(fill=gdp2020), color = "white", lwd = 0.05) + scale_fill_viridis_c(option = "magma", name = "GDP (log)") + theme_map() + theme(legend.direction="horizontal", legend.position=c(-0.0, -0.15))+ labs(subtitle = "(b) Spatial distribution of GDP per capita") p1 / p2

p <- ggplot(df_map) + geom_sf(aes(fill=gdp2019, label=district), lwd = 0.05) + scale_fill_viridis_c(option = "magma", name = "GDP (log)") + theme_bw() + labs(x = "Longitude", y="Latitude", title = "GDP across 118 districts in Java, Indonesia, 2019 (constant price 2010)") p_inter <- ggplotly(p) p_inter