Is it possible to exclude areas of maps made with sf and tmap?

Question

I have created a map of Europe using a shape file downloaded from Eurostat. The file also includes several European administered areas such as French Guiana. I would like to exclude these areas. Is that possible?

Here is a picture of the map, and here is my r code:

require(readxl)
library(sf)
library(tmap)
library(tmaptools)

options(scipen = 999)

data1 <- read_excel("National Data.1-kopi.xlsx")
mydata <- subset(data1, TIME == 2007)


mymap <- st_read("NUTS_RG_20M_2021_3035.shp", stringsAsFactors = FALSE)
str(mymap)

map_and_data <- inner_join(mymap, mydata)

tm_shape(map_and_data) +
  tm_polygons("POP25", id = "NUTS_ID", palette = "Blues")

Sure, I think, just fishing around here `which(my_nuts$LEVL_CODE == 0)` vs `plot(my_nuts$geometry[which(my_nuts$LEVL_CODE == 0)[1:12]])`, you want to knock out guyana, curocao, etc, you can also do your join to data, then set smaller bounding box to plot, assess...so playing around with leaving out LEVL_CODE(s) and you'll likely end up with something like `==0)[c(1:8, 14:23,26)]`... as a first approach. — Chris, Mar 23 '22 at 13:41
To arrive at a bounding box useful to your purposes, as proposed below. — Chris, Mar 23 '22 at 14:18

score 3 · Answer 1 · answered Mar 23 '22 at 13:45

Yes, you can. There is an option in tm_shape called bbox that is the one you should use to specify the limits of your map.

One important thing to remember (the tricky part) is that the box should be in the same unit system than the shape. If your shape is in longitude/latitude, use degrees, but in your case (crs=EPSG:3035), that is quite common for Europe, I use c(2377294, 1313597, 7453440, 5628510):

library(sf)
library(tmap)
library(tmaptools)
# Use giscoR, same source of map
library(giscoR)

mymap <- gisco_get_nuts(
  year = 2021,
  resolution = 20, epsg = 3035, nuts_level = 0
)


# Lets just pretend this is your data
set.seed(1241)
POP25 <- sample(1:50000, nrow(mymap))
mydata <- data.frame(
  NUTS_ID = mymap$NUTS_ID,
  POP25 = POP25
)

map_and_data <- inner_join(mymap, mydata)

# You need to pass the bbox argument. For Europe EPSG:3035
# this work
box <- c(2377294, 1313597, 7453440, 5628510)
class(box) <- "bbox"


tm_shape(map_and_data, bbox = box) +
  tm_polygons("POP25", id = "NUTS_ID", palette = "Blues")

score 3 · Answer 2 · answered Mar 23 '22 at 13:45

I feel your pain, as I was facing exactly the same issue when preparing data for a blog post https://www.jla-data.net/eng/breadbaskets-of-europe/

At the end I resolved the issue in a hacky way - I filtered out problematic areas exceeding a certain distance from Bern (which makes a nice center for "civilized Europe" - I was happy to leave the bulk of Russia out, but needed the part around Königsberg in).

So consider this piece of code:

library(sf)
library(tmap)
library(dplyr)

# download shapefiles from GISCO API
europe <- giscoR::gisco_get_countries(resolution = "10",
                                      region = "Europe")


# Overseas France, Azores + Canaries, and don't get me started on East Prussia!
troublemakers <- subset(europe, CNTR_ID %in% c("FR", "PT", "ES", "RU"))

buffer <- tidygeocoder::geo("Bern", quiet = T) %>% 
  st_as_sf(coords = c("long", "lat"), crs = 4326) %>% 
  st_buffer(units::as_units(1750, "km")) %>% 
  st_geometry() # don't need the address column anymore

# remove the bits too far from Bern (roughly midway between Cadiz and Königsberg)
troublemakers <- troublemakers %>% 
  st_intersection(buffer) 

# administrative regions in equal area projection
europe <- europe %>% 
  filter(!CNTR_ID %in% c("RU", "SJ", "FR", "PT", "ES")) %>% 
  rbind(troublemakers) %>% 
  st_transform(3035) # equal area Lambert / "European Albers"

# apply it in a {tmap} fashion...
tm_shape(europe) +
  tm_polygons()

I liked 'bread baskets' very much when I saw it a while ago. — Chris, Mar 23 '22 at 14:21
@Chris thanks! I wrote the post in late 2021, not fully realizing how relevant the issue may turn in light of the current situation in Ukraine (a _major_ wheat producing area of Europe) — Jindra Lacko, Mar 23 '22 at 15:09
Yes, after invasion I was looking at Ukraine's many exports, including pharma clinical trial data, saw your's, noting the 'consulting', and if you have a chance to look at [things go wrong when you don't understand `data.table` sufficiently and attempt a forced marriage with `sf`](https://stackoverflow.com/questions/71413790/r-st-sample-same-multipolygon-for-different-point-classes), as that guy seems in search of a solution. — Chris, Mar 23 '22 at 15:16

Is it possible to exclude areas of maps made with sf and tmap?

2 Answers2