Following an email exchange with Rafael Pereira, and the excellent StackOverflow answer from Robin Lovelace, I wrote up a function, points_to_line
, that is designed to resemble the “Points to Line” tool available in ArcGIS. I borrow the approach slightly as well from the book Applied Spatial Data Analysis with R.
The function takes five arguments:
data
(required): a data frame with XY coordinateslong
(required): the column containing the X coordinateslat
(required): the column containing the Y coordinatesid_field
(optional): field that identifies unique lines; not needed if there is only one linesort_field
(optional): the field that identifies the sequence the points should be in.The function depends on the sp and maptools packages. It will return a SpatialLines
object; if you want a SpatialLinesDataFrame
, you’ll need to add in the data you want afterwards.
library(sp)
library(maptools)
points_to_line <- function(data, long, lat, id_field = NULL, sort_field = NULL) {
# Convert to SpatialPointsDataFrame
coordinates(data) <- c(long, lat)
# If there is a sort field...
if (!is.null(sort_field)) {
if (!is.null(id_field)) {
data <- data[order(data[[id_field]], data[[sort_field]]), ]
} else {
data <- data[order(data[[sort_field]]), ]
}
}
# If there is only one path...
if (is.null(id_field)) {
lines <- SpatialLines(list(Lines(list(Line(data)), "id")))
return(lines)
# Now, if we have multiple lines...
} else if (!is.null(id_field)) {
# Split into a list by ID field
paths <- sp::split(data, data[[id_field]])
sp_lines <- SpatialLines(list(Lines(list(Line(paths[[1]])), "line1")))
# I like for loops, what can I say...
for (p in 2:length(paths)) {
id <- paste0("line", as.character(p))
l <- SpatialLines(list(Lines(list(Line(paths[[p]])), id)))
sp_lines <- spRbind(sp_lines, l)
}
return(sp_lines)
}
}
Now let’s put the function to work with the GTFS dataset of transit routes in Las Vegas Rafael sent me, and plot the result on a Leaflet map. To reproduce, download the data file from my miscdata GitHub repository.
library(leaflet)
dat <- read.csv('shapes.txt')
v_lines <- points_to_line(data = dat,
long = "shape_pt_lon",
lat = "shape_pt_lat",
id_field = "shape_id",
sort_field = "shape_pt_sequence")
leaflet(data = v_lines) %>%
addTiles() %>%
addPolylines()
The function can certainly be improved, but it works! Please let me know if you have any comments or suggestions.