1 Background of project

With the steady growth in spatial analysis with R, there have been attempts to improve map visualization comparable or better than traditional graphical user interface (GUI) based geographical information systems (GIS) platforms such as QGIS.In that process, several packages have been developed in R including tmap, Mapdeck, sf among others. The list has not left out ggplot2 that is the bread-and-butter of data visualization in R. When sf package was built for wrangling spatial data in a manner similar to regular data frames in R, the spatial science community appreciated the analytic capabilities but were not fully satisfied with visualization of the outputs. This is the niche that has been filled by the development of mapsf package see. In this short tutorial, we are going to see how to use mapsf to visualize several vector map types, at some point I shall make one for rasters and combination of the two data types/structures. You can find the code generating the html from my page on GitHub.

2 About mapsf and installation

mapsf is a package for easier creation of maps of simple features (sf) in R. For more information on the package, see.

You can install and load the mapsf package by un-commenting the following codes:

# install.packages('mapsf') or the development version using remotes::install_github("riatelab/mapsf")
# library(mapsf)

Since it is mapping simple features, loading the library will show Loading required package: sf and of course link to GEOS, GDAL, and PROJ.

3 Loading data to map with mapsf

The package is coming with an inbuilt data-set for Martinique Island. We can access the data-set and save it to an object called study_area using the following simple code:

study_area <- mf_get_mtq() # Most of the functions in mapsf package start with mf_

3.1 Initializing map development

We will now initialize the mapping process. We can think of this as designing the layout pane on which to display the map in QGIS, for example. In this case, we use mf_init() function as follows:

mf_init(study_area, expandBB = rep(0, 4), theme = 'iceberg') # This is the layout to throw map onto

# You can check other possible themes using ?mf_theme

3.2 Previewing the data to map

Before doing the mapping itself, it is important to preview the data-set so that we understand the available attributes (fields) that can be mapped. We can achieve that by using the head() function just like in the case of normal data.frames in R.

head(study_area) # multipolygon, WGS84 projection zone 20N. Six features and seven fields plus bounding box. geom holds the coordinates and not regarded as a field in the seven fields.
## Simple feature collection with 6 features and 7 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: 695444 ymin: 1598818 xmax: 717731 ymax: 1645182
## Projected CRS: WGS 84 / UTM zone 20N
##   INSEE_COM              STATUS            LIBGEO  POP   MED CHOM  ACT
## 1     97201 Simple municipality L'Ajoupa-Bouillon 1902 13633  254  801
## 2     97202 Simple municipality Les Anses-d'Arlet 3737 14558  425 1659
## 3     97203 Simple municipality      Basse-Pointe 3357 14456  400 1395
## 4     97204 Simple municipality         Le Carbet 3683 18847  285 1568
## 5     97205 Simple municipality       Case-Pilote 4458 21005  347 2096
## 6     97206 Simple municipality        Le Diamant 5976 19704  430 2654
##                             geom
## 1 MULTIPOLYGON (((699261 1637...
## 2 MULTIPOLYGON (((709840 1599...
## 3 MULTIPOLYGON (((697602 1638...
## 4 MULTIPOLYGON (((702229 1628...
## 5 MULTIPOLYGON (((698805 1621...
## 6 MULTIPOLYGON (((709840 1599...

We can as well go ahead and get the column names plus data types of all the attributes in the sf object.

glimpse(study_area)
## Rows: 34
## Columns: 8
## $ INSEE_COM <chr> "97201", "97202", "97203", "97204", "97205", "97206", "97207~
## $ STATUS    <chr> "Simple municipality", "Simple municipality", "Simple munici~
## $ LIBGEO    <chr> "L'Ajoupa-Bouillon", "Les Anses-d'Arlet", "Basse-Pointe", "L~
## $ POP       <dbl> 1902, 3737, 3357, 3683, 4458, 5976, 17792, 790, 82502, 17540~
## $ MED       <dbl> 13633, 14558, 14456, 18847, 21005, 19704, 19341, 15151, 1692~
## $ CHOM      <dbl> 254, 425, 400, 285, 347, 430, 1546, 84, 10046, 2342, 78, 142~
## $ ACT       <dbl> 801, 1659, 1395, 1568, 2096, 2654, 8352, 343, 37615, 8178, 2~
## $ geom      <MULTIPOLYGON [m]> MULTIPOLYGON (((699261 1637..., MULTIPOLYGON ((~

In this case, one must be understanding the data fully to make sense of what the abbreviated column names mean for example POP, MED, CHOM, and ACT. However, the geom column, that is sometimes called geometry, is a unique attribute which stores a list of all the coordinates of the features in the data. When dropped, we end up with a normal data.frame object in R.

4 Creating map shadow

There is a shadow feature in mapsf which can be plotted just before placing the actual map on top. Actually, it is the shadow of the actual map. We can go ahead and plot the shadow using mf_shadow() as follows:

mf_init(study_area, 
        expandBB = rep(0, 4), 
        theme = 'dark')
mf_shadow(study_area, 
          col = 'purple', 
          cex = 2, 
          add = TRUE) # Now we have shadow on the layout added

4.1 Adding map on the shadow

The next step now is to plot the map itself on top of its shadow. Sounds like an overlay on top of the same layer. However, the shadow will be slightly offset diagonally to depict ‘shadowy’ plot.

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE) 
mf_map(study_area, type = 'base', add = TRUE) # Awesome

5 The map layout development

In as much as the map layout looks fine, it is still way below what we can achieve with mapsf. Let us use the mf_layout() function to create a layout more suitable for map creation. The layout will capture north arrow, scale, Source, author, package and version used, and title of the map.

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_layout(title = 'Martinique',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))

6 Adding one of the fields/attributes to the map using proportional circles

Now we want to show a map of the population since POP is a field in the study_area data per region. We want to represent this as graduated symbol so that we have larger symbols for those regions with higher population values.

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_layout(title = 'Martinique',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(study_area, add = TRUE)
mf_map(x = study_area, var = 'POP', type = 'prop', inches = 0.25, 
       col = 'purple', leg_pos = 'topright', leg_title = 'Total population'
)

7 Adding choropleth map using one attribute

Awesome so far. Here we make a step to the choropleth maps. Since regions are polygons, we can calculate area of each polygon using the st_area() function from sf package. Then we can use the population column and the created area values to create a new column called population density. This is given by dividing the population by area.

study_area$POPDENS <- 1e6 * study_area$POP / st_area(study_area)
mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_layout(title = 'Martinique',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(x = study_area, var = "POPDENS", type = "choro", breaks = "geom", nbreaks = 5,
  col = "Mint", border = "white", lwd = 0.5, leg_pos = "topright", 
  leg_title = "Population Density\n(people per km2)", add = TRUE
)

8 Mapping categorical variable on the map

We can also plot a map which is showing the character variable in the fields such as the STATUS field. Then we color the map based on individual categories, so all regions having same category will have the same color. Here I changed the shadow color to yellow because magenta has been used in mapping and may bring confusion.

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'yellow', cex = 2, add = TRUE)
mf_layout(title = 'Martinique',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(x = study_area, var = "STATUS", type = "typo", pal = c("purple", "magenta", "cyan"), lwd = .5, val_order = c("Prefecture", "Sub-prefecture", 
"Simple municipality"), leg_pos = "topright", leg_title = "", add = TRUE)

9 Labeling some regions of the map

We may need to label some of the regions. This is achievable by using the mf_label() function as follows:

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_layout(title = 'Martinique',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(x = study_area, var = "STATUS", type = "typo", pal = c("purple", "magenta", "cyan"), lwd = .5, val_order = c("Prefecture", "Sub-prefecture", 
"Simple municipality"), leg_pos = "topright", leg_title = "", add = TRUE)
mf_label(x = study_area[study_area$STATUS != "Simple municipality", ], var = "LIBGEO", cex = 0.9, halo = TRUE, r = 0.15)

10 Combining choropleth with proportional map

Both the proportional map can be combined with choropleth coloration. The size of circles and the color density represent different variables. This is achievable with the following code:

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_layout(title = 'Martinique',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(study_area, add = TRUE)
mf_map(x = study_area, var = c("POP", "MED"), type = "prop_choro", border = "grey50",
  lwd = 1, leg_pos = c("topright", "right"), leg_title = c("Population","Median\nIncome\n(in euros)"), breaks = "equal", nbreaks = 4, 
  pal = "Greens", leg_val_rnd = c(0, -2), leg_frame = c(TRUE, TRUE)
)

11 Making proportional map and typology in one

Proportional map and typology can be mapped on the same map layout using the following code:

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_layout(title = 'Martinique with proportional and typology',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(study_area, add = TRUE)
mf_map(x = study_area, var = c("POP","STATUS"), type = "prop_typo", symbol = "square",
  border = "white", lwd = .5, leg_pos = c("right",  "topright"), 
  leg_title = c("Population", "Administrative\nStatus"),
  val_order = c("Prefecture", "Sub-prefecture", "Simple municipality"))

Other functions such as mf_arrow, mf_scale,and mf_credits can be used to alter the position of the map aesthetics.

11.1 Adding graticule on the map

It may sometimes, not always, be necessary to add graticule (longitudes and latitudes) on the map. That can be achieved simply with the following:

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')

mf_theme(mar = c(2, 2, 1.2, 0)+5)
plot(st_geometry(study_area), border = NA, col = NA, graticule = st_crs(4326),
     axes = TRUE, bg = mf_theme(mar = c(2, 2, 1.2, 0)+2)$bg,
     lon = seq(-62, -60, by = 0.2), lat = seq(14, 15, by = 0.2))
mf_layout(title = 'Martinique with graticules',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_shadow(study_area, col = 'magenta', cex = 1, add = TRUE)
mf_map(study_area, add = TRUE)

12 Adding inset map

Again we may need to show a map as an inset from a larger map. For example, a district within a country. Here we indicate inset from the whole globe. The point of the map from the globe is marked with red triangle.

mf_init(study_area, expandBB = rep(0, 4), theme = 'dark')
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_layout(title = 'Martinique with inset shown',
          credits = paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_map(x = study_area, var = "POPDENS", type = "choro", breaks = "geom", nbreaks = 5,
  col = "Mint", border = "white", lwd = 0.5, leg_pos = "left", 
  leg_title = "Population Density\n(people per km2)", add = TRUE)
mf_inset_on(x = 'worldmap', pos = 'topright')
mf_worldmap(study_area, col = "red")
mf_inset_off()

13 Fun with background

If you have an interesting image which you want to make appear in the background of your map, then you can add it with mf_background function.

mf_init(study_area, theme = 'dark')
mf_background(system.file("img/background.jpg", package = "mapsf"))
mf_shadow(study_area, col = 'magenta', cex = 2, add = TRUE)
mf_map(x = study_area, var = "POPDENS", type = "choro", breaks = "geom", nbreaks = 5,
  col = "Mint", border = "white", lwd = 0.5, leg_pos = "left", 
  leg_title = "Population Density\n(people per km2)", add = TRUE)
mf_credits(paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_title("Martinique map with a background")
mf_inset_on(x = 'worldmap', pos = 'topright')
mf_worldmap(study_area, col = "red")
mf_inset_off()

14 Exporting a map out to png for publication

mf_export(x = study_area, filename = "mtq.png", width = 600, 
          theme = "dark", expandBB = c(0,0,0,.3)) 
mf_shadow(study_area, col = "magenta", add = TRUE)
mf_map(x = study_area, var = "MED", type = "choro",
       pal = "Dark Mint", breaks = "quantile", nbreaks = 6, 
       leg_title = "Median Income\n(euros)", leg_val_rnd = -2, add = TRUE)
mf_inset_on(x = "worldmap", pos = "right")
mf_worldmap(study_area, col = "red")
mf_inset_off()
mf_title("Wealth in Martinique, 2015")
mf_credits(paste0('Sources: IGN, 2018\nAuthor: Wyclife Agumba Oluoch\n',
          'Package: mapsf ',
          packageVersion('mapsf')))
mf_scale(size = 5)
mf_arrow('topleft')
dev.off()
## png 
##   2

You will find the exported image in your working directory. Here is my exported image The Exported Map.

References

Allaire, JJ, Yihui Xie, Jonathan McPherson, Javier Luraschi, Kevin Ushey, Aron Atkins, Hadley Wickham, Joe Cheng, Winston Chang, and Richard Iannone. 2021. Rmarkdown: Dynamic Documents for r. https://CRAN.R-project.org/package=rmarkdown.
Giraud, Timothée. 2021. Mapsf: Thematic Cartography. https://CRAN.R-project.org/package=mapsf.
Pebesma, Edzer. 2018. “Simple Features for R: Standardized Support for Spatial Vector Data.” The R Journal 10 (1): 439–46. https://doi.org/10.32614/RJ-2018-009.
———. 2021. Sf: Simple Features for r. https://CRAN.R-project.org/package=sf.
R Core Team. 2021. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing. https://www.R-project.org/.
Wickham, Hadley, Romain François, Lionel Henry, and Kirill Müller. 2021. Dplyr: A Grammar of Data Manipulation. https://CRAN.R-project.org/package=dplyr.
Xie, Yihui. 2014. “Knitr: A Comprehensive Tool for Reproducible Research in R.” In Implementing Reproducible Computational Research, edited by Victoria Stodden, Friedrich Leisch, and Roger D. Peng. Chapman; Hall/CRC. http://www.crcpress.com/product/isbn/9781466561595.
———. 2015. Dynamic Documents with R and Knitr. 2nd ed. Boca Raton, Florida: Chapman; Hall/CRC. https://yihui.org/knitr/.
———. 2021. Knitr: A General-Purpose Package for Dynamic Report Generation in r. https://yihui.org/knitr/.
Xie, Yihui, J. J. Allaire, and Garrett Grolemund. 2018. R Markdown: The Definitive Guide. Boca Raton, Florida: Chapman; Hall/CRC. https://bookdown.org/yihui/rmarkdown.
Xie, Yihui, Christophe Dervieux, and Emily Riederer. 2020. R Markdown Cookbook. Boca Raton, Florida: Chapman; Hall/CRC. https://bookdown.org/yihui/rmarkdown-cookbook.