leaflet_tutorial.Rmd

---
title: "Leaflet Maps"
author: "Insert Name"
date: "27/11/2021"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning = FALSE, message = FALSE)
```

# Package Installation

```{r}
# Use install.packages("package name") to install any of the packages below
# To install the rgeoboundaries package you will need to run the following two lines of code
# install.packages("remotes")
# remotes::install_gitlab("dickoa/rgeoboundaries")

library(tidyverse) ## For plotting and data wrangling.

library(leaflet) ## For leaflet interactive maps

library(sf) ## For spatial data

library(RColorBrewer) ## For colour palettes

library(htmltools) ## For html

library(leafsync) ## For placing plots side by side

library(kableExtra) ## Table output

library(stringr) ## String manipulation

library(rgeoboundaries) ## Administrative boundaries
```


# Birth Rates in North Carolina

* We will use the North Carolina (`nc`) data from the `sf` package.

* Let's load in the data for North Carolina using the function `st_read`.

```{r nc-data, eval = TRUE, messages = FALSE}

nc_df <- st_read(system.file("shape/nc.shp", package="sf"))

```

* `st` = spatial type and `.shp` is a common shape file format (e.g. GIS). 

## Rename the columns and view the data North Carolina

* Number of births for counties in North Carolina in 1974 

* Rename our columns to country, births, and geometry.


```{r nc-data2, eval = TRUE, messages = FALSE}

nc <- nc_df %>%
        select("NAME", "BIR74", "BIR79", "geometry") %>%
        rename("county" = "NAME", "births1974" = "BIR74", "births1979" = "BIR79")

```


## Let's inspect the data

* Let's load in the data for North Carolina (nc)

```{r nc-view, eval = TRUE}

head(nc)


```

## And inspect the structure


```{r nc-str, eval = TRUE}

str(nc)

```

## Building a map in ggplot2

```{r first-map1a, eval=TRUE}
ggplot(nc) #<<

```

```{r first-map1e, eval=TRUE}
ggplot(nc) +
  geom_sf(aes(fill = births1974)) +
  labs(title = "Births per county in 1974",
       x = "Longitude",
       y = "Latitude", 
       fill = "Births") +
  scale_y_continuous(breaks = 34:36) #<<
```

## Accessing Administrative Boundaries

```{r access admin 0 boundaries, fig.alt="Map of Nigeria"}

# One country
nigeria <- geoboundaries(country = "Nigeria") #<<

ggplot(data = nigeria) +
  geom_sf()
```


```{r access admin 0 boundareis for nigeria and chad, fig.alt="Map of Nigeria and Chad: admin level 0."}

# Two countries
nigeria_chad <- geoboundaries(country = c("Nigeria", "Chad")) #<<

ggplot(data = nigeria_chad) +
  geom_sf()
```

```{r access admin level 1 boundaries, fig.alt="Map of Nigeria: admin level 1."}

# Admin Level 1
nigeria <- geoboundaries(country = "Nigeria", 
                         adm_lvl = "adm1") #<<

ggplot(data = nigeria) +
  geom_sf()
```

`r icons::icon_style(icons::fontawesome("lightbulb", style = "solid"), scale = 1)` Find out more administrative boundaries and other open data sources tutorials at [rspatialdata.github.io](https://rspatialdata.github.io)


```{r View Nigeria data}
nigeria %>%
  kableExtra::kable()
```

### Joining spreadsheet data to spatial data


* Let's say we have population data we want to match. 

```{r nigeria-pop}

# Two countries
nigeria_pop <- read_csv(file = "data/nga_admpop_adm1_2020.csv") #<<

head(nigeria_pop)
```

`r icons::icon_style(icons::fontawesome("lightbulb", style = "solid"), scale = 2)` In what column do we find spatial information we can use to match this data to the other dataset?

* We can match the data based on the names of the administrative boundaries.

```{r nigeria-pop-spatial-join}

# Joining Nigeria spatial data to the population data by the admin name columns

nigeria_pop <- nigeria_pop %>%
  mutate(ADM1_NAME = str_to_title(ADM1_NAME, locale = "en"))

# how well do they match up?
anti_join(x = nigeria, 
          y = nigeria_pop,
          by = c("shapeName" = "ADM1_NAME"))
```


Often one of the trickiest part is to get the spelling of regions to match up. We can check where mismatches occur with an `anti_join`. The NA here is from Akwa Lbom and Akwa lbom. 

* We can fix mismatches using `case_when`

```{r nigeria-pop-spatial-join-2}

# Joining Nigeria spatial data to the population data by the administrative data column

nigeria_pop <- nigeria_pop %>%
  mutate(ADM1_NAME = str_to_title(ADM1_NAME, locale = "en"),
         ADM1_NAME_SC = case_when(ADM1_NAME == "Akwa Ibom" ~ "Akwa lbom",
                               TRUE ~ ADM1_NAME))

nigeria_df <- left_join(x = nigeria, 
                        y = nigeria_pop,
                        by = c("shapeName" = "ADM1_NAME_SC"))

```

`r icons::icon_style(icons::fontawesome("info-circle", style = "solid"), scale = 1)` We can use `case_when` to change the name to match. 

* Now we can create our map

```{r first-nigeria-map3, fig.alt="Map of Nigeria with admin 1 areas filled by population"}

ggplot(data = nigeria_df) +
  geom_sf(aes(fill = T_TL)) +
  labs(x = "Longitude",
       y = "Latitude",
       fill = "Total Population")

```


## Coordinate reference system

* Every location on earth is specified by a longitude and latitude. 

* The Coordinate Reference system (CRS) determines how the data will be projected onto a map. 

* We can check the CRS using `st_crs`:

```{r checking the CRS}

st_crs(nc)

```

* The CRS is specified in the attributes `epsg` and `proj4string`. 


```{r checking the CRS Nigeria, eval = FALSE}

### Fill in the blank to check the nigeria data set

st_crs(____)

```


## Transforming coordinate reference system


* You can transform a coordinate reference system using the `st_transform()`. 

* But what is a sensible coordinate reference system to assign? 

* Well, a good place to start is with the one that leaflet uses for plotting the world:  EPSG 4326.

```{r changing the CRS using st_transform}

nc <- st_transform(nc, "+init=epsg:4326")

st_crs(nc)

```


```{r changing the CRS using st_transform nigeria, eval = FALSE}

## Fill in the blanks for the nigeria dataset
nigeria_df <- st_transform(nigeria_df, "+init=epsg:4326")

______(nc)

```


# Our first leaflet map

* Every plot starts with `leaflet()`

```{r first-leaflet-map1a, eval=TRUE}

leaflet(data = nc) 

```


* Layers are added using `%>%`

```{r first-leaflet-map1b, eval=TRUE}
leaflet(data = nc) %>%
  addTiles() #<<
```


* N.B. Layers are added with `%>%` in `leaflet` and `+` in `ggplot`. `%>%` also is used in the `tidyverse` packages. 


* We can set the view using `setView()`

```{r first-leaflet-map1c, eval=TRUE}
leaflet(data = nc) %>%
  addTiles() %>%
  setView(lng = -80, #<<
          lat = 34.5, #<<
          zoom = 5) #<<
```

* Add different background map using addProviderTiles

```{r first-leaflet-map1d, eval=TRUE}
leaflet(data = nc) %>%
  addProviderTiles(providers$Esri.OceanBasemap) %>%
  setView(lng = -80, 
          lat = 34.5, 
          zoom = 5)
```

* N.B. The different provider tiles come with different licensing.

* Add polygons using `addPolygons()`

```{r first-leaflet-map1e, eval=TRUE}
leaflet(data = nc) %>%
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = -80, 
          lat = 34.5, 
          zoom = 5) %>%
  addPolygons() #<<
```

* Fill in the blanks to create a leaflet map with the Nigeria data. What longitude and latitude do you want to set the view with?


```{r nigeria-first-leaflet, eval = FALSE}

leaflet(data = ) %>%
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = _____, 
          lat = _____, 
          zoom = 5) %>%
  addPolygons() #<<

```


## Creating a colour palette

* What type of data are we showing? 

* Who is our audience? 

* Tools

  * [ColorBrewer2.org](https://colorbrewer2.org/#type=sequential&scheme=BuGn&n=3)
  * [Color Picker](https://www.colorcodepicker.com/)


## Creating a colour palette

* `RColorBrewer` includes **sequential** colour palettes (e.g. number of people).

```{r pal-sequential, eval = TRUE}
display.brewer.all(type = "seq")
```

## Creating a colour palette

* `RColorBrewer` includes **diverging** colour palettes (e.g. to show distinct categories).

```{r pal-diverging, eval = TRUE}
display.brewer.all(type = "div")
```


* First we will define the colour palette and bins for the plot.

```{r Defining colour bins, eval = TRUE}

summary(nc$births1974)
summary(nc$births1979)

bins <- seq(0, 35000, 5000) 
```


* Then we can define the colours for the palette:


```{r Defining the palette, eval = TRUE}
pal74 <- colorBin("OrRd", domain = nc$births1974, bins = bins) 

pal79 <- colorBin("OrRd", domain = nc$births1979, bins = bins)
```


* What might make sense as bins for showing our Nigeria population?

```{r bins for nigeria plot, eval = FALSE}

summary(nigeria_df$TL)

bins <- seq(from = ___, to = ____, by = ____)


pal <- colorBin("OrRd", domain = ________, bins = bins) 

```



* Customising `addPolygons()`

```{r first-leaflet-map1f, eval=TRUE}
leaflet(data = nc) %>%
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = -80, 
          lat = 34.5, 
          zoom = 6) %>%
  addPolygons(
    fillColor = ~pal74(nc$births1974),
    fillOpacity = 0.7, #<<
    color = "white", #<<
    opacity = 1, #<<
    weight = 2 #<<
  ) 
```




```{r nigeria leaflet 2, eval = FALSE}

# Try further customizing your Nigeria plot

leaflet(data = ) %>%
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = _____, 
          lat = _____, 
          zoom = 5) %>%
  addPolygons(
    fillColor = ~pal(_______),
    fillOpacity = 0.7, #<<
    color = "white", #<<
    opacity = 1, #<<
    weight = 2 #<<
  ) #<<

```


* Customising `addPolygons()`

```{r first-leaflet-map1g, eval=TRUE}
leaflet(data = nc) %>%
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = -80, 
          lat = 34.5, 
          zoom = 6) %>%
  addPolygons(
    fillColor = ~pal74(nc$births1974), #<<
    fillOpacity = 1, #<<
    color = "blue", #<<
    opacity = 0.7, #<<
    weight = 1 #<<
  ) 
```

# What can you customise in addPolygons()

```{r what can you customise with addPolygons, eval = FALSE}

?addPolygons()

```

* `color:` stroke color
* `weight:` stroke width in pixels
* `opacity:` stroke opacity
* `fillColor:` fill color
* `fillOpacity:` fill opacity
* `highlightOptions:` Options for highlighting the shape on mouse over. 

* Let's assign our plot to an object.

```{r first-leaflet-map1h, eval=TRUE}
m1 <- leaflet(data = nc) %>% #<<
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = -80, 
          lat = 34.5, 
          zoom = 6)

m1 %>%
  addPolygons(
    fillColor = ~pal74(nc$births1974), 
    fillOpacity = 0.7,
    opacity = 1,
    color = "white", 
    weight = 2) 


```

* Let's add some `highlightOptions`

```{r first-leaflet-map1j, eval=TRUE}
m1 %>%
  addPolygons(
      fillColor = ~pal74(nc$births1974), 
      fillOpacity = 0.7, 
      color = "white", 
      opacity = 1,
      weight = 2,
    highlight = highlightOptions( #<<
        weight = 3, #<<
        color = "blue", #<<
        fillOpacity = 1, #<<
        bringToFront = TRUE)) #<< 
```

## Let's add some labels!

`sprintf`: returns a character vector containing a formatted combination of text and variable values.


```{r Make our labels}
labels <- sprintf("<strong>%s</strong><br/>%g births", 
                  nc$county, nc$births1974) %>% lapply(htmltools::HTML)

head(labels, 1)
```

html - markup language for the web
* `<strong>` = bold; `<br/>` = new line

PHP - Hypertext Preprocessor
* `%s` = place holder for a character string; `%g` = general format place holder for a number

```{r first-leaflet-map1k, eval=TRUE}
(m1 <- m1 %>%
  addPolygons(data = nc,
      fillColor = ~pal74(nc$births1974),
      fillOpacity = 0.7,
      color = "white",
      opacity = 1,
      weight = 2,
      highlight = highlightOptions(
        weight = 3,
        color = "blue",
        fillOpacity = 1,
        bringToFront = TRUE),
      label = labels)) #<< 

```

* Let's add a legend

```{r first-leaflet-map1l, eval=TRUE}
m1 <- m1 %>%
  addLegend( #<<
    position = "bottomright", #<<
    pal = pal74, #<<
    values = ~nc$births1974, #<<
    title = "Births by county in 1974", #<<
    opacity = 1) #<<

m1
```

## Let's create a second map

* Let's create a second map of births in 1979.

* First we'll need to create a new set of labels

```{r labels for second map}

labels79 <- sprintf(
  "<strong>%s</strong><br/>%g births",
  nc$county, nc$births1979
) %>% lapply(htmltools::HTML)


```

## Let's create a second map

```{r second map}

m2 <- leaflet(data = nc) %>%
  addProviderTiles(providers$Stamen.Terrain) %>%
  setView(lng = -80, lat = 34.5, zoom = 6) %>%
  addPolygons(
      fillColor = ~pal79(nc$births1979),
      fillOpacity = 0.7,
      color = "white",
      opacity = 1,
      weight = 2,
      highlight = highlightOptions(
        weight = 3,
        color = "blue",
        fillOpacity = 1,
        bringToFront = TRUE),
      label = labels79)

```

```{r leaflet-map2a, eval=TRUE}
(m2 <- m2 %>%
  addLegend(
    position = "bottomright",
    pal = pal79,
    values = ~nc$births1979,
    title = "Births by country in 1979",
    opacity = 1))

```

## Placing two maps side by side

```{r leaflet-map3}

leafsync::sync(m1, m2, ncol = 2, sync = "all")
```

# Leaflet maps with points

```{r leaflet-map-points, echo = TRUE}

work <- data.frame(
  "location" = c("Valparaíso, Chile", 
                 "Curitiba, Brasil", 
                 "Sable Island, Nova Scotia", 
                 "Greifswald, Germany",
                 "Arusha, Tanzania",
                 "Kigali, Rwanda",
                 "Kingston, Jamaica",
                 "Asunción, Paraguay",
                 "East Kilbride, Scotland"), 
  "institute" = c("Instituto de Fomento Pesquero", 
                  "Universidade Federal do Paraná", 
                  "Dalhousie University", 
                  "Friedrich Loeffler Institut",
                  "Nelson Mandela African Institute of Science and Technology",
                  "National Institute of Statistics Rwanda",
                  "Caribbean National Statistical Offices",
                  "El Ministerio de Salud Pública y Bienestar Social y el Ministerio de Educación y Ciencias - Paraguay",
                  "Foreign Commonwealth Development Office"), 
  "work" = c("Chilean Pink Cusk Eel", 
             "Fox rabies", 
             "Grey seals", 
             "Fox rabies",
             "Teaching",
             "Teaching",
             "Teaching",
             "Teaching",
             "Teaching"), 
  "lat" = c(-33.0472, 
            -25.4290, 
            43.9337, 
            54.0865,
            -3.3995,
            -1.9415,
            18.0179,
            -25.2637,
            55.760869), 
  "lon" = c(-71.6127, 
            -49.2671, 
            -59.9149, 
            13.3923,
            36.7968,
            30.0574,
            -76.8099,
            -57.5759,
            -4.22407), 
  "icon" = c("fish", 
             "disease", 
             "gps", 
             "disease",
             "training",
             "training",
             "training",
             "training", 
             "training"))
```

```{r data-preview} 

work %>%
  kbl() %>%
  kable_paper("hover")

```


```{r work-leaflet, eval=TRUE}
leaflet(work) %>% 
      addProviderTiles(providers$Stamen.Watercolor) %>%
      addProviderTiles(providers$Stamen.TerrainLabels) %>%
      addCircleMarkers(~lon, ~lat)
```

* **Your Turn*** What else can you change about addCircleMarkers? Hint: Type `??addControl` Try adding: `clusterOptions = markerClusterOptions()` to your map. 

* Add labels

```{r point-labels}


labels <- sprintf(
  "<strong>%s</strong>",
  work$institute) %>% lapply(htmltools::HTML)

```

```{r}
(work_map <- leaflet(work) %>% 
                addProviderTiles(providers$Stamen.Watercolor) %>%
                addProviderTiles(providers$Stamen.TerrainLabels) %>%
                addCircleMarkers(~lon, ~lat, popup = ~labels))
```