Create sf object from data frame with multiple coordinates

I want to create an sf object from a data frame that contains multiple coordinates under different columns for each row. In the repex below, each ID contains start and end coordinates for both latitude and longitude.

dfr <- structure(list(ID = c("1001A", "1002A", "1003A", "1004A", "1005A", 
"1006A", "1007A", "1008A", "1009A", "1010A"), StartLat = c(33.53418, 
33.60399, 33.40693, 33.64672, 33.57127, 33.42848, 33.54936, 33.49554, 
33.5056, 33.61696), StartLong = c(-112.09114, -111.92731, -112.02982, 
-111.92548, -112.04899, -112.0998, -112.09123, -111.9687, -112.05629, 
-111.98657), EndLat = c(33.53488, 33.60401, 33.40687, 33.64776, 
33.57125, 33.42853, 33.54893, 33.49647, 33.5056, 33.61654), EndLong = c(-112.09114, 
-111.93097, -112.03429, -111.93031, -112.04807, -112.09929, -112.09122, 
-111.97105, -112.0541, -111.98657)), row.names = c(3028L, 8618L, 
6322L, 1171L, 691L, 6590L, 2008L, 4552L, 2894L, 1909L), class = "data.frame")

I tried using the sf package's st_as_sf function, but it yields an error:

dfr_sf <- st_as_sf(dfr, coords = c(c("StartLong", "EndLong"), c("StartLat", "EndLat")),
                   crs = " proj=longlat  datum=WGS84")

Error in points_rcpp(as.matrix(cc), dim) : 
  dim(pts)[2] == nchar(gdim) is not TRUE

The start and end coordinates in each row define a road segment. All the coordinates should go under the geometry column in the final sf object so that these can be plotted as polylines.

CodePudding user response：

Please see if dfr_line is what you want. I think the idea is we need to create a point sf object first, and the convert it to linestring.

library(tidyverse)
library(sf)

dfr2 <- dfr %>%
  pivot_longer(cols = -ID, names_to = "Coord", values_to = "Value") %>%
  extract(Coord, into = c("Type", "Coord"), regex = "(Start|End)(Long|Lat)") %>%
  pivot_wider(names_from = "Coord", values_from = "Value")

dfr_point <- dfr2 %>%
  st_as_sf(coords = c("Long", "Lat"), crs = " proj=longlat  datum=WGS84")

dfr_line <- dfr_point %>%
  select(-Type) %>%
  group_by(ID) %>%
  summarize() %>%
  st_cast("LINESTRING")

CodePudding user response：

Please find an alternative solution with data.table and sf libraries. I hope I am not wrong: the lines appear in Arizona.

library(data.table)
library(sf)
library(dplyr)

# Reshape the data.table 'dfr' into long format
dfr2 <- melt(setDT(dfr), measure= patterns("Long", "Lat"), 
          value.name= c("Longitude", "Latitude"), na.rm=TRUE)[,variable:=NULL][]


# function 'rename_geom_col' to rename the geometry column of the 'sf' object
rename_geom_col <- function(x, new_name){
  current_name  <-  attr(x, "sf_column")
  names(x)[names(x)==current_name]  <- new_name
  st_geometry(x) <- new_name
  return(x)
}

# Convert the data.table 'dfr2' into the linestring 'sf' object 'results'
results <- dfr2 %>% 
  split(., by = "ID", keep.by = FALSE) %>% 
  lapply(., as.matrix) %>% 
  lapply(.,st_linestring) %>%
  st_as_sfc() %>% 
  st_as_sf %>% 
  st_set_crs(4326) %>% 
  mutate(ID = dfr2[1:(nrow(dfr2)/2),"ID"]) %>% 
  relocate(x, .after = last_col()) %>% 
  rename_geom_col(., "geom")


results
#> Simple feature collection with 10 features and 1 field
#> Geometry type: LINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: -112.0998 ymin: 33.40687 xmax: -111.9255 ymax: 33.64776
#> Geodetic CRS:  WGS 84
#>       ID                           geom
#> 1  1001A LINESTRING (-112.0911 33.53...
#> 2  1002A LINESTRING (-111.9273 33.60...
#> 3  1003A LINESTRING (-112.0298 33.40...
#> 4  1004A LINESTRING (-111.9255 33.64...
#> 5  1005A LINESTRING (-112.049 33.571...
#> 6  1006A LINESTRING (-112.0998 33.42...
#> 7  1007A LINESTRING (-112.0912 33.54...
#> 8  1008A LINESTRING (-111.9687 33.49...
#> 9  1009A LINESTRING (-112.0563 33.50...
#> 10 1010A LINESTRING (-111.9866 33.61...

^{Created on 2021-11-16 by the reprex package (v2.0.1)}