5

I am trying to create a shapefile from a .csv file in R. I can successfully create a working shapefile with this code, but when I try to view the attribute table for it in arcMap 10.4, the table only contains columns for objectID, latitude, and longitude. I need the info in the other columns for symbolizing the data correctly.

Is there a way to specify which columns of the data to keep in the shapefile?

(I know how to make a data frame into a shapefile within arcMap, but I would like to automate this step in R if possible.)

library(maptools)
library(rgdal)
library(sp)

#data
site <- c("a","b","c","d")
prop_c <- c(0.88,0.48,0.15,0.47)
prop_b <- c(0.17,0.18,0.09,0.08)
minus_c <- 1-prop_c
minus_b <- 1-prop_b
lat <- c(44.22,38.38,33.35,43.48)
long <- c(-124.45, -123.70, -124.40, -124.05)
MyData <- cbind.data.frame(site, prop_c, prop_b, minus_c, minus_b, lat, long)

#convert data to shapefile    

WGScoor<-  MyData #data to convert
coordinates(WGScoor)=~long+lat #column names of the lat long cols
proj4string(WGScoor)<- CRS("++proj=longlat +datum=WGS84") # set coordinate system to WGS
WGScoor.df <- SpatialPointsDataFrame(WGScoor, data.frame(id=1:length(WGScoor)))
LLcoor<-spTransform(WGScoor.df,CRS("+proj=longlat"))
LLcoor.df=SpatialPointsDataFrame(LLcoor, data.frame(id=1:length(LLcoor)))
writeOGR(LLcoor.df, dsn=getwd(),layer="MyShapefile",driver="ESRI Shapefile")
#successfully creates shapefile that opens in arcGIS,
#but shapefile does not include other columns in original data frame

2 Answers 2

14

You are overcomplicating things.

> WGScoor<-  MyData
> coordinates(WGScoor)=~long+lat
> proj4string(WGScoor)<- CRS("+proj=longlat +datum=WGS84")

at this point you have something you can save as a shapefile with all the columns intact, but you seem to want to project out of WGS84, so let's do that:

> LLcoor<-spTransform(WGScoor,CRS("+proj=longlat"))

and let's save this:

> raster::shapefile(LLcoor, "MyShapefile.shp")

I'm using here the shapefile function from the raster package. Its just a nice wrapper to the writeOGR function you used, but makes things a bit easier.

Now let's see what columns are in the resulting shapefile via the command line:

$ ogrinfo -so -al MyShapefile.shp 
INFO: Open of `MyShapefile.shp'
      using driver `ESRI Shapefile' successful.

Layer name: MyShapefile
Geometry: Point
Feature Count: 4
Extent: (-124.450000, 33.350000) - (-123.700000, 44.220000)
Layer SRS WKT:
GEOGCS["GCS_WGS_1984",
    DATUM["unknown",
        SPHEROID["WGS84",6378137,298.257223563]],
    PRIMEM["Greenwich",0],
    UNIT["Degree",0.017453292519943295]]
site: String (80.0)
prop_c: Real (24.15)
prop_b: Real (24.15)
minus_c: Real (24.15)
minus_b: Real (24.15)

which shows they are all there. I'm not sure your transform from "+proj=longlat +datum=WGS84" (and you had an extra + in your string) to "+proj=longlat" actually does anything, but I've left it.

3
  • Thank you! That really helps clear up the steps of the code, and now the shapefile saves with the columns intact.
    – Mon Mo
    Commented Oct 14, 2016 at 23:10
  • @Spacedman What is the explanation of this CRS("+proj=longlat +datum=WGS84")? I cannot understand
    – MK Huda
    Commented Nov 7, 2022 at 0:50
  • It uses the sp package functions to convert a string "+proj=longlat +datum=WGS84" to a coordinate reference system object. CRS("+proj=longlat +datum=WGS84") - note this Q is 6 years old and if you are coding new things now try and use the sf package. If you still don't understand something try posting as a new question.
    – Spacedman
    Commented Nov 7, 2022 at 8:18
1

When you use the function SpatialPointsDataFrame(), the "data.frame" argument is where you get all the data you want in your shape. Now, you're putting only an ID.

Try:

LLcoor.df2=SpatialPointsDataFrame(LLcoor, data.frame(id=1:length(LLcoor), MyData))

You can also do later on:

LLcoor.df2@data <- MyData

However this is risky as you can get order problems. Make sure the other of your shape and your data frame is teh same before replacing it that way.

Hope it helps,

1
  • This method works too. It creates a few additional columns in the attribute table when opened in ArcGIS, but all the essential columns are there.
    – Mon Mo
    Commented Oct 14, 2016 at 23:14

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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