R script to chart and map Paris tree data

paris_trees.R

library(tidyverse)
library(lubridate)
library(stringr)
library(hrbrthemes)
library(sf)
library(leaflet)

# R script to chart and map Paris tree data

# data available here https://opendata.paris.fr/explore/dataset/les-arbres/ ----------------

trees <- read_csv2("les-arbres.csv") %>% 
  separate(geo_point_2d, c("lat", "lon"), sep = ",") %>% 
  mutate(lat = as.numeric(lat), lon = as.numeric(lon),
         CIRCONFERENCEENCM = as.numeric(CIRCONFERENCEENCM),
         `HAUTEUR (m)` = as.numeric(`HAUTEUR (m)`)) %>% 
  mutate(year = year(DATEPLANTATION),
         century_planted = case_when(
           year < 1700 ~ "17th",
           year < 1800 ~ "18th",
           year < 1900 ~ "19th",
           year < 2000 ~ "20th",
           year < 2018 ~ "21st",
           is.na(year) ~ "Unknown"
         ))

levels <- map_chr(1:20, paste0, "E")
levels[1] <- "1ER"

arr_sums <- trees %>%
  group_by(ARRONDISSEMENT) %>% 
  summarise(Total = n()) %>% 
  filter(str_detect(ARRONDISSEMENT, "PARIS")) %>% 
  mutate(nom = str_replace(ARRONDISSEMENT, "ARRDT", "ARRONDISSEMENT"),
         nom_num = factor(paste0(str_extract(nom, "[[:digit:]]+"), ifelse(str_extract(nom, "[[:digit:]]+") == 1, "ER", "E")), levels = levels)) %>% 
  select(-ARRONDISSEMENT)


# shp data available here ---------------------------------------------------
# https://www.europeandataportal.eu/data/en/dataset/decoupage-administratif-communal-francais-issu-d-openstreetmap/resource/eb4f8187-43b6-4bf9-a3cd-4e2963846f64

boundaries <- st_read("arr_shps/arrondissements-municipaux-20160128.shp") %>% 
  filter(grepl("Paris", nom)) %>%
  mutate(nom = str_to_upper(nom)) %>% 
  left_join(arr_sums, by = "nom") %>% 
  mutate(surf_ha = surf_ha / 100, 
         `Per km2` = round(Total / surf_ha), 0)
         
         
# total trees per arrondissement chart (that looks like trees) -------------------------------

total_trees <- ggplot(boundaries, aes(nom_num, Total)) +
  geom_col(fill = "Sienna", width = .5) +
  geom_text(aes(label=scales::comma(Total), y = Total + 300), vjust=0,  position = position_dodge(0.9), colour = "MediumSeaGreen",
            size = 6, family = "Iosevka") + #nudge_y=500,
  labs(title = "Nombre d'arbres ? Paris par Arrondissement", subtitle = "(Number of trees in Paris by district)", x = "Arrondissement", y = NULL,
       caption = "Chart: @paulcampbell91 Source: Mairie de Paris / Direction des Espaces Verts et de l'Environnement") +
  theme_ipsum(base_family = "Iosevka") + 
  theme(axis.text.y = element_blank(), axis.ticks = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        plot.background = element_rect(fill = "#f5f5f5", colour = "#f5f5f5"), plot.caption = element_text(margin = margin(t = 30)))


# per km2 chart ----------------------------------------------------------------------------------------

persqkm <- ggplot(boundaries, aes(nom_num, `Per km2`)) +
  geom_col(fill = "Sienna", width = .5) +
  geom_text(aes(label=scales::comma(`Per km2`), y = `Per km2` + 50), vjust=0,  position = position_dodge(0.9), colour = "MediumSeaGreen",
            size = 6, family = "Iosevka") +
  labs(title = "Arbres par kilom?tre carr? ? Paris par Arrondissement", subtitle = "(Trees per square kilometre in Paris by district)", x = "Arrondissement", y = NULL,
       caption = "Chart: @paulcampbell91 Source: Mairie de Paris / Direction des Espaces Verts et de l'Environnement") +
  theme_ipsum(base_family = "Iosevka") + #grid = FALSE,
  theme(axis.text.y = element_blank(), axis.ticks = element_blank(), panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        plot.background = element_rect(fill = "#f5f5f5", colour = "#f5f5f5"), plot.caption = element_text(margin = margin(t = 30)))


# ggplot2 map --------------------------------------------------------------------

centroids <- st_centroid(boundaries) %>% 
  mutate(x = sapply(geometry, "[[", 1), y = sapply(geometry, "[[", 2)) %>% 
  select(nom, x, y)

p <- ggplot(boundaries) +
  geom_sf(aes(fill = trees_per_km2)) +
  geom_text(data = centroids, aes(x, y, label = str_extract(nom, "[[:digit:]]+"))) +
  scale_fill_viridis_c()
  
  
# leaflet map --------------------------------------------------------------------

pal <- colorNumeric(palette = "viridis", domain = boundaries$trees_per_km2)
tooltip <- paste0(boundaries$nom, "<br>Trees: ", boundaries$n,
                  "<br>Per Square KM: ", round(boundaries$trees_per_km2, 0)) %>% 
  lapply(htmltools::HTML)

leaflet(boundaries) %>% 
  addProviderTiles(provider = "CartoDB.Positron") %>%
  addPolygons(fill = ~trees_per_km2, fillColor = ~pal(trees_per_km2), fillOpacity = .4,
              label = tooltip, stroke = FALSE, smoothFactor = 0) %>% 
  addLegend("topright", 
            pal = pal, 
            values = ~trees_per_km2,
            title = "Trees per km2",
            opacity = 1)
            
# fin