4

I'm having kind of a peculiar issue that I'm hoping someone can shed some light on. I'm actually creating my map in R, but I'm using QGIS as a tool to check my mapping process (R newb here). So let's dig into it.

I'm trying to plot rivers in my area of interest (Tanzania), but the source shapefile (diva-gis.org) includes tons of unnecessary water lines. For the sake of reproducability, I'll be including all edited DBFs and shapefiles. Here is the original DBF file and here is the original shapefile. It looks like this. enter image description here

As you can see in the DBF, many of the rivers are unknown/unnamed and listed as "UNK" which isn't very useful for me. In QGIS, I am simply able to open the attribute table, sort by name, and remove all UNK values. Alternatively, I can filter out those values through the SQL expression:

"NAM" != "UNK"

Both of these yield the same results. Here is the QGIS edited DBF and here is the QGIS edited shapefile. This is what the map looks like: enter image description here

Then things get quirky. Again, I'm actually creating this map in R. Using the original source data, I edit the DBF with the "foreign" package, overwrite the DBF, and then map using ggplot. I'm getting a totally different result than the one pictured above. Here's a little snippet of example code (not for the full map, but this is essentially what this portion of the map code will look like):

library(rgdal)
library(foreign)
library(ggplot2)

#set your appropriate working directory
setwd("D:/Mapping-R/Returns-Practice")

#read dbf, sql expression to remove all UNK values in NAM
water.lines <- read.dbf("original-waterlines.dbf")
water.lines <- water.lines[water.lines$NAM != "UNK",]

#overwrite dbf
write.dbf(water.lines, "TZA_water_lines_dcw.dbf")

#read shapefile, fortify to create readable data frame
water.lines.shp <- readOGR(dsn = "D:/Mapping-R/Returns-Practice", 
                  layer = "original-waterlines")
water.lines.shp <- fortify(water.lines.shp)

#plot
ggplot() + geom_polygon(data = water.lines.shp, 
             aes(long, lat, group = group), 
             color = "slategray4")

The result of this code looks like this: enter image description here

Kind of stretched out and funky, but clearly does not follow the same line pattern as the QGIS-edited image. This has been driving me absolutely nuts. I actually went through and checked the R attributes against the QGIS attributes and they matched up perfectly, but the image result is completely different. After playing around with it some more, I figured that my code was to blame and that's where I messed up. But then I opened the original DBF and manually removed all UNK values and plugged it back into QGIS. Here is the resulting DBF and here is its corresponding shapefile, which matches the R image/plot and looks like this: enter image description here

What gives?! The attribute tables are identical but the images are different. I feel like there has to be something really obvious that I'm missing but I've been playing with this for literally hours and can't figure it out.

  • Looks like different projections used for views, geographic and projected in qgis – FelixIP Dec 4 '15 at 19:52
  • 1
    1. You don't need the foreign package. Just use water.lines.shp <- subset(water.lines.sho, NAM != "UNK") after readOGR(). 2. add coord_map() to your plot. – rcs Dec 4 '15 at 20:27
  • You should have received an error if the .shp feature count was not equal to the .dbf record count. I'd consider that a bug in R. – Vince Dec 4 '15 at 20:55
4

Here is how I'd do it in R.

library(rgdal)
# read in the shapefile
tz.wl <- readOGR(".", "TZA_water_lines_dcw")

# now select 
tz.wl.selection <- tz.wl[tz.wl$NAM != 'UNK',]

# you can write this out and compare with your Qgis file:
writeOGR(tz.wl.selection, ".", "TZA_water_lines_selection", driver = "ESRI Shapefile")
  • Thank you, this appeared to work! This is actually where I started out, but I couldn't get past the @data issue. I originally had gotten to the selection point, but instead of using "tz.wl.selection <- tz.wl[tz.wl$NAM != "UNK",]", I used"tz.wl.selection<- tz.wl$NAM[...]". This fixed that issue, and now I have the new shapefile. Thanks again! – Lauren Dec 4 '15 at 21:33
4

Don't mess with Shapefile DBFs!

The connection between geometries and attributes gets messed up easily, especially if you start removing or reordering entries in the DBF. This will never work.

Use a different file format or use a library that supports Shapefiles as a whole and not just DBF.

  • I googled to hell and back how to remove the UNK values for the NAM column but that's the only solution I found! any other suggestions, or should I ask another question for that? – Lauren Dec 4 '15 at 19:54
  • 1
    Sorry, I don't do R. Better open a new question. – underdark Dec 4 '15 at 19:55

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.