Note
Go to the end to download the full example code
Graphs from geographic points#
This example shows how to build a graph from a set of points using PySAL and geopandas. In this example, we’ll use the famous set of cholera cases at the Broad Street Pump, recorded by John Snow in 1853. The methods shown here can also work directly with polygonal data using their centroids as representative points.
Traceback (most recent call last):
File "/builddir/build/BUILD/networkx-networkx-3.2.1/examples/geospatial/plot_points.py", line 54, in <module>
add_basemap(facet)
File "/usr/lib/python3.12/site-packages/contextily/plotting.py", line 134, in add_basemap
image, extent = bounds2img(
^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/contextily/tile.py", line 262, in bounds2img
arrays = Parallel(n_jobs=n_connections, prefer=preferred_backend)(
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/joblib/parallel.py", line 1918, in __call__
return output if self.return_generator else list(output)
^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/joblib/parallel.py", line 1847, in _get_sequential_output
res = func(*args, **kwargs)
^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/joblib/memory.py", line 573, in __call__
return self._cached_call(args, kwargs, shelving=False)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/joblib/memory.py", line 530, in _cached_call
return self._call(call_id, args, kwargs, shelving)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/joblib/memory.py", line 762, in _call
output = self.func(*args, **kwargs)
^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/contextily/tile.py", line 291, in _fetch_tile
request = _retryer(tile_url, wait, max_retries)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/contextily/tile.py", line 434, in _retryer
request = requests.get(tile_url, headers={"user-agent": USER_AGENT})
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/requests/api.py", line 73, in get
return request("get", url, params=params, **kwargs)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/requests/api.py", line 59, in request
return session.request(method=method, url=url, **kwargs)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/requests/sessions.py", line 589, in request
resp = self.send(prep, **send_kwargs)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/requests/sessions.py", line 703, in send
r = adapter.send(request, **kwargs)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
File "/usr/lib/python3.12/site-packages/requests/adapters.py", line 519, in send
raise ConnectionError(e, request=request)
requests.exceptions.ConnectionError: HTTPSConnectionPool(host='a.tile.openstreetmap.fr', port=443): Max retries exceeded with url: /hot/17/65484/43579.png (Caused by NewConnectionError('<urllib3.connection.HTTPSConnection object at 0xffff57f23530>: Failed to establish a new connection: [Errno -3] Temporary failure in name resolution'))
from libpysal import weights, examples
from contextily import add_basemap
import matplotlib.pyplot as plt
import networkx as nx
import numpy as np
import geopandas
# read in example data from a geopackage file. Geopackages
# are a format for storing geographic data that is backed
# by sqlite. geopandas reads data relying on the fiona package,
# providing a high-level pandas-style interface to geographic data.
cases = geopandas.read_file("cholera_cases.gpkg")
# construct the array of coordinates for the centroid
coordinates = np.column_stack((cases.geometry.x, cases.geometry.y))
# construct two different kinds of graphs:
## 3-nearest neighbor graph, meaning that points are connected
## to the three closest other points. This means every point
## will have exactly three neighbors.
knn3 = weights.KNN.from_dataframe(cases, k=3)
## The 50-meter distance band graph will connect all pairs of points
## that are within 50 meters from one another. This means that points
## may have different numbers of neighbors.
dist = weights.DistanceBand.from_array(coordinates, threshold=50)
# Then, we can convert the graph to networkx object using the
# .to_networkx() method.
knn_graph = knn3.to_networkx()
dist_graph = dist.to_networkx()
# To plot with networkx, we need to merge the nodes back to
# their positions in order to plot in networkx
positions = dict(zip(knn_graph.nodes, coordinates))
# plot with a nice basemap
f, ax = plt.subplots(1, 2, figsize=(8, 4))
for i, facet in enumerate(ax):
cases.plot(marker=".", color="orangered", ax=facet)
add_basemap(facet)
facet.set_title(("KNN-3", "50-meter Distance Band")[i])
facet.axis("off")
nx.draw(knn_graph, positions, ax=ax[0], node_size=5, node_color="b")
nx.draw(dist_graph, positions, ax=ax[1], node_size=5, node_color="b")
plt.show()
Total running time of the script: (0 minutes 0.169 seconds)