It's easy to run into projection issues when you using spatial data from R packages without paying attention. In this quick post I will outline a few easy solutions for fixing common projection issues.
dggridr
Let's start with the dggridR package, which generates discrete global grids as SpatialPolygons
or data frames. Here I'm generating an Icosahedron Snyder Equal Area (ISEA) grid of resolution 5. I'm using a Robinson projection, but EPSG:4326
will have the same issues.
library(dggridR)
dggs <- dgconstruct(res = 5)
grid <- dgearthgrid(dggs, frame = FALSE) %>%
st_as_sf()
ggplot() +
geom_sf(data = grid) +
coord_sf(crs = "ESRI:54030") +
theme_void()
There's clearly an issue with polygons crossing the date line. Using an orthographic projection we can take a better look at the polygons causing the artefacts on our map.
ggplot() +
geom_sf(data = grid) +
geom_sf(data = grid[1153,], color = "#00ccff", fill = "#ccf5ff") +
geom_sf(data = st_sfc(st_linestring(rbind(c(180, 0), c(180, 90))), crs = 4326), color = "#cc3300") +
coord_sf(crs = "+proj=ortho +lat_0=90 +lon_0=0") +
theme_void()
To fix this issue we can use st_wrap_dateline(options = c("WRAPDATELINE=YES", "DATELINEOFFSET=90"))
, which will split polygons that cross the date line. Note the DATELINEOFFSET
option which has value 10 by default, and will determine the range around the date line within which polygons will be split. Let's reproduce the previous map, but now with the st_wrap_dateline()
applied to the highlighted grid cell:
ggplot() +
geom_sf(data = grid) +
geom_sf(data = st_sfc(st_linestring(rbind(c(180, 0), c(180, 90))), crs = 4326), color = "#cc3300") +
geom_sf(data = grid[1153,] %>% st_wrap_dateline(options = c("WRAPDATELINE=YES", "DATELINEOFFSET=90")), color = "#00ccff", fill = "#ccf5ff") +
coord_sf(crs = "+proj=ortho +lat_0=90 +lon_0=0") +
theme_void()
After applying this to the entire grid, our Robinson projection map will look a whole lot better:
grid_fixed <- grid %>%
st_wrap_dateline(options = c("WRAPDATELINE=YES", "DATELINEOFFSET=270"))
ggplot() +
geom_sf(data = grid_fixed) +
coord_sf(crs = "ESRI:54030") +
theme_void()
We could opt to also use st_segmentize()
to make the map a little more accurate (i.e. with slightly curved cell boundaries instead of straight ones), but that would make the rendering considerably slower.
ggplot2::map_data() and rnaturalearth
The ggplot2
package has a function map_data()
which wraps the maps
package and provides shapes of countries and states ready for use in ggplot2
maps. However, if you play close attention you will see that there's an issue with Antarctica when using the WGS84 coordinate reference system (see bottom left).
world <- map_data(map = "world") %>%
filter(region == "Antarctica")
ggplot() +
geom_polygon(data = world, aes(x = long, y = lat, group = group), color ="#595959", fill ="#e5e5e5") +
coord_sf(crs = 4326) +
theme_void()
An easy solution is to just switch to the rnaturalearth package:
world <- ne_countries(scale = "large", returnclass = "sf") %>%
filter(admin == "Antarctica")
ggplot() +
geom_sf(data = world) +
coord_sf() +
theme_void()
However, in contrast to the ggplot2::map_data()
Antarctica shape, the one from rnaturalearth
does not work well with an orthographic projection:
ggplot() +
geom_sf(data = world) +
coord_sf(crs = "+proj=ortho +lat_0=-90 +lon_0=0") +
theme_void()
The solution here is to apply a minimal buffer using st_buffer()
after transforming the shape to the orthographic projection:
world_fixed <- world %>%
st_transform("+proj=ortho +lat_0=-90 +lon_0=0") %>%
st_buffer(1)
ggplot() +
geom_sf(data = world_fixed) +
coord_sf(crs = "+proj=ortho +lat_0=-90 +lon_0=0") +
theme_void()