5

Using the shapefile available here I am trying two merge the polygons of Sudan and South Sudan, so that I get the borders of Sudan in 2010. My code to make the shapefile available in R is

library(rgdal)
africa    <- readOGR(dsn =  "Data/Shapes", layer = "AfricanCountries")
class(africa)
[1] "SpatialPolygonsDataFrame" attr(,"package") [1] "sp"

I tried different packages and solutions like raster::intersect, rgeos::gIntersect or maptools::unionSpatialPolygons. I always end up with a spatial object that lost the data belonging two both polygons or no union at all.

So far I have the following solution which does not seem to be very handy:

# Seperate Polygons and Data
sp_africa  <- SpatialPolygons(africa@polygons, proj4string = CRS(proj4string(africa)))
dat_africa <- africa@data

# Delete South Sudan from dataframe
dat_africa <- dat_africa[-which(dat_africa$COUNTRY == "South Sudan"),]

# Add row.names to data according to polygon ID'S
rownames(dat_africa) <- dat_africa$OBJECTID

# Get all ID's
allIDs <- africa$OBJECTID

# Get the ID for Sudan, Primary Land (so we only merge the
# Sudan main land (no islands) with the South Sudan main land
sudanID <- africa[which(africa$COUNTRY == "Sudan" & africa$Land_Type == "Primary land"),]$OBJECTID

# Change ID of South Sudan to that of Sudan
allIDs[which(africa$COUNTRY == "South Sudan")] <- sudanID

# Now unite polygons and afterwards merge polygons with data
tmp     <- unionSpatialPolygons(sp_africa, IDs = allIDs)
africa2 <- SpatialPolygonsDataFrame(tmp, data = dat_africa)

If there is an easier, more handy way I would be glad to know about it.

Martin Schmelzer
  • 23,283
  • 6
  • 73
  • 98

1 Answers1

6

you can use aggregate from the raster package.

nsudan <- africa[africa$COUNTRY == "Sudan",]
ssudan <- africa[africa$COUNTRY == "South Sudan",]

plot(nsudan, axes = T, ylim = c(0, 25))
plot(ssudan, add = T)

north and south sudan

sudan <- aggregate(rbind(ssudan, nsudan))
plot(sudan)

sudan

since you create a new SpatialPolygons object, it will be difficult to keep the data from all features. You probably should append the old information from nsudan

# remove Sudan and South Sudan
africa <- africa[!africa$COUNTRY %in% c("Sudan", "South Sudan"),]

# Adjust the Polygon IDs
africa <- spChFIDs(africa, as.character(africa$OBJECTID))
sudan <- spChFIDs(sudan, as.character(max(africa$OBJECTID) + 1))

library(maptools)
africaNew <- spRbind(sudan, africa)
plot(africaNew)

africa new

loki
  • 9,816
  • 7
  • 56
  • 82
  • In the end, I want to have the same SpatialPolygonsDataFrame just without the border between Sudan and South Sudan. Would that work with aggregate too? – Martin Schmelzer Nov 16 '16 at 13:08
  • I also get the error `Error in match.fun(FUN) : argument "FUN" missing` when trying to aggregate both countries. – Martin Schmelzer Nov 16 '16 at 13:12
  • Not sure why, but `aggregate` is not available in the `sp` package allthough it is inside the documentation: `Error: object ‘aggregate’ is not exported by 'namespace:sp'`. I am using version 1.2-3. – Martin Schmelzer Nov 16 '16 at 13:17
  • I edited my question in that sense, that I added the solution I have so far. It works but it does not seem to be the most handy approach. I will see if I can make it any shorter by using `aggregate`. – Martin Schmelzer Nov 16 '16 at 14:41
  • Egypt seems to be missing in your final map. Also, the border between Sudan and South Sudan is still there. – Martin Schmelzer Nov 17 '16 at 09:15
  • Your code now still has the problem of non-unique ID's ;) – Martin Schmelzer Nov 17 '16 at 10:27