Today I will illustrate how to create an interactive map. Imagine you have a dataset where each observation is a country or area. It could be useful for communication to present a map summarizing the main findings in each paper you have reviewed. Instead of a table, it can be more interesting to your audience to see a map.
This application can be quite useful to create maps summarizing the studies considered in a meta analysis, for example. Not only do we illustrate the spatial spread of our studies, we can also give some more detailed information about each study in such a way that the audience can interact with the map.
In order to fiddle with maps, there are a few packages that need to be installed beforehand. These are: dplyr, tidyverse, leaflet, rgeos, rworldmap and leafpop.
install.packages("dplyr", "tidyverse", "leaflet", "rgeos", "rworldmap", "leafpop")
After installing the packages, remember to use the library function:
library(dplyr) library(tidyverse) library(leaflet) library(rgeos) library(rworldmap) library(leafpop)
We will use publicly available data from De Groot et al. (2012). The authors reviewed hundreds of studies that used different valuation methods to estimate the value of a variety of ecosystem services. In the Supplementary material provided by the authors, they share several tables that summarize all the studies and values obtained.
Let us map all the valuation studies that focused on the value of coral reefs as reviewed by De Groot et al. (2012) (see Table 2.2 in the Supplementary Material). Coral reefs are unique ecosystems that provide a variety of ecosystem services. I copy pasted the corresponding table into an excel file, which you can download here: Groot studies.
I converted the excel file into a csv file and then imported it to R:
coral_reefs <- read.csv("Groot studies.csv", sep=";")
I can start looking into my data. First, I will convert the Value vector into a numeric vector and then use the summarize function to see the descriptive statistics.
> coral_reefs$Value <- as.numeric(coral_reefs$Value..2007...ha.yr.) > summary(coral_reefs$Value) Min. 1st Qu. Median Mean 3rd Qu. Max. 1.00 19.25 42.50 41.39 62.75 84.00
De Groot et al. (2012) find a total of 94 estimates of the value of coral reefs. As reported above, the mean value of coral reefs is 47.39 $ per hectare and per year, ranging from 1$ to 84$.
I also want to take a look at which valuation methods the reviewed papers have used.
> table(coral_reefs$Valuation.Method) AC CV DMP FI / PF GV Other PES RC TC 6 12 62 3 1 1 1 4 4
These are the meanings of the acronyms for the valuation methods:
- AC – Avoided Cost;
- CV – Contingent Valuation;
- DMP – Direct market pricing;
- FI /PF – Factor Income / Production Function;
- GV – Group Valuation;
- PES – Payment for Ecosystem services;
- RC – Replacement cost;
- TC – Travel Cost.
The majority of the papers included in their review used the direct market pricing method, and a minority used the contingent valuation method.
The data as it is, however, does not yet include coordinates to match each observation to a point in the map. In order to get coordinates for each country, I will use the code described here. To implement this code, the packages rgeos and rworldmap are needed (we already installed them before). This is the part of the code that we need:
wmap <- getMap(resolution="high") centroids <- gCentroid(wmap, byid=TRUE) df <- as.data.frame(centroids)
The code above produces the df dataframe. This dataframe now records the latitude and longitude for each country. In order to look up the coordinates of the centroid of the United States, I can type:
> df["United States of America",1:2] x y United States of America -112.4617 45.67955
where x is the longitude and y is the latitude.
However, we need to merge the df and coral_reefs dataframes. I want the coordinates from the df dataframe to be matched with the country in the coral_reefs dataframe. I will use the merge function, as follows:
df$Country <- row.names(df) df.coral = merge(df[, c("Country", "x", "y")], coral_reefs[, c("Country", "Value", "Valuation.Method", "Reference")])
Unfortunately not all countries could be merged because the countrys name in the df dataframe do not match the country names in the coral_reefs dataframe. The United States is an example. Since this is just an illustration, I will ignore this problem in this blog post. If you want country names to match 100%, you should change the names manually, either in the df or in the coral_reefs dataframe. We are left with 74 matches, compared with the 94 observations in the coral_reefs dataframe.
We now have all we need to generate our map. This tutorial helped me to get started on my map. I then personalized it by, for example, changing the basemap to providers$CartoDB.Positron. I also used the popupTable function to define which information shows up when I click on a marker. You can copy and paste the code below to reproduce the map. If you are using your own data to generate your map, make sure to replace “df.coral” and variable names with the ones from your own data.
df.coral$longitude % addProviderTiles(providers$Esri.WorldImagery, group = "World Imagery") %>% addProviderTiles(providers$CartoDB.Positron) %>% addLayersControl(baseGroups = c("Toner Lite", "World Imagery")) %>% addMarkers(label = df.coral$Value, popup = popupTable(df.coral, zcol = c("Value", "Valuation.Method", "Reference"))) %>% addMiniMap( toggleDisplay = TRUE, tiles = providers$Stamen.TonerLite )
The code above generates the following map. Click on it to be redirected to an interactive map.
If you click on the blue markers, the paper reference, value (per hectare and year) and the valuation method will appear. I defined which variables should appear in the “zcol” option when writing the code for the map (see above)
This type of map can also be a neat way of summarizing one’s work (e.g. articles and case studies you have focused on) or even international experience (which places you have been affiliated with).
De Groot, R., Brander, L., Van Der Ploeg, S., Costanza, R., Bernard, F., Braat, L., … & Hussain, S. (2012). Global estimates of the value of ecosystems and their services in monetary units. Ecosystem services, 1(1), 50-61.