Motivation

Reference: R Charts - Chart 1 - Chart 2

Total rainfall at some monitoring stations 2022

Level: Vietnam Province

Download Total rainfall (mm) in 2022.xlsx

Bubble map

# Clear workspace: 
rm(list = ls())

# Setwd
setwd("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Rain_22")

library(giscoR)
library(dplyr)
library(sf)
library(ggplot2)

# Load tidyverse package: 

library(tidyverse)

# Vietnam Map (include Spratly Islands and Paracel Islands) : 

link <- "https://data.opendevelopmentmekong.net/dataset/999c96d8-fae0-4b82-9a2b-e481f6f50e12/resource/2818c2c5-e9c3-440b-a9b8-3029d7298065/download/diaphantinhenglish.geojson?fbclid=IwAR1coUVLkuEoJRsgaH81q6ocz1nVeGBirqpKRBN8WWxXQIJREUL1buFi1eE"

vn_spatial <- sf::st_read(link)

ggplot() + 
  geom_sf(data = vn_spatial)

# Import rain data

library(rio)
rain_22 <- import("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Rain_22/Rain_22.xlsx")

# Combine the two datasets: 
rain_map_df <- right_join(vn_spatial, rain_22, by = c ("Name" = "Province"))

# Replace blank values with 0
rain_map_df <- rain_map_df %>%
  mutate(Year = ifelse(is.na(Year), 0, Year))

# Countries centroids (Multi to Point)
rain_point_df <- st_centroid(rain_map_df, of_largest_polygon = TRUE)

# Mapping

library(giscoR)
library(dplyr)
library(sf)
library(ggplot2)

ggplot() +
  geom_sf(data = vn_spatial, fill = "white") +
  geom_sf(data = rain_map_df, fill = "grey90")+
  geom_sf(
    data = rain_point_df,
    pch = 21,
    aes(size = Year, fill = Year),
    col = "grey90", alpha = 0.9)+
  # Adjust scale size and legend of scale
  scale_size(
    range = c(1, 9),
    guide = guide_legend(
      direction = "vertical",
      nrow = 4,
      label.position = "right"))+
  guides(fill = guide_legend(title = "")) +
  # Add title, caption
  labs(title = "Total rainfall at some monitoring stations 2022",
       caption = "Unit: mm\nSource: GSO",
       size = "") +
  theme_void() +
  theme(legend.position = c(0.9,0.7))+
  # Adjust background
  theme(plot.background = element_rect(fill = "#f0f0f0", color = "#f0f0f0")) +
  theme(panel.background = element_rect(fill = "#f0f0f0", color = "#f0f0f0")) +
  # Adjust margin
  theme(plot.margin = unit(c(0.1, 0.5, 0.1, 0.5), "cm"))+
  # Adjust title, caption, legend
  theme(plot.title = element_text(size = 11, face = "bold", color = "grey30"))+
  theme(plot.caption = element_text(size = 6, face = "italic", color = "grey50"))+
  theme(legend.text = element_text(size = 8, color = "grey50"))

  ggsave("rain_22.png", width = 3.8, height = 4.5,dpi = 300,units = c("in"))

Choropleth map

# Clear workspace: 
rm(list = ls())

# Setwd
setwd("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Rain_22")

library(giscoR)
library(dplyr)
library(sf)
library(ggplot2)

# Load tidyverse package: 

library(tidyverse)

# Vietnam Map (include Spratly Islands and Paracel Islands) : 

link <- "https://data.opendevelopmentmekong.net/dataset/999c96d8-fae0-4b82-9a2b-e481f6f50e12/resource/2818c2c5-e9c3-440b-a9b8-3029d7298065/download/diaphantinhenglish.geojson?fbclid=IwAR1coUVLkuEoJRsgaH81q6ocz1nVeGBirqpKRBN8WWxXQIJREUL1buFi1eE"

vn_spatial <- sf::st_read(link)

ggplot() + 
  geom_sf(data = vn_spatial)

# Import rain data

library(rio)
rain_22 <- import("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Rain_22/Rain_22.xlsx")

# Combine the two datasets: 
rain_map_df <- right_join(vn_spatial, rain_22, by = c ("Name" = "Province"))

# Replace blank values with 0
rain_map_df <- rain_map_df %>%
  mutate(Year = ifelse(is.na(Year), 0, Year))

# Countries centroids (Multi to Point)
rain_point_df <- st_centroid(rain_map_df, of_largest_polygon = TRUE)

# Mapping

library(giscoR)
library(dplyr)
library(sf)
library(ggplot2)

ggplot() +
    geom_sf(data = vn_spatial, fill = "grey95") +
    geom_sf(data = rain_map_df, aes(fill = Year))+
    labs(title = "Total rainfall at some monitoring stations 2022",
         caption = "Unit: mm   Source: GSO")+
    # Adjust theme
    theme_minimal()+
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank())+
    # Adjust fill color and legend format
    scale_fill_gradientn(colours = hcl.colors(5, "YlGnBu", rev = TRUE, alpha = 0.9),
                         guide = guide_colourbar(direction = "horizontal",
                                                 barheight = unit(2, units = "mm"),
                                                 barwidth = unit(30, units = "mm"),
                                                 title.hjust = 0.5,
                                                 label.hjust = 0.5, 
                                                 title.position = "top"))+
    # Adjust legend position
    theme(legend.position = c(0.8,0.85))+
    # Adjust margin
    theme(plot.margin = unit(c(0.5, 0.5, 0.5, 0.5), "cm"))+
    # Adjust background
    theme(plot.background = element_rect(fill = "white", color = "white")) +
    theme(panel.background = element_rect(fill = "white", color = "white")) +
    # Adjust title, sub, caption
    theme(plot.title = element_text(color = "grey20", size = 30, face = "bold")) + 
    theme(plot.subtitle = element_text(size = 22, color = "gray40")) + 
    theme(plot.caption = element_text(size = 18, colour = "grey40", face = "italic")) + 
    theme(legend.text = element_text(color = "grey40", size = 18)) + 
    theme(legend.title = element_text(color = "grey20", size = 18))
  
  
  ggsave("rain_22.1.png", width = 3.8, height = 4.5,dpi = 300,units = c("in"))

Forest areas in Mekong Delta Vietnam 2022

Level: Provinces in the Mekong Delta

Download Forest areas in Mekong Delta Vietnam.xlsx

Bubble map

# Clear workspace: 
rm(list = ls())

# Setwd
setwd("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Forest_22")

# Bubble Map

# Get geospatial data for Viet Nam: 

library(raster)
vietnam_province <- getData("GADM", country = "Vietnam", level = 1)

detach(package:raster)
library(tidyverse)
vietnam_df <- vietnam_province %>% fortify(region = "NAME_1")

library(stringi)

vietnam_df <- vietnam_df %>% 
  mutate(id_prov = stri_trans_general(id, "Latin-ASCII")) %>% 
  mutate(id_prov = case_when(str_detect(id_prov, "Ba Ria") ~ "BRVT", 
                             str_detect(id_prov, "Ho Chi Minh") ~ "TP.HCM", 
                             str_detect(id_prov, "Thua Thien Hue") ~ "TT-Hue", 
                             TRUE ~ id_prov))
# Get forest data

  library(rio)
  forest_22 <- import("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Forest_22/Forest_22.xlsx")
  
  # https://www.gps-latitude-longitude.com/address-to-longitude-latitude-gps-coordinates
  
  area_forest_22_b <- forest_22 %>% 
    select("Province", "Area", "Lat","Long")
  
  area_forest_22_b <- area_forest_22_b %>% 
    mutate(id_prov = stri_trans_general(Province, "Latin-ASCII"))
  
# Join data frame
  
  mekong_df <- vietnam_df %>%
    semi_join(area_forest_22 %>% select(id_prov), by = c("id_prov" = "id_prov"))
  
  mekong_df <- left_join(mekong_df, area_forest_22, by = "id_prov")
  
  library(summarytools)
  freq(mekong_df$id_prov)
  
  library(showtext)
  showtext.auto() 
  my_font <- "Roboto Condensed"
  font_add_google(name = my_font, family = my_font)
  
# Mapping
  
  ggplot() + 
    geom_polygon(data = mekong_df, aes(long, lat, group = group), fill = "grey90", color = "white", linewidth = 0.1) +
    coord_map("albers", lat0 = 30, lat1 = 40) +
    geom_point(data = area_forest_22_b, aes(Long, Lat, size = Area, color = Area))+
    # Adjust scale size and legend of scale
    scale_size(
      range = c(3, 10),
      guide = guide_legend(
        direction = "vertical",
        nrow = 1,
        label.position = "bottom"))+
    # Adjust scale color
    scale_color_gradientn(colours = hcl.colors(5, "ag_GrnYl", rev = TRUE, alpha = 0.9)) +
    guides(color = guide_legend(title = ""))+
    # Add title, sub, and caption
    labs(title = "Forest areas (thousands of hectares) in Mekong Delta Vietnam",
         subtitle = "Total: 246,700 hec    Year: 2022",
         caption = "Data Source: GSO",
         size = "")+
    theme_void()+
    # Adjust background
    theme(plot.background = element_rect(fill = "white", color = "white")) +
    theme(panel.background = element_rect(fill = "white", color = "white")) +
    # Legend position
    theme(legend.position = c(0.8,0.2))+
    # Adjust title, sub, caption
    theme(plot.title = element_text(family = my_font, color = "grey20", size = 30, face = "bold")) + 
    theme(plot.subtitle = element_text(family = my_font, size = 22, color = "gray40")) + 
    theme(plot.caption = element_text(family = my_font, size = 18, colour = "grey40", face = "italic"))+
    theme(legend.text = element_text(family = my_font, color = "grey40", size = 18)) + 
    # Adjust plot margin
    theme(plot.margin = unit(c(0.2, 0.5, 0.2, 0.2), "cm"))
  
  ggsave("forest_22b.png", width = 3.8, height = 3,dpi = 300,units = c("in"))

Choropleth map

Option 1: Viridis package

# Clear workspace: 
rm(list = ls())

# Setwd
setwd("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Forest_22")

# Choropleth Map

# Get geospatial data for Viet Nam: 

library(raster)
vietnam_province <- getData("GADM", country = "Vietnam", level = 1)

detach(package:raster)
library(tidyverse)
vietnam_df <- vietnam_province %>% fortify(region = "NAME_1")

library(stringi)

vietnam_df <- vietnam_df %>% 
  mutate(id_prov = stri_trans_general(id, "Latin-ASCII")) %>% 
  mutate(id_prov = case_when(str_detect(id_prov, "Ba Ria") ~ "BRVT", 
                             str_detect(id_prov, "Ho Chi Minh") ~ "TP.HCM", 
                             str_detect(id_prov, "Thua Thien Hue") ~ "TT-Hue", 
                             TRUE ~ id_prov))
# Get forest data

library(rio)
forest_22 <- import("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Forest_22/Forest_22.xlsx")

# https://www.gps-latitude-longitude.com/address-to-longitude-latitude-gps-coordinates

names(forest_22) <- c("Province", "Area", "Natural", "Planted", "Percent", "Lat", "Long")

area_forest_22 <- forest_22 %>% 
  select("Province", "Area")

area_forest_22 <- area_forest_22 %>% 
  mutate(id_prov = stri_trans_general(Province, "Latin-ASCII"))

# Join data frame

mekong_df <- vietnam_df %>%
  semi_join(area_forest_22 %>% select(id_prov), by = c("id_prov" = "id_prov"))

mekong_df <- left_join(mekong_df, area_forest_22, by = "id_prov")

library(summarytools)
freq(mekong_df$id_prov)

library(showtext)
showtext.auto() 
my_font <- "Roboto Condensed"
font_add_google(name = my_font, family = my_font)

# Mapping

library(viridis)

ggplot() + 
  geom_polygon(data = mekong_df, aes(long, lat, group = group, fill = Area), color = "white", linewidth = 0.1) +
  coord_map("albers", lat0 = 30, lat1 = 40)+
  labs(title = "Forest areas (thousands of hectares) in Mekong Delta Vietnam",
       subtitle = "Total: 246,700 hec    Year: 2022",
       caption = "Data Source: GSO") +
  theme(axis.line = element_blank(),
        axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks = element_blank(),
        axis.title.x = element_blank(),
        axis.title.y = element_blank(),
        panel.grid = element_blank(),
        plot.background = element_rect(fill = "white"),
        panel.background = element_rect(fill = "white"),
        legend.background = element_rect(fill = NA),
        panel.border = element_blank())+
  theme(plot.title = element_text(family = my_font, color = "grey20", size = 30, face = "bold")) + 
  theme(plot.subtitle = element_text(family = my_font, size = 22, color = "gray40")) + 
  theme(plot.caption = element_text(family = my_font, size = 18, colour = "grey40", face = "italic")) + 
  theme(legend.text = element_text(family = my_font, color = "grey40", size = 18)) + 
  theme(legend.title = element_text(family = my_font, color = "grey20", size = 18)) + 
  theme(plot.margin = unit(c(0.2, 0.5, 0.2, 0.2), "cm"))+
  theme(legend.position = c(0.8,0.2)) + 
  theme(plot.caption = element_text(hjust = 0))+
  scale_fill_viridis(direction = 1, 
                     option = "D", 
                     name = "Area", 
                     guide = guide_colourbar(direction = "horizontal",
                                             barheight = unit(2, units = "mm"),
                                             barwidth = unit(30, units = "mm"),
                                             title.hjust = 0.5,
                                             label.hjust = 0.5, 
                                             title.position = "top"),
                     na.value = "grey90")

  ggsave("forest_22.png", width = 3.8, height = 3,dpi = 300,units = c("in"))

Option 2

# Clear workspace: 
rm(list = ls())

# Setwd
setwd("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Forest_22")

# Choropleth Map

# Get geospatial data for Viet Nam: 

library(raster)
vietnam_province <- getData("GADM", country = "Vietnam", level = 1)

detach(package:raster)
library(tidyverse)
vietnam_df <- vietnam_province %>% fortify(region = "NAME_1")

library(stringi)

vietnam_df <- vietnam_df %>% 
  mutate(id_prov = stri_trans_general(id, "Latin-ASCII")) %>% 
  mutate(id_prov = case_when(str_detect(id_prov, "Ba Ria") ~ "BRVT", 
                             str_detect(id_prov, "Ho Chi Minh") ~ "TP.HCM", 
                             str_detect(id_prov, "Thua Thien Hue") ~ "TT-Hue", 
                             TRUE ~ id_prov))
# Get forest data

library(rio)
forest_22 <- import("D:/0 - My documents/TOOLS/R/The Economist/Bubble Map/Forest_22/Forest_22.xlsx")

# https://www.gps-latitude-longitude.com/address-to-longitude-latitude-gps-coordinates

names(forest_22) <- c("Province", "Area", "Natural", "Planted", "Percent", "Lat", "Long")

area_forest_22 <- forest_22 %>% 
  select("Province", "Area")

area_forest_22 <- area_forest_22 %>% 
  mutate(id_prov = stri_trans_general(Province, "Latin-ASCII"))

# Join data frame

mekong_df <- vietnam_df %>%
  semi_join(area_forest_22 %>% select(id_prov), by = c("id_prov" = "id_prov"))

mekong_df <- left_join(mekong_df, area_forest_22, by = "id_prov")

library(summarytools)
freq(mekong_df$id_prov)

library(showtext)
showtext.auto() 
my_font <- "Roboto Condensed"
font_add_google(name = my_font, family = my_font)

# Mapping

ggplot() + 
  geom_polygon(data = mekong_df, aes(long, lat, group = group, fill = Area), color = "white", linewidth = 0.1) +
  coord_map("albers", lat0 = 30, lat1 = 40)+
  labs(title = "Forest areas (thousands of hectares) in Mekong Delta Vietnam",
       subtitle = "Total: 246,700 hec    Year: 2022",
       caption = "Data Source: GSO")+
  # Adjust theme
  theme_minimal()+
  theme(axis.title.x = element_blank(),
        axis.title.y = element_blank())+
  # Adjust fill color and legend format
  scale_fill_gradientn(colours = hcl.colors(5, "ag_GrnYl", rev = TRUE, alpha = 0.9),
                       guide = guide_colourbar(direction = "horizontal",
                                               barheight = unit(2, units = "mm"),
                                               barwidth = unit(30, units = "mm"),
                                               title.hjust = 0.5,
                                               label.hjust = 0.5, 
                                               title.position = "top"),
                       na.value = "grey90")+
  # Adjust legend position
  theme(legend.position = c(0.8,0.2))+
  # Adjust plot margin
  theme(plot.margin = unit(c(0.2, 0.5, 0.2, 0.2), "cm"))+
  # Adjust background
  theme(plot.background = element_rect(fill = "white", color = "white")) +
  theme(panel.background = element_rect(fill = "white", color = "white")) +
  # Adjust title, sub, caption
  theme(plot.title = element_text(family = my_font, color = "grey20", size = 30, face = "bold")) + 
  theme(plot.subtitle = element_text(family = my_font, size = 22, color = "gray40")) + 
  theme(plot.caption = element_text(family = my_font, size = 18, colour = "grey40", face = "italic")) + 
  theme(legend.text = element_text(family = my_font, color = "grey40", size = 18)) + 
  theme(legend.title = element_text(family = my_font, color = "grey20", size = 18))

ggsave("forest_22.1.png", width = 3.8, height = 3,dpi = 300,units = c("in"))