read and open shapefile data and plot the ecoregion world map takes too long time

Question

Here is the link to download shapefile data: https://www.facebook.com/flx/warn/?u=https%3A%2F%2Fwww.caribbeanmarineatlas.net%2Fgeoserver%2Fows%3Fservice%3DWFS%26version%3D1.0.0%26request%3DGetFeature%26typename%3Dgeonode%253Atm_ecoregions_2017%26outputFormat%3DSHAPE-ZIP%26srs%3DEPSG%253A4326%26format_options%3Dcharset%253AUTF-8%26fbclid%3DIwAR1pnvtJcBMBuh6Qt4K4JmRNzWbhNDmtO3Je5YtUKK4IwsYHQC7Sk78Tz28&h=AT0d9mHza5n0nmdg9SDxLgozigI5l6Yw2LhgU5pDLDfRAvNlRI3K-3Ga1wTAAh41Fd5VK7NyhrqJQjqyxQaJfQKhgnC3v_kMkFkKdJ_VXohgIULZ0U-blYl9oa918aJkIZEDk9big0dUaBZHgU7n42MdPGTdEuBlbGYJ3nnZdjHDH-dkaWMr0pZE0MzfY5XwlukVHy4g3h6jFESlG3oNPXDnzOb9G6ruqbC-ovRKXcM-rzKJaqJwDaRGRwY9RtG_9h-7gcuq

With the shp. file, firstly I read and open it using readOGR in rgdal, and then fortify it in ggplot2 package to get the long and lat, and then plot it in ggplots. I think the steps are good, but the problem is it runs very slowly. I checked the fortified dataframe and found there is more than 690000 rows... Anyone knows how to aggregate or simplify the dataset?

Mycode is like this:

library(ggplot2)
library(rgdal)
library(ggspatial)
library(raster)
library(cowplot)
library(egg)
library(rnaturalearth)
library(sf)
eco.regions <- readOGR(paste0(getwd(),"/tm_ecoregions_2017/"), "tm_ecoregions_2017")
eco.df <- fortify(eco.regions, region ="BIOME_NAME")


ggplot() +geom_polygon(data = eco.d, aes(x = long, y = lat,group=group, fill=id),
           alpha = 1, linetype = 1)+geom_sf(data = world, fill= NA,    lwd=1,color="white") 

I want to get a world map filled with 14 biomes, but it just takes too long time, and I can't make it finally. Please do help!

Answer 1

First, if you just want to make a nice map then I'd use QGIS. Without much learning, you can load the data in and style it and get this:

QGIS is really quick at this because it takes advantage of spatial indexing and other tricks that R doesn't seem to do.

However, to do efficiently in R, don't use the old sp package and friends, use the sf package:

 eco.regions = st_read("./tm_ecoregions_2017.shp")
 plot(eco.regions[,"BIOME_NAME"])

and after a couple of minutes you get this:

It looks horrible because there's so much outlining. You can change the line style and palette with options. To speed things up you may also dissolve adjacent polygons, generalise the lines to simplify the geometry, or drop the tiniest island features.

I tried doing some of this in R, but there are some topological problems with the map data. Doing a dissolve in QGIS was pretty quick and results in a layer with one feature per biome type, which draws quickly and could be saved to a GeoPackage if you need to do some analysis in R.

I don't doubt this could be done in R using st_union and friends, but sometimes its not the best tool for the job.