Adjacency neighbourhood
In [1]:
Copied!
import geopandas as gpd
import numpy as np
from shapely.geometry import Point
from srai.constants import WGS84_CRS
from srai.neighbourhoods import AdjacencyNeighbourhood
from srai.plotting.folium_wrapper import plot_all_neighbourhood, plot_neighbours, plot_regions
from srai.regionalizers import (
AdministrativeBoundaryRegionalizer,
VoronoiRegionalizer,
geocode_to_region_gdf,
)
import geopandas as gpd
import numpy as np
from shapely.geometry import Point
from srai.constants import WGS84_CRS
from srai.neighbourhoods import AdjacencyNeighbourhood
from srai.plotting.folium_wrapper import plot_all_neighbourhood, plot_neighbours, plot_regions
from srai.regionalizers import (
AdministrativeBoundaryRegionalizer,
VoronoiRegionalizer,
geocode_to_region_gdf,
)
Adjacency Neighbourhood¶
It can generate neighbourhoods for all geodataframes with touching geometries.
Real boundaries example - Italy¶
In [2]:
Copied!
it_gdf = geocode_to_region_gdf(query=["R365331"], by_osmid=True)
plot_regions(it_gdf)
it_gdf = geocode_to_region_gdf(query=["R365331"], by_osmid=True)
plot_regions(it_gdf)
Out[2]:
Make this Notebook Trusted to load map: File -> Trust Notebook
In [3]:
Copied!
regionalizer = AdministrativeBoundaryRegionalizer(admin_level=4)
it_regions_gdf = regionalizer.transform(it_gdf)
regionalizer = AdministrativeBoundaryRegionalizer(admin_level=4)
it_regions_gdf = regionalizer.transform(it_gdf)
Loading boundaries: 0: 0%| | 0/6 [00:00<?, ?it/s]
Loading boundaries: 1: 0%| | 0/6 [00:00<?, ?it/s]
Loading boundaries: 2: 0%| | 0/6 [00:00<?, ?it/s]
Loading boundaries: 3: 0%| | 0/6 [00:00<?, ?it/s]
Loading boundaries: 4: 0%| | 0/6 [00:00<?, ?it/s]
Loading boundaries: 5: 0%| | 0/6 [00:00<?, ?it/s]
Loading boundaries: 6: 0%| | 0/6 [00:01<?, ?it/s]
Loading boundaries: 7: 0%| | 0/6 [00:01<?, ?it/s]
Loading boundaries: 8: 0%| | 0/6 [00:01<?, ?it/s]
Loading boundaries: 9: 0%| | 0/6 [00:01<?, ?it/s]
Loading boundaries: 10: 0%| | 0/6 [00:02<?, ?it/s]
Loading boundaries: 11: 0%| | 0/6 [00:02<?, ?it/s]
Loading boundaries: 12: 0%| | 0/6 [00:02<?, ?it/s]
Loading boundaries: 13: 0%| | 0/6 [00:03<?, ?it/s]
Loading boundaries: 14: 0%| | 0/6 [00:03<?, ?it/s]
Loading boundaries: 15: 0%| | 0/6 [00:03<?, ?it/s]
Loading boundaries: 16: 0%| | 0/6 [00:04<?, ?it/s]
Loading boundaries: 17: 0%| | 0/6 [00:04<?, ?it/s]
Loading boundaries: 18: 0%| | 0/6 [00:04<?, ?it/s]
Loading boundaries: 19: 0%| | 0/6 [00:05<?, ?it/s]
Loading boundaries: 20: 0%| | 0/6 [00:05<?, ?it/s]
Loading boundaries: 21: 0%| | 0/6 [00:05<?, ?it/s]
Loading boundaries: 22: 0%| | 0/6 [00:06<?, ?it/s]
Loading boundaries: 23: 0%| | 0/6 [00:06<?, ?it/s]
Loading boundaries: 24: 0%| | 0/6 [00:06<?, ?it/s]
Loading boundaries: 25: 0%| | 0/6 [00:07<?, ?it/s]
Loading boundaries: 26: 0%| | 0/6 [00:07<?, ?it/s]
Loading boundaries: 27: 0%| | 0/6 [00:08<?, ?it/s]
Loading boundaries: 28: 0%| | 0/6 [00:08<?, ?it/s]
Loading boundaries: 29: 0%| | 0/6 [00:09<?, ?it/s]
Loading boundaries: 30: 0%| | 0/6 [00:09<?, ?it/s]
Loading boundaries: 31: 0%| | 0/6 [00:10<?, ?it/s]
Loading boundaries: 32: 0%| | 0/6 [00:10<?, ?it/s]
Loading boundaries: 33: 0%| | 0/6 [00:10<?, ?it/s]
Loading boundaries: 34: 0%| | 0/6 [00:11<?, ?it/s]
Loading boundaries: 35: 0%| | 0/6 [00:11<?, ?it/s]
Loading boundaries: 36: 0%| | 0/6 [00:11<?, ?it/s]
Loading boundaries: 37: 0%| | 0/6 [00:12<?, ?it/s]
Loading boundaries: 38: 0%| | 0/6 [00:12<?, ?it/s]
Loading boundaries: 39: 0%| | 0/6 [00:13<?, ?it/s]
Loading boundaries: 40: 0%| | 0/6 [00:13<?, ?it/s]
Loading boundaries: 41: 0%| | 0/6 [00:13<?, ?it/s]
Loading boundaries: 42: 0%| | 0/6 [00:14<?, ?it/s]
Loading boundaries: 43: 0%| | 0/6 [00:14<?, ?it/s]
Loading boundaries: 44: 0%| | 0/6 [00:15<?, ?it/s]
Loading boundaries: 45: 0%| | 0/6 [00:15<?, ?it/s]
Loading boundaries: 46: 0%| | 0/6 [00:15<?, ?it/s]
Loading boundaries: 47: 0%| | 0/6 [00:16<?, ?it/s]
Loading boundaries: 48: 0%| | 0/6 [00:16<?, ?it/s]
Loading boundaries: 49: 0%| | 0/6 [00:17<?, ?it/s]
Loading boundaries: 50: 0%| | 0/6 [00:17<?, ?it/s]
Loading boundaries: 51: 0%| | 0/6 [00:17<?, ?it/s]
Loading boundaries: 52: 0%| | 0/6 [00:18<?, ?it/s]
Loading boundaries: 53: 0%| | 0/6 [00:18<?, ?it/s]
Loading boundaries: 54: 0%| | 0/6 [00:19<?, ?it/s]
Loading boundaries: 55: 0%| | 0/6 [00:19<?, ?it/s]
Loading boundaries: 56: 0%| | 0/6 [00:19<?, ?it/s]
Loading boundaries: 57: 0%| | 0/6 [00:20<?, ?it/s]
Loading boundaries: 58: 0%| | 0/6 [00:20<?, ?it/s]
Loading boundaries: 59: 0%| | 0/6 [00:20<?, ?it/s]
Loading boundaries: 60: 0%| | 0/6 [00:21<?, ?it/s]
Loading boundaries: 61: 0%| | 0/6 [00:21<?, ?it/s]
Loading boundaries: 62: 0%| | 0/6 [00:22<?, ?it/s]
Loading boundaries: 63: 0%| | 0/6 [00:22<?, ?it/s]
Loading boundaries: 64: 0%| | 0/6 [00:23<?, ?it/s]
Loading boundaries: 65: 0%| | 0/6 [00:23<?, ?it/s]
Loading boundaries: 66: 0%| | 0/6 [00:24<?, ?it/s]
Loading boundaries: 67: 0%| | 0/6 [00:25<?, ?it/s]
Loading boundaries: 68: 0%| | 0/6 [00:25<?, ?it/s]
Loading boundaries: 68: 17%|█▋ | 1/6 [00:25<02:09, 25.85s/it]
Loading boundaries: 68: 33%|███▎ | 2/6 [00:26<00:43, 10.82s/it]
Loading boundaries: 68: 50%|█████ | 3/6 [00:26<00:18, 6.00s/it]
Loading boundaries: 68: 67%|██████▋ | 4/6 [00:26<00:07, 3.74s/it]
Loading boundaries: 68: 100%|██████████| 6/6 [00:26<00:00, 1.83s/it]
Loading boundaries: 68: 100%|██████████| 6/6 [00:26<00:00, 4.49s/it]
/opt/hostedtoolcache/Python/3.10.14/x64/lib/python3.10/site-packages/srai/regionalizers/administrative_boundary_regionalizer.py:168: FutureWarning: `unary_union` returned None due to all-None GeoSeries. In future, `unary_union` will return 'GEOMETRYCOLLECTION EMPTY' instead. ].geometry.unary_union
In [4]:
Copied!
plot_regions(it_regions_gdf)
plot_regions(it_regions_gdf)
Out[4]:
Make this Notebook Trusted to load map: File -> Trust Notebook
In [5]:
Copied!
neighbourhood = AdjacencyNeighbourhood(it_regions_gdf)
neighbourhood = AdjacencyNeighbourhood(it_regions_gdf)
Nearest neighbours¶
In [6]:
Copied!
region_id = "Lazio"
neighbours = neighbourhood.get_neighbours(region_id)
neighbours
region_id = "Lazio"
neighbours = neighbourhood.get_neighbours(region_id)
neighbours
Out[6]:
{'Abruzzo', 'Campania', 'Marche', 'Molise', 'Tuscany', 'Umbria'}
In [7]:
Copied!
plot_neighbours(it_regions_gdf, region_id, neighbours)
plot_neighbours(it_regions_gdf, region_id, neighbours)
Out[7]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Neighbours at a distance¶
In [8]:
Copied!
region_id = "Basilicata"
neighbours = neighbourhood.get_neighbours_at_distance(region_id, 2)
neighbours
region_id = "Basilicata"
neighbours = neighbourhood.get_neighbours_at_distance(region_id, 2)
neighbours
Out[8]:
{'Lazio', 'Molise'}
In [9]:
Copied!
plot_neighbours(it_regions_gdf, region_id, neighbours)
plot_neighbours(it_regions_gdf, region_id, neighbours)
Out[9]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Regions without neighbours¶
In [10]:
Copied!
region_id = "Sardinia"
neighbours = neighbourhood.get_neighbours(region_id)
neighbours
region_id = "Sardinia"
neighbours = neighbourhood.get_neighbours(region_id)
neighbours
Out[10]:
set()
In [11]:
Copied!
plot_neighbours(it_regions_gdf, region_id, neighbours)
plot_neighbours(it_regions_gdf, region_id, neighbours)
Out[11]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Plotting all neighbourhood¶
In [12]:
Copied!
region_id = "Campania"
plot_all_neighbourhood(it_regions_gdf, region_id, neighbourhood)
region_id = "Campania"
plot_all_neighbourhood(it_regions_gdf, region_id, neighbourhood)
Out[12]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Voronoi example - Australia¶
In [13]:
Copied!
au_gdf = geocode_to_region_gdf(query=["R80500"], by_osmid=True)
plot_regions(au_gdf)
au_gdf = geocode_to_region_gdf(query=["R80500"], by_osmid=True)
plot_regions(au_gdf)
Out[13]:
Make this Notebook Trusted to load map: File -> Trust Notebook
In [14]:
Copied!
def generate_random_points(shape, n_points=500):
"""Generates random points."""
minx, miny, maxx, maxy = shape.bounds
pts = []
rng = np.random.default_rng()
while len(pts) < 4:
randx = rng.uniform(minx, maxx, n_points)
randy = rng.uniform(miny, maxy, n_points)
coords = np.vstack((randx, randy)).T
# use only the points inside the geographic area
pts = [p for p in list(map(Point, coords)) if p.within(shape)]
del coords # not used any more
return pts
def generate_random_points(shape, n_points=500):
"""Generates random points."""
minx, miny, maxx, maxy = shape.bounds
pts = []
rng = np.random.default_rng()
while len(pts) < 4:
randx = rng.uniform(minx, maxx, n_points)
randy = rng.uniform(miny, maxy, n_points)
coords = np.vstack((randx, randy)).T
# use only the points inside the geographic area
pts = [p for p in list(map(Point, coords)) if p.within(shape)]
del coords # not used any more
return pts
In [15]:
Copied!
pts = generate_random_points(au_gdf.geometry[0])
au_seeds_gdf = gpd.GeoDataFrame(
{"geometry": pts},
index=list(range(len(pts))),
crs=WGS84_CRS,
)
pts = generate_random_points(au_gdf.geometry[0])
au_seeds_gdf = gpd.GeoDataFrame(
{"geometry": pts},
index=list(range(len(pts))),
crs=WGS84_CRS,
)
/opt/hostedtoolcache/Python/3.10.14/x64/lib/python3.10/site-packages/geopandas/geoseries.py:641: FutureWarning: Series.__getitem__ treating keys as positions is deprecated. In a future version, integer keys will always be treated as labels (consistent with DataFrame behavior). To access a value by position, use `ser.iloc[pos]` val = getattr(super(), mtd)(*args, **kwargs)
In [16]:
Copied!
vr = VoronoiRegionalizer(seeds=au_seeds_gdf)
au_result_gdf = vr.transform(gdf=au_gdf)
vr = VoronoiRegionalizer(seeds=au_seeds_gdf)
au_result_gdf = vr.transform(gdf=au_gdf)
Generating spherical polygons: 0%| | 0/88 [00:00<?, ?it/s]
Generating spherical polygons: 2%|▏ | 2/88 [00:00<00:05, 14.72it/s]
Generating spherical polygons: 20%|██ | 18/88 [00:00<00:01, 67.69it/s]
Generating spherical polygons: 32%|███▏ | 28/88 [00:00<00:00, 79.08it/s]
Generating spherical polygons: 52%|█████▏ | 46/88 [00:00<00:00, 95.83it/s]
Generating spherical polygons: 64%|██████▎ | 56/88 [00:00<00:00, 68.07it/s]
Generating spherical polygons: 73%|███████▎ | 64/88 [00:01<00:00, 45.65it/s]
Generating spherical polygons: 80%|███████▉ | 70/88 [00:01<00:00, 47.74it/s]
Generating spherical polygons: 86%|████████▋ | 76/88 [00:01<00:00, 44.33it/s]
Generating spherical polygons: 93%|█████████▎| 82/88 [00:01<00:00, 33.83it/s]
Generating spherical polygons: 100%|██████████| 88/88 [00:01<00:00, 51.06it/s]
Interpolating edges: 0%| | 0/341 [00:00<?, ?it/s]
Interpolating edges: 0%| | 1/341 [00:00<00:51, 6.54it/s]
Interpolating edges: 45%|████▌ | 154/341 [00:00<00:00, 739.14it/s]
Interpolating edges: 92%|█████████▏| 313/341 [00:00<00:00, 1082.79it/s]
Interpolating edges: 100%|██████████| 341/341 [00:00<00:00, 910.02it/s]
Generating polygons: 0%| | 0/88 [00:00<?, ?it/s]
Generating polygons: 70%|███████ | 62/88 [00:00<00:00, 569.21it/s]
Generating polygons: 100%|██████████| 88/88 [00:00<00:00, 513.20it/s]
In [17]:
Copied!
folium_map = plot_regions(au_result_gdf, tiles_style="CartoDB positron")
au_seeds_gdf.explore(
m=folium_map,
style_kwds=dict(color="#444", opacity=1, fillColor="#f2f2f2", fillOpacity=1),
marker_kwds=dict(radius=3),
)
folium_map = plot_regions(au_result_gdf, tiles_style="CartoDB positron")
au_seeds_gdf.explore(
m=folium_map,
style_kwds=dict(color="#444", opacity=1, fillColor="#f2f2f2", fillOpacity=1),
marker_kwds=dict(radius=3),
)
Out[17]:
Make this Notebook Trusted to load map: File -> Trust Notebook
In [18]:
Copied!
neighbourhood = AdjacencyNeighbourhood(regions_gdf=au_result_gdf)
neighbourhood = AdjacencyNeighbourhood(regions_gdf=au_result_gdf)
Nearest neighbours¶
In [19]:
Copied!
region_id = 0
neighbours = neighbourhood.get_neighbours(region_id)
neighbours
region_id = 0
neighbours = neighbourhood.get_neighbours(region_id)
neighbours
Out[19]:
{3, 12, 24, 26, 83}
In [20]:
Copied!
plot_neighbours(au_result_gdf, region_id, neighbours)
plot_neighbours(au_result_gdf, region_id, neighbours)
Out[20]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Neighbours at a distance¶
In [21]:
Copied!
region_id = 0
neighbours = neighbourhood.get_neighbours_at_distance(region_id, 3)
neighbours
region_id = 0
neighbours = neighbourhood.get_neighbours_at_distance(region_id, 3)
neighbours
Out[21]:
{6, 7, 13, 15, 19, 23, 33, 40, 45, 50, 59, 63, 64, 67, 70, 85, 87}
In [22]:
Copied!
plot_neighbours(au_result_gdf, region_id, neighbours)
plot_neighbours(au_result_gdf, region_id, neighbours)
Out[22]:
Make this Notebook Trusted to load map: File -> Trust Notebook
Plotting all neighbourhood¶
In [23]:
Copied!
region_id = 0
plot_all_neighbourhood(au_result_gdf, region_id, neighbourhood)
region_id = 0
plot_all_neighbourhood(au_result_gdf, region_id, neighbourhood)
Out[23]:
Make this Notebook Trusted to load map: File -> Trust Notebook