Spatial Weight Matrix

[mohdalattar]

I propose the implementation of a feature for constructing spatial weight matrix, as existing methods insufficiently consider the nonplanarity of street networks. This means that bridges are currently treated as adjacent to all streets they cross over. While OSMnx effectively accounts for nonplanarity in calculating street centrality, I am not aware of any method within it that addresses this aspect for spatial weight matrix.

[martinfleis]

Can you expand on that a bit with some code examples? I am not sure I fully understand the use case that is not yet covered.

When measuring street network attributes using functions in momepy like closeness, you need a networkx.Graph object. And our gdf_to_nx conversion function preserves information about non-planarity of networks as it is based on end points of each segment, not an intersection between geometries.

If you are talking about libpysal spatial weights objects representing weights matrix, then it depends on how you generate it but if you use Queen contiguity, the resulting matrix encodes the same as gdf_to_nx, which means that non-planarity is preserved.

[mohdalattar]

My code is

import osmnx as ox
import pysal
import networkx as nx
import libpysal as lps
import matplotlib.pyplot as plt
from collections import defaultdict

place_name = "San Francisco, California, USA"
graph = ox.graph_from_place(place_name, network_type='drive')

# Convert the graph to a GeoDataFrame
nodes, edges = ox.graph_to_gdfs(graph)


# Add a StreetID column to edges
edges['StreetID'] = range(len(edges))

# Initialize an empty dictionary to hold edge adjacencies
edge_adjacencies = defaultdict(set)

# Iterate over edges to find shared nodes and use StreetID
for street_id, (u, v, key) in enumerate(graph.edges(keys=True)):
    # For each node in the edge, add the StreetID to the adjacency set of the node
    edge_adjacencies[u].add(street_id)
    edge_adjacencies[v].add(street_id)

# Create an adjacency list for edges using StreetID
edge_adj_list = {}
for street_id in edges['StreetID']:
    u, v, key = edges.iloc[street_id].name
    shared_u = edge_adjacencies[u]
    shared_v = edge_adjacencies[v]
    # Intersection gives edges that share a node with the current edge
    adjacent_edges = (shared_u | shared_v) - {street_id}
    edge_adj_list[street_id] = list(adjacent_edges)

# Convert the adjacency list to a libpysal W object
w = lps.weights.W(edge_adj_list)

and to verify the results I used

# Specify the StreetID of interest
specific_street_id = 900

# Initialize plot
fig, ax = ox.plot_graph(graph, show=False, close=False, edge_color='black', node_size=0)

# Coordinates for zooming
coords = []

# Plot the specific StreetID edge
u, v, _ = edges.loc[edges['StreetID'] == specific_street_id].index[0]
x1, y1 = graph.nodes[u]['x'], graph.nodes[u]['y']
x2, y2 = graph.nodes[v]['x'], graph.nodes[v]['y']
coords.extend([(x1, y1), (x2, y2)])
ax.plot([x1, x2], [y1, y2], color='blue', linewidth=2, label='StreetID 0')

# Plot each neighbor edge
if specific_street_id in w.neighbors:
    for neighbor_id in w.neighbors[specific_street_id]:
        nu, nv, _ = edges.loc[edges['StreetID'] == neighbor_id].index[0]
        nx1, ny1 = graph.nodes[nu]['x'], graph.nodes[nu]['y']
        nx2, ny2 = graph.nodes[nv]['x'], graph.nodes[nv]['y']
        coords.extend([(nx1, ny1), (nx2, ny2)])
        ax.plot([nx1, nx2], [ny1, ny2], color='green', linestyle='--', linewidth=1, label=f'Neighbor {neighbor_id}')

# Set the axes limits to zoom into the specific StreetID and its neighbors
if coords:
    xs, ys = zip(*coords)
    ax.set_xlim(min(xs) - 0.001, max(xs) + 0.001)  # Adding a small buffer around the edges
    ax.set_ylim(min(ys) - 0.001, max(ys) + 0.001)

# Set plot title and legend
ax.set_title('StreetID 0 and its Neighbors')
ax.legend()

# Show plot
plt.show()

However, I'm uncertain about its accuracy since it's behaving somewhat unpredictably.

I'm kinda happy with the bridge neighbors since it ignores the streets underneath. [ StreetID = 0 ]

But confused with the results below [ StreetID = 900 ]

[martinfleis]

Your confusion comes from the fact that edges in OSM model are directional. It means that a bi-directional street is present twice. Once as an edge going from A to B and then as an edge going from B to A. Notice the overlaps in the dashed line in the second plot. It is a way of representing MultiDiGraph as GeoDataFrames.