Output shapefile for the igraph network in R

Question

Hello I have a network in R using the igraph library

Vertices: 616 
Edges: 6270 
Directed: TRUE 
No graph attributes.
Vertex attributes: name, Lat, Lon.
Edge attributes: V3.

How can I generate two shapefiles for the Vertices and the Edges using the Lat, Lon info in the vertex?

Amen

Answer 1

You can do this using the sp and maptools packages. There are handy functions writePointsShape() and writeLinesShape() in maptools that will write to the ESRI shapefile format.

Before doing this, it is necessary to extract the lat/lon information from the graph vertices and put it into a SpatialPoints object for the vertices, and a SpatialLinesDataFrame object for the edges.

This code produces a very simple igraph object for the following example:

library(igraph)

## Produce a ring graph with 4 vertices
x <- graph.ring(4)

## Add lat/lon information to vertices
V(x)$lat <- c(50, 50, 51, 51)
V(x)$lon <- c(40, 41, 41, 40)

Now, create the SpatialPoints object for the vertices

library(sp)
library(maptools)

## Create SpatialPoints object containing coordinates
xV <- SpatialPoints(cbind(V(x)$lon, V(x)$lat))

## Write vertices to a shapefile
writePointsShape(xV, fn="vertices")

Finally, create the SpatialLinesDataFrame object for the edges. This is a little messy, but I am yet to find a quick way to produce a SpatialLines object given coordinates.

## Create SpatialLinesDataFrame object describing edges
edges <- get.edgelist(x)+1
edges <- cbind(edgeNum=1:nrow(edges), v1=edges[,1], v2=edges[,2])
xE <- apply(edges, 1, function(i) Lines(Line(cbind(c(V(x)$lon[i["v1"]], V(x)$lon[i["v2"]]), c(V(x)$lat[i["v1"]], V(x)$lat[i["v2"]]))), ID=as.character(i["edgeNum"])))
xE <- SpatialLinesDataFrame(SpatialLines(xE), data=data.frame(edgeNum=1:nrow(edges)))

## Write edges to a shapefile
writeLinesShape(xE, fn="edges")

Nerine

Had transposed Lat and Lon. Corrected now. Repeat after me: (Lon=X, Lat=Y) :) – Nerine

Error in validityMethod(object) : coords cannot contain missing values. That is how I get for the "apply" function for the line. Do you have any suggestions that I can check the values. Sorry I am a new hand in R... – Amen

I am assuming that you are running this on your 616 vertex network, not the example 4 vertex ring given. Sounds like there are some vertices that don't have lat/lon information. You can check using print(yourGraph, vertex=T) – Nerine

Finally I made it work! Thanks! The last question remaining is to distinguish the directions of the edges of the graph.... – Amen

Not sure how you would represent a directed edge (i.e. arc) in a shapefile. Are arcs/arrows part of the ESRI shapefile specification? – Nerine