Step by step on how to map with R.
Basic knowledge in GIS.
R and Rstudio
For this tuto, R version 4.2.2 Patched (2022-11-10 r83330) is used.
R Package needed for this first tuto is rgdal which provides functions for reading and writing vector and raster spatial data.
FYI, here is the website for all spatial packages available : https://cran.r-project.org/web/views/Spatial.html
PS : In this tuto, we use the Madagascar 23 regions shapefile from this link : https://www.unsalb.org/data/mdg?fbclid=IwAR0bVCV5CJdmKy833ESaY9bbERWo6AbI9Jdv5EasQFpnAkcl_IYoRZ6NeMQ
library(rgdal)
Next, assign the new spatial object within a variable named mada_regions , and use the function readOGR in order to load the shapefile
mada_regions <- readOGR(dsn = "data/mada_23_regions.shp")
str(mada_regions)
sapply(mada_regions@data, class)
plot(mada_regions, col = "lightgrey")
In here, we use the Shape_Area column to check the value.
big_r <- mada_regions$Shape_Area == max(mada_regions$Shape_Area)
plot(mada_regions[big_r, ], col = "orange", add = TRUE)
title("Madagascar, the 23 regions")
To conclude, this is just a simple manipulation to begin with. We will see for the next coming tuto, thematic maps with some spatial analysis!
If you would like to further your skills in mapping with R, there are lots of resources out there, not to mention https://github.com/Robinlovelace/Creating-maps-in-R
In here, we load a csv which contains the statistic data related to number of habitants in Madagascar that I get from the internet (means it is not an official dataset, just for the practice)
stat_pop <- read.csv("pop_region.csv",
stringsAsFactors = FALSE)
Now, we compare the ADM2NM column in mada_regions to region column in stat_pop to see which rows match.
mada_regions$ADM2NM %in% stat_pop$region
The result indicates that all names match
[1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE
[20] TRUE TRUE TRUE TRUE
Join the spatial and non-spatial datasets (the csv) by loading the dplyr package
library(dplyr)
mada_regions@data <- left_join(mada_regions@data, stat_pop, by = c('ADM2NM' = 'region'))
In here, we use tmap package
library(tmap)
tm_shape(mada_regions) + tm_polygons(col = "pop", palette="Reds", title="Number of Population", style="jenks") + tm_layout(bg.color = "lightblue")+tm_scale_bar(breaks = c(0, 100, 200), text.size = 0.7)
- tm_shape takes the input shapefile.
- tm_polygons speficies the column associated with the geometry to get mapped, the color, the title of the legend, and the style for the classification (we use jenks style here)
- tm_layout for background
- tm_scale_bar for the scale bar
- You can check various alternative in here https://geocompr.robinlovelace.net/adv-map.html