For an analysis in R I want to create a "new" border from a shapefile. The Shapefile comes with a dummy variable (zone_dummy) indicating the regions laying east and west from the border (see example below). Now I wanted to create a geomatric path along the border of the green and red municipalities:

enter image description here

Here is my shapefile and below an example code:

rm(list = ls(all = TRUE))


area <- readShapePoly("VZ250_GEM.shp")


zone<-data.frame(c("Weiden i.d. OPf.", "Bad Kissingen", "Schweinfurt", "    Mellrichstadt", "Brückenau", "Rhön-Grabfeld",
               "Hofheim", "Ebern",  "Haßberge", "Coburg", "Hof", "Selb", "Kulmbach",
               "Marktredwitz", "Bayreuth", "Bamberg","Staffelstein", "Lichtenfels", "Kronach",
               "Stadtsteinach", "Naila", "Münchberg", "Rehau", "Wunsiedel", "Weiden",
               "Tirschenreuth", "Kemnath", "Neustadt a.d. Waldnaab", "Bad Neustadt a.d. Saale", "Neustadt b. Coburg", "Vohenstrauß", "Nabburg", 
               "Overviechtach", "Waldmünchen", "Cham", "Roding",
               "Deggendorf", "Passau", "Kötzling" ,"Viechtach", "Regen", "Bogen", "Grafenau",
               "Wolfstein", "Wegscheid", "Wunsiedel i. Fichtelgebirge", "Schwandorf", "Freyung-Grafenau"
#create dummy to mark municipalities in zone
#merge with original data
kreise.data.zone<-merge(kreise.data, zone, by="GEN_K", all=T)
## Replace NA's with 0
kreise.data.zone[is.na(kreise.data.zone)] <- 0
##Merge assembled Shape-file Data Frame to Original Shapefile 
area.zone<-merge(area, kreise.data.zone)
#we only care about bavaria:
bayern<-subset(area.zone, area.zone@data$NUTS1_CODE=="DE2")
#plot a map of bavaria with zone: 
my.palette <- brewer.pal(n = 3, name = "Set1")
spplot(bayern, "zone_dummy", cuts = 1 , col.regions = my.palette, col="white", colorkey=FALSE)

My question is now: How can I create a border along the zone (=sharp border marking the path along the red and green municipalities)?

The idea is to (1) create this border in the shapefile (2) estimate the distance of each municipality in Bavaria to this border.


Compute the intersection of both geometries (e.g. sf::st_intersection), here is some R code using the sf package:


zone <- c("Weiden i.d. OPf.", "Bad Kissingen", "Schweinfurt", 
          "Mellrichstadt", "Brückenau", "Rhön-Grabfeld", "Hofheim", "Ebern",
          "Haßberge", "Coburg", "Hof", "Selb", "Kulmbach", "Marktredwitz",
          "Bayreuth", "Bamberg","Staffelstein", "Lichtenfels", "Kronach",
          "Stadtsteinach", "Naila", "Münchberg", "Rehau", "Wunsiedel",
          "Weiden", "Tirschenreuth", "Kemnath", "Neustadt a.d. Waldnaab",
          "Bad Neustadt a.d. Saale", "Neustadt b. Coburg", "Vohenstrauß",
          "Nabburg", "Overviechtach", "Waldmünchen", "Cham", "Roding",
          "Deggendorf", "Passau", "Kötzling" ,"Viechtach", "Regen", "Bogen",
          "Grafenau", "Wolfstein", "Wegscheid", "Wunsiedel i. Fichtelgebirge",
          "Schwandorf", "Freyung-Grafenau")

shp <- st_read("VZ250_GEM.shp") %>%
       filter(NUTS1_CODE=="DE2") %>%
       mutate(zone=factor(if_else(GEN_K %in% zone, 1, 0)))

border <- st_intersection(st_geometry(filter(shp, zone==1)),
                          st_geometry(filter(shp, zone==0)))

ggplot() +
  geom_sf(data=shp, aes(fill=zone)) +
  geom_sf(data=border, colour="white", size=2)



Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.