Morphological Graphs from Overture Maps and OpenStreetMap with City2Graph¶

# Core libraries
import numpy as np
import geopandas as gpd
import matplotlib.pyplot as plt
import contextily as cx
from shapely import Point

# city2graph for morphological graph analysis
import city2graph as c2g

# Configure matplotlib for high-quality visualizations
plt.rcParams['figure.figsize'] = (14, 10)
plt.rcParams['figure.dpi'] = 100
plt.rcParams['font.size'] = 11
plt.style.use('ggplot')

print("Dependencies loaded successfully")
print(f"city2graph version: {c2g.__version__}")
print(f"geopandas version: {gpd.__version__}")

# Download data from Overture Maps (uncomment to download fresh data)
c2g.load_overture_data(
    place_name="Liverpool, UK",
    types=["segment", "building", "connector"],
    output_dir="./data/overturemaps",
    prefix="liverpool_",
    save_to_file=True,
    return_data=False,
)

# Load the downloaded GeoJSON files
buildings_gdf = gpd.read_file("./data/overturemaps/liverpool_building.geojson")
segments_gdf = gpd.read_file("./data/overturemaps/liverpool_segment.geojson") 
connectors_gdf = gpd.read_file("./data/overturemaps/liverpool_connector.geojson")

# Convert to British National Grid (EPSG:27700) for accurate distance calculations
buildings_gdf = buildings_gdf.to_crs(epsg=27700)
segments_gdf = segments_gdf.to_crs(epsg=27700)
connectors_gdf = connectors_gdf.to_crs(epsg=27700)

print("Data loaded successfully")
print("Dataset summary:")
print(f"   - Buildings: {len(buildings_gdf):,}")
print(f"   - Road segments: {len(segments_gdf):,}")
print(f"   - Connectors: {len(connectors_gdf):,}")
print(f"   - CRS: {buildings_gdf.crs}")

# Filter to keep only road segments (excluding pedestrian paths, railways, etc.)
segments_gdf = segments_gdf[segments_gdf["subtype"] == "road"].copy()

# Process segments to handle bridges/tunnels and create proper spatial barriers
segments_gdf = c2g.process_overture_segments(
    segments_gdf=segments_gdf,
    get_barriers=True,
    connectors_gdf=connectors_gdf
)

print(f"Processed {len(segments_gdf)} road segments")
print("Barrier geometries created for tessellation")

# Check the geometry types in the barrier_geometry column
geometry_types = segments_gdf["barrier_geometry"].geom_type.value_counts()
print("\nBarrier geometry types:")
for geom_type, count in geometry_types.items():
    print(f"   - {geom_type}: {count:,}")

segments_G = c2g.segments_to_graph(segments_gdf, as_nx=True)

c2g.plot_graph(segments_G)

# Define center point for the analysis area (Liverpool city centre)
center_point = gpd.GeoSeries([Point(-2.9879004, 53.4062724)], crs='EPSG:4326').to_crs(epsg=27700)

# Create the morphological graph
morpho_nodes, morpho_edges = c2g.morphological_graph(
    buildings_gdf=buildings_gdf,
    segments_gdf=segments_gdf,
    center_point=center_point,
    distance=500,                    # Analysis radius in meters
    clipping_buffer=300,            # Buffer for edge effects
    primary_barrier_col='barrier_geometry',
    contiguity="queen",             # Adjacency rule for tessellation
    keep_buildings=True,            # Preserve building geometries
)

print(f"Node types: {list(morpho_nodes.keys())}")
print(f"Edge types: {list(morpho_edges.keys())}")
print(f"Private spaces: {len(morpho_nodes['private']):,}")
print(f"Public spaces: {len(morpho_nodes['public']):,}")

for edge_type, edge_gdf in morpho_edges.items():
    print(f"   • {edge_type}: {len(edge_gdf):,} connections")

morpho_nodes["private"].head()

morpho_nodes["public"].head()

morpho_edges[('public', 'connected_to', 'public')].head()

morpho_edges[('private', 'faced_to', 'public')].head()

morpho_edges[('private', 'touched_to', 'private')].head()

morpho_nodes["public"]

fig, ax = plt.subplots(figsize=(14, 12), facecolor='#f9f9f9')

morpho_nodes["private"]["building_geometry"].plot(ax=ax, color='#e0e0e0', edgecolor='#c0c0c0', linewidth=0.3, alpha=0.7)
morpho_nodes["public"]["segment_geometry"].plot(ax=ax, color='#404040', linewidth=0.7, alpha=0.6)

c2g.plot_graph(
    nodes=morpho_nodes,
    edges=morpho_edges,
    subplots=False,
    ax=ax,
    bgcolor="#f9f9f9",
    labelcolor="#000000",
    node_color={'private': 'red', 'public': 'blue'},
    edge_color={
        ('private', 'touched_to', 'private'): '#B22222',
        ('public', 'connected_to', 'public'): '#0000FF',
        ('private', 'faced_to', 'public'): '#7B68EE',
    },
    markersize=10,
    linewidth=0.2,
)

# Check for PyTorch Geometric availability
import torch

try:
    import torch_geometric
    HAS_TORCH_GEOMETRIC = True
    print("PyTorch Geometric imported successfully")
    print(f"   PyTorch version: {torch.__version__}")
    print(f"   PyTorch Geometric version: {torch_geometric.__version__}")
    print("Ready for graph representation learning")
except ImportError:
    HAS_TORCH_GEOMETRIC = False
    print("PyTorch Geometric is not installed")
    print("Install with: pip install torch-geometric")
    print("Graph conversion examples will continue with available backends")

# Compute spatial features for private spaces (tessellation cells)
morpho_nodes["private"]["area"] = morpho_nodes["private"].geometry.area
morpho_nodes["private"]["perimeter"] = morpho_nodes["private"].geometry.length
morpho_nodes["private"]["compactness"] = (
    4 * np.pi * morpho_nodes["private"]["area"] / (morpho_nodes["private"]["perimeter"] ** 2)
)

print("Computed spatial features for private spaces:")
feature_stats = morpho_nodes["private"][['area', 'perimeter', 'compactness']].describe()
print(feature_stats.round(2))

# Set proper index names for graph conversion
morpho_nodes["private"].index.name = "tess_id"
morpho_edges[('private', 'touched_to', 'private')].index.names = ['from_private_id', 'to_private_id']

print("\nFeatures computed and indices configured")

if HAS_TORCH_GEOMETRIC:
    # Create homogeneous graph of private spaces
    private_graph = c2g.gdf_to_pyg(
        nodes=morpho_nodes["private"],
        edges=morpho_edges[('private', 'touched_to', 'private')],
        node_feature_cols=['area', 'perimeter', 'compactness'],
        edge_feature_cols=None
    )
    
    print("Created homogeneous graph:")
    print(f"   - Nodes: {private_graph.num_nodes:,}")
    print(f"   - Edges: {private_graph.num_edges:,}")
    print(f"   - Node features: {private_graph.x.shape[1]} dimensions")
    print(f"   - Feature tensor shape: {private_graph.x.shape}")
    print("\nGraph is ready for training")
    
else:
    print("Skipping PyTorch Geometric conversion (not installed)")
    private_graph = None

# Define center point for the analysis area (Liverpool city centre)
center_point = gpd.GeoSeries([Point(-2.9879004, 53.4062724)], crs='EPSG:4326').to_crs(epsg=27700)

# Create the morphological graph
print("Creating morphological graph...")
morpho_nodes, morpho_edges = c2g.morphological_graph(
    buildings_gdf=buildings_gdf,
    segments_gdf=segments_gdf,
    center_point=center_point,
    distance=500,                    # Analysis radius in meters
    clipping_buffer=300,            # Buffer for edge effects
    primary_barrier_col='barrier_geometry',
    contiguity="queen",             # Adjacency rule for tessellation
    keep_buildings=True,            # Preserve building geometries
)

print("Morphological graph created successfully")
print("Network summary:")
print(f"   - Node types: {list(morpho_nodes.keys())}")
print(f"   - Edge types: {list(morpho_edges.keys())}")
print(f"   - Private spaces: {len(morpho_nodes['private']):,}")
print(f"   - Public spaces: {len(morpho_nodes['public']):,}")

for edge_type, edge_gdf in morpho_edges.items():
    print(f"   - {edge_type}: {len(edge_gdf):,} connections")

# Recompute node features after rebuilding morpho_nodes
morpho_nodes["private"]["area"] = morpho_nodes["private"].geometry.area
morpho_nodes["private"]["perimeter"] = morpho_nodes["private"].geometry.length
morpho_nodes["private"]["compactness"] = (
    4 * np.pi * morpho_nodes["private"]["area"] / (morpho_nodes["private"]["perimeter"] ** 2)
)
morpho_nodes["public"]["length"] = morpho_nodes["public"].geometry.length

# Keep index names explicit so downstream conversions preserve structure
morpho_nodes["private"].index.name = "private_id"
morpho_nodes["public"].index.name = "public_id"
morpho_edges[('private', 'touched_to', 'private')].index.names = ['from_private_id', 'to_private_id']
morpho_edges[('public', 'connected_to', 'public')].index.names = ['from_public_id', 'to_public_id']
morpho_edges[('private', 'faced_to', 'public')].index.names = ['from_private_id', 'to_public_id']

print("Computed spatial features for private and public nodes")
print("Private feature columns:", ["area", "perimeter", "compactness"])
print("Public feature columns:", ["length"])

if HAS_TORCH_GEOMETRIC:
    # Create heterogeneous graph with all morphological relationships
    hetero_graph = c2g.gdf_to_pyg(
        morpho_nodes,
        morpho_edges,
        node_feature_cols={
            "private": ['area', 'perimeter', 'compactness'],
            "public": ['length']
        }
    )
    
    print("Created heterogeneous morphological graph:")
    print("Node types and counts:")
    for node_type in hetero_graph.node_types:
        num_nodes = hetero_graph[node_type].x.shape[0]
        num_features = hetero_graph[node_type].x.shape[1]
        print(f"   - {node_type}: {num_nodes:,} nodes, {num_features} features")
    
    print("Edge types and counts:")
    for edge_type in hetero_graph.edge_types:
        num_edges = hetero_graph[edge_type].edge_index.shape[1]
        print(f"   - {edge_type}: {num_edges:,} edges")
    
    print("\nComplete morphological graph is ready for heterogeneous GNN workflows")
    
else:
    print("Skipping PyTorch Geometric conversion (not installed)")
    hetero_graph = None

hetero_graph

# Convert morphological graph to PyG and NetworkX
morpho_graph = c2g.gdf_to_pyg(
    morpho_nodes,
    morpho_edges,
    node_feature_cols={
        "private": ['area', 'perimeter', 'compactness'],
        "public": ['length']
    },
)

# Convert morphological graph to NetworkX for traditional graph analysis
G_morpho = c2g.gdf_to_nx(morpho_nodes, morpho_edges, keep_geom=True)

print("Converted to NetworkX:")
print(f"   - Total nodes: {G_morpho.number_of_nodes():,}")
print(f"   - Total edges: {G_morpho.number_of_edges():,}")
print(f"   - Graph type: {'Directed' if G_morpho.is_directed() else 'Undirected'}")

# Check heterogeneity and inspect attributes
has_node_types = any('node_type' in data for _, data in G_morpho.nodes(data=True))
print(f"   - Heterogeneous: {'Yes' if has_node_types else 'No'}")

sample_node_id, sample_node_data = next(iter(G_morpho.nodes(data=True)))
non_private_keys = [k for k in sample_node_data.keys() if not k.startswith('_')]
print(f"   - Sample node id: {sample_node_id}")
print(f"   - Sample node attributes: {non_private_keys}")

if has_node_types:
    from collections import Counter

    node_types = [data.get('node_type', 'unknown') for _, data in G_morpho.nodes(data=True)]
    type_counts = Counter(node_types)
    print("   - Node type distribution:")
    for node_type, count in type_counts.items():
        print(f"     - {node_type}: {count:,}")

print("\nNetworkX graph is ready for analysis")

import osmnx as ox

# Download and create a street network graph for Liverpool, UK
G_osmnx = ox.graph_from_place("Liverpool, UK", network_type="all")

print("OSMnx graph created for Liverpool, UK")
print(f"   - Nodes: {len(G_osmnx.nodes):,}")
print(f"   - Edges: {len(G_osmnx.edges):,}")

# Download POIs (amenities, shops, etc.) in Liverpool, UK
tags = {"building": True}
osm_buildings = ox.features_from_place("Liverpool, UK", tags)

print(f"Retrieved {len(osm_buildings):,} buildings from OSMnx")
osm_buildings.head()

osm_data

osm_street_nodes.head()

osm_street_edges.head()

# Define center point for the analysis area (Liverpool city centre)
center_point = gpd.GeoSeries([Point(-2.9879004, 53.4062724)], crs='EPSG:4326').to_crs(epsg=27700)

# Create the morphological graph
print("Creating morphological graph...")
morpho_nodes, morpho_edges = c2g.morphological_graph(
    buildings_gdf=osm_buildings,
    segments_gdf=osm_street_edges,
    center_point=center_point,
    distance=500,                    # Analysis radius in meters
    clipping_buffer=300,            # Buffer for edge effects
    primary_barrier_col='barrier_geometry',
    contiguity="queen",             # Adjacency rule for tessellation
    keep_buildings=True,            # Preserve building geometries
    keep_segments=True,             # Preserve street segment geometries
)

print(f"   - Node types: {list(morpho_nodes.keys())}")
print(f"   - Edge types: {list(morpho_edges.keys())}")
print(f"   - Private spaces: {len(morpho_nodes['private']):,}")
print(f"   - Public spaces: {len(morpho_nodes['public']):,}")

for edge_type, edge_gdf in morpho_edges.items():
    print(f"   - {edge_type}: {len(edge_gdf):,} connections")

fig, ax = plt.subplots(figsize=(14, 12), facecolor='#f9f9f9')

morpho_nodes["private"]["building_geometry"].plot(ax=ax, color='#e0e0e0', edgecolor='#c0c0c0', linewidth=0.3, alpha=0.7)
morpho_nodes["public"]["segment_geometry"].plot(ax=ax, color='#404040', linewidth=0.7, alpha=0.6)

c2g.plot_graph(
    nodes=morpho_nodes,
    edges=morpho_edges,
    subplots=False,
    ax=ax,
    bgcolor="#f9f9f9",
    labelcolor="#000000",
    node_color={'private': 'red', 'public': 'blue'},
    edge_color={
        ('private', 'touched_to', 'private'): '#B22222',
        ('public', 'connected_to', 'public'): '#0000FF',
        ('private', 'faced_to', 'public'): '#7B68EE',
    },
    markersize=10,
    linewidth=0.2,
)

osm_hetero_data