Shapefile Mapping in R – Solving Partial Plot Issues with ggplot2

Question

I would like to plot shapefiles from the EUFORGEN website (http://www.euforgen.org/) in R (European beech and Scots pine). At first I had problems, because they are files with holes. I was able to solve this problem for beech using sebdalgarno’s hint.

For pine it is a bit more difficult, because with the same method, suddenly only single fragments (isolated population) are displayed but not the core population (see picture).

I could imagine that I could "adjust" something when reading in the file, so that the whole layer is displayed.

Any ideas?

Edit (referring to the comment from @Spacedman)

I'm sorry, but I don't know how to read a shape-file via r-script from the internet. If someone knows how to do that, I can change that with pleasure. Therefore here again the link to the cloud

Cloud_Beech(Fagus sylvatica) and Pine(Pinus sylvestris)

or to the zip-file from the Webpage
EUFORGEN.org

Here is the map as it is almost correct. The darker small areas (e.g., in Norway and Sweden) are holes in the polygon, but are not displayed correctly with geom_polygon().

library("rnaturalearth")
library("rgdal")
library("ggplot2")


gg_world <- ne_countries(scale = "medium", returnclass = "sp")  # world map as sp class

Pinus <- readOGR("./chorological_maps_dataset/Pinus sylvestris/shapefiles/Pinus_sylvestris_plg_clip.shp")


ggplot()+
  geom_polygon(data = gg_world, mapping=aes(x = long, y = lat, group = group), fill = "white", color = "grey50")+
  # Species shape file
  geom_polygon(data = Pinus, aes(x = long, y = lat, group = group), colour = "black", fill = "darkblue", alpha = 0.2)+
  coord_cartesian(xlim=c(-15, 40), ylim=c(35, 75), default = F)+
  # Layout
  theme_bw()+
  # axes
  scale_x_continuous(expand = c(0, 0), name = "Longitude (°E)")+
  scale_y_continuous(expand = c(0, 0), name = "Latitude (°N)")

So I continued with sebdalgarno's suggestion and used geom_sf() and as also above, with the hint from Spacedman to use a world map in sf format. I used ne_countries from the rnaturalearth packages

gg_world <- ne_countries(scale = "medium", returnclass = "sf")
Pinus_st <- st_read("./chorological_maps_dataset/Pinus sylvestris/shapefiles/Pinus_sylvestris_plg_clip.shp")


ggplot() + 
 
  
  geom_sf(data = Pinus_st, colour = NA, fill = "darkblue", alpha = 0.5)+
  #geom_sf(data = Fagus_st, colour = NA, fill = "red", alpha = 0.5)+
  geom_sf(data = gg_world, fill = NA, color = "grey50")+
  
  coord_sf(xlim=c(-15, 40), ylim=c(35, 75), expand = F)+
  
  # Layout
  theme_bw()+
  
  # axes
  scale_x_continuous(expand = c(0, 0), name = "Longitude (°E)")+
  scale_y_continuous(expand = c(0, 0), name = "Latitude (°N)")

and get this

A new step-by-step build of the code then showed that the map section boundary was causing the error. As was also the case here Link. If one takes out the alpha command, it works.

gg_world <- ne_countries(scale = "medium", returnclass = "sf")
Pinus_st <- st_read("./chorological_maps_dataset/Pinus sylvestris/shapefiles/Pinus_sylvestris_plg_clip.shp")


ggplot() + 
 
  
  geom_sf(data = Pinus_st, colour = NA, fill = "darkblue")+
  geom_sf(data = gg_world, fill = NA, color = "grey50")+
  coord_sf(xlim=c(-15, 40), ylim=c(35, 75), expand = F)+
  # Layout
  theme_bw()+
  # axes
  scale_x_continuous(expand = c(0, 0), name = "Longitude (°E)")+
  scale_y_continuous(expand = c(0, 0), name = "Latitude (°N)")

Answer 1

Something is buggy somewhere. Here's how I'd map that shapefile with tmap:

library(sf)
library(tmap)
library(tmaptools)

data(World)

Pinus = st_read("./chorological_maps_dataset/Pinus sylvestris/shapefiles/Pinus_sylvestris_plg_clip.shp")

xlim=c(-15, 40)
ylim=c(35, 75)


map = tm_shape(World, bbox=bb(xlim=xlim, ylim=ylim)) +
    tm_graticules() + 
    tm_polygons("white") +
    tm_shape(Pinus) + tm_polygons("darkblue", alpha=0.2) 

print(map) # or just "map" if interactive

Adjust colours, labels etc as per requirements, but this is close enough:

I see the holes in Scandinavia which I guess were the real problems.

Or with ggplot, using the same Pinus read as before, as an sf object:

ggplot() + geom_sf(data=World) + geom_sf(data = Pinus, fill="darkblue", alpha=0.2) +
    coord_sf(xlim=xlim, ylim=ylim)

I prefer tmap because its easy to change the mode and make an interactive browser map, but I give the ggplot map as well anyway.

And for completion here's the base graphics version:

plot(st_geometry(World), xlim=xlim, ylim=ylim)
plot(Pinus, add=TRUE, col="#00008040")