Replies: 2 comments 1 reply
-
It's possible to convert the WKT to sfc geometry using {wk} like so. I think having an library(sf)
#> Linking to GEOS 3.11.0, GDAL 3.5.3, PROJ 9.1.0; sf_use_s2() is TRUE
library(dplyr)
library(igraph)
library(sfnetworks)
gurl <- "https://raw.githubusercontent.com/gdsbook/book/master/data/cache/yoyogi_park_graph.graphml"
g <- igraph::read_graph(gurl, format = "graphml")
tg <- tidygraph::as_tbl_graph(g)
nodes <- tg |>
as_tibble() |>
mutate(across(c(x, y), as.numeric)) |>
st_as_sf(coords = c("x", "y"), crs = 4326)
# extract edges as tibble
edges <- tg |>
activate(edges) |>
as_tibble() |>
mutate(geometry = st_as_sfc(wk::wkt(na_if(geometry, "")))) |>
st_as_sf(crs = 4326)
st_geometry(edges)
#> Geometry set for 287 features (with 85 geometries empty)
#> Geometry type: LINESTRING
#> Dimension: XY
#> Bounding box: xmin: 139.6909 ymin: 35.66601 xmax: 139.7014 ymax: 35.67426
#> Geodetic CRS: WGS 84
#> First 5 geometries:
#> LINESTRING (139.6943 35.67009, 139.6942 35.6701...
#> LINESTRING (139.6943 35.67009, 139.6946 35.6700...
#> LINESTRING (139.6943 35.67009, 139.6943 35.67, ...
#> LINESTRING EMPTY
#> LINESTRING (139.6995 35.66972, 139.6994 35.6698... Created on 2023-09-19 with reprex v2.0.2 |
Beta Was this translation helpful? Give feedback.
-
library(sf)
library(tidygraph)
library(dplyr)
library(igraph)
library(tidyr)
library(sfnetworks)
gurl <- "https://raw.githubusercontent.com/gdsbook/book/master/data/cache/yoyogi_park_graph.graphml"
g <- igraph::read_graph(gurl, format = "graphml")
# covert to tbl_graph
tg = tidygraph::as_tbl_graph(g) So after playing a bit more with this example I found a couple of ways someone could create the I first will extract network components # create geometry column from nodes X and Y
nodes_coords = as_tibble(tg) |>
transmute(across(c(x, y), as.numeric)) |>
st_as_sf(coords = c("x", "y"), crs = 4326) |>
st_geometry()
# extract only nodes as sf
nodes_sf = tg |>
as_tibble() |>
mutate(geometry = nodes_coords) |>
st_as_sf(crs = 4326)
# add geometry to the `tbl_graph` nodes
nodes_tg = tg |>
activate(nodes) |>
mutate(geometry = nodes_coords)
# extract only edges as tibble
edges = tg |>
activate(edges) |>
as_tibble()
# create an sf object from the edges including empty geometries
edges_empty = edges |>
mutate(geometry = st_as_sfc(wk::wkt(na_if(geometry, "")))) |>
st_as_sf(crs = 4326)
# create an sf object from the edges excluding empty geometries
edges_filt = edges |>
mutate(geometry = na_if(geometry, "")) |>
drop_na() |>
mutate(geometry = st_as_sfc(wk::wkt(geometry))) |>
st_as_sf(crs = 4326) Now we can create sfnetwork objects in a couple of ways # using the nodes as sf objects and the empty geometries
# we need to use force=TRUE otherwise building the network fails
sfn_empty = sfnetwork(nodes_sf, edges_empty, force = TRUE)
# using the nodes as sf objects and the filterd geometries
sfn_filt = sfnetwork(nodes_sf, edges_filt, force = FALSE)
#> Checking if spatial network structure is valid...
#> Spatial network structure is valid
# using the tbl_graph object which contains nodes and edges
# but the edges don't have an explicit geometry
sfn_tg = as_sfnetwork(nodes_tg)
#> Checking if spatial network structure is valid...
#> Spatial network structure is valid
# using only the nodes as sf objects and letting sfnetworks
# build the edges
sfn_nodes = as_sfnetwork(nodes_sf)
par(mar = c(1, 1, 1, 1), mfrow = c(2,2))
plot(sfn_empty, main = "Includes EMPTY geoms")
plot(sfn_filt, main = "Excludes EMPTY geoms")
plot(sfn_tg, main = "Built based on tidygraph with added geom column")
plot(sfn_nodes, main = "Built only from node geom information") Now the first two networks look quite OK, mainly when we compare to how the Geographic Data Science with Python book renders it: But you will see there are some edges missing, which correspond, I assume to those empty geometries. We can see that when checking the number of nodes and edges # edges
ecount(sfn_empty)
#> [1] 287
ecount(sfn_filt)
#> [1] 202
# nodes
vcount(sfn_empty)
#> [1] 106
vcount(sfn_filt)
#> [1] 106 I am not sure how osmnx builds the network but with sfnetworks some extra steps would be needed to get all the information from the GRAPHML file. This is one way I solved this based on the sfnetwork that includes the EMPTY edge geometries # Filter the network to contain edges with a geometry
# and assign to a new network
sfn_explicit_edges = sfn_empty |>
activate(edges) |>
filter(!st_is_empty(geometry))
# Filter the network for EMPTY edge geometry
# drop the geometry column for the edges
# remove isolated nodes
# and finally convert the edges to spatially
# explicit with a morpher
sfn_empty_edges = sfn_empty |>
activate(edges) |>
filter(st_is_empty(geometry)) |>
st_drop_geometry() |>
activate("nodes") |>
filter(!node_is_isolated()) |>
convert(to_spatial_explicit)
# Finally, join both networks
sfn_complete = st_network_join(sfn_empty_edges, sfn_explicit_edges)
vcount(sfn_complete)
#> [1] 106
ecount(sfn_complete)
#> [1] 287
par(mar = c(1, 1, 1, 1), mfrow = c(1,1))
plot(sfn_complete) After all this exploration, a couple of ideas popped as to how to handle this cases:
@luukvdmeer do you think this is something we could support? Otherwise, keeping this discussion and maybe turning it into a vignette could be worth for these cases? |
Beta Was this translation helpful? Give feedback.
-
From #256 @JosiahParry shared a dataset from the Geographic Data Science with Python book.
If I take the way the graph was created using tidygraph, you will notice there is actually a geometry column for the edges.
Josiah built the sfnetwork from the nodes by parsing the x and y column from the nodes into an sf geometry object and copying it back in:
This works, BUT it basically ignores the original geometry stored on the edges. Why? Because the geometry in the edges is not a list-column, so it is not recognized. When calling
g_sf |> tidygraph::convert(to_spatial_explicit)
we are basically creating straight lines between the nodes coordinates, losing the original topology stored in the edges.So instead of building the sfnetwork from the nodes as Josiah did, I tried to build it from the edges but got into problems. In escence, the geometry column has empty rows, and sf does not support this (see r-spatial/sf#1034).
So one way to solve this would be to drop all empty geometries as:
But this will break the original node data since we will be dropping information. Another thing would be to create straight lines between the nodes surrounding the edges with empty geometries since those coordinates are available from the original node data.
I am not sure what the rest of the analysis looks like for your example Josiah but I guess at least in terms of plotting you are not getting the same results as in the Python book.
Also, how do they handle the empty geometries? I would be curious to know that.
Beta Was this translation helpful? Give feedback.
All reactions