Access S2S via Database#
from typing import Dict
import folium as flm # Comment out if you have not installed, or run pip install folium matplotlib mapclassify
import geopandas as gpd
import numpy as np
import pandas as pd
from dotenv import load_dotenv
from geojson_pydantic import Feature, Polygon
from lonboard import Map, ScatterplotLayer
from shapely import from_geojson
from space2stats import StatsTable
Expected environment variables to be set from dotenv file
PGHOST=
PGPORT=
PGDATABASE=
PGUSER=
PGPASSWORD=
PGTABLENAME=space2stats
load_dotenv("../../db.env")
True
Connect to the database and print available fields
with StatsTable.connect() as stats_table:
fields = stats_table.fields()
fields
['pop_flood_pct',
'ogc_fid',
'sum_pop_f_0_2020',
'sum_pop_f_10_2020',
'sum_pop_f_15_2020',
'sum_pop_f_1_2020',
'sum_pop_f_20_2020',
'sum_pop_f_25_2020',
'sum_pop_f_30_2020',
'sum_pop_f_35_2020',
'sum_pop_f_40_2020',
'sum_pop_f_45_2020',
'sum_pop_f_50_2020',
'sum_pop_f_55_2020',
'sum_pop_f_5_2020',
'sum_pop_f_60_2020',
'sum_pop_f_65_2020',
'sum_pop_f_70_2020',
'sum_pop_f_75_2020',
'sum_pop_f_80_2020',
'sum_pop_m_0_2020',
'sum_pop_m_10_2020',
'sum_pop_m_15_2020',
'sum_pop_m_1_2020',
'sum_pop_m_20_2020',
'sum_pop_m_25_2020',
'sum_pop_m_30_2020',
'sum_pop_m_35_2020',
'sum_pop_m_40_2020',
'sum_pop_m_45_2020',
'sum_pop_m_50_2020',
'sum_pop_m_55_2020',
'sum_pop_m_5_2020',
'sum_pop_m_60_2020',
'sum_pop_m_65_2020',
'sum_pop_m_70_2020',
'sum_pop_m_75_2020',
'sum_pop_m_80_2020',
'sum_pop_m_2020',
'sum_pop_f_2020',
'sum_pop_2020',
'pop',
'pop_flood']
Define a geometry, in this case we use a bounding box surrounding Kenya and verify the GeoJSON Feature
AOIModel = Feature[Polygon, Dict]
# ~burundi
minx, miny, maxx, maxy = 29.038924, -4.468958, 30.850461, -2.310523
aoi = {
"type": "Feature",
"geometry": {
"type": "Polygon",
"coordinates": [
[
[minx, maxy],
[minx, miny],
[maxx, miny],
[maxx, maxy],
[minx, maxy],
]
],
},
"properties": {"name": "Updated AOI"},
}
feat = AOIModel(**aoi)
Connect to the database and retrieve hexagons that intersect with the bounding box
with StatsTable.connect() as stats_table:
data = stats_table.summaries(
aoi=feat,
spatial_join_method="touches",
fields=["sum_pop_2020"],
geometry="point",
)
df = pd.DataFrame(data)
df.head()
| hex_id | geometry | sum_pop_2020 | |
|---|---|---|---|
| 0 | 866ad8087ffffff | {"type":"Point","coordinates":[30.081627708207... | 25345.874701 |
| 1 | 866ad808fffffff | {"type":"Point","coordinates":[30.139814038919... | 11870.360712 |
| 2 | 866ad8097ffffff | {"type":"Point","coordinates":[30.032345934224... | 19034.332476 |
| 3 | 866ad809fffffff | {"type":"Point","coordinates":[30.090538449149... | 14700.551092 |
| 4 | 866ad80c7ffffff | {"type":"Point","coordinates":[30.247261402176... | 12067.935215 |
Retrieve all fields and polygon geometries for hexagons
with StatsTable.connect() as stats_table:
data = stats_table.summaries(
aoi=feat, spatial_join_method="touches", fields=fields, geometry="polygon"
)
df = pd.DataFrame(data)
df.head()
| hex_id | geometry | pop_flood_pct | ogc_fid | sum_pop_f_0_2020 | sum_pop_f_10_2020 | sum_pop_f_15_2020 | sum_pop_f_1_2020 | sum_pop_f_20_2020 | sum_pop_f_25_2020 | ... | sum_pop_m_60_2020 | sum_pop_m_65_2020 | sum_pop_m_70_2020 | sum_pop_m_75_2020 | sum_pop_m_80_2020 | sum_pop_m_2020 | sum_pop_f_2020 | sum_pop_2020 | pop | pop_flood | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 866ad8087ffffff | {"type":"Polygon","coordinates":[[[30.10399360... | 0.050492 | 6194996 | 318.054565 | 1747.088623 | 1923.321899 | 1251.022461 | 990.854004 | 1019.086853 | ... | 206.265076 | 122.173080 | 48.124065 | 19.658331 | 24.460270 | 12000.383184 | 13345.491516 | 25345.874701 | 26163.988 | 1321.075323 |
| 1 | 866ad808fffffff | {"type":"Polygon","coordinates":[[[30.16218780... | 0.107920 | 6194997 | 146.099442 | 825.470764 | 822.941284 | 574.661377 | 440.004395 | 454.970947 | ... | 122.356934 | 70.530174 | 25.380550 | 14.198433 | 14.668016 | 5801.703337 | 6068.657375 | 11870.360712 | 12280.807 | 1325.350181 |
| 2 | 866ad8097ffffff | {"type":"Polygon","coordinates":[[[30.05471035... | 0.039941 | 6194998 | 276.175537 | 1284.143066 | 1353.945312 | 1069.982666 | 753.500305 | 765.620483 | ... | 140.353241 | 89.866302 | 38.299133 | 19.963312 | 19.103374 | 9014.848974 | 10019.483501 | 19034.332476 | 18063.918 | 721.491072 |
| 3 | 866ad809fffffff | {"type":"Polygon","coordinates":[[[30.11291078... | 0.035242 | 6194999 | 190.571701 | 1002.784302 | 1023.486633 | 749.586670 | 556.957886 | 575.871033 | ... | 142.405670 | 81.206604 | 30.017090 | 15.441221 | 16.896353 | 7089.448963 | 7611.102129 | 14700.551092 | 14131.070 | 498.011870 |
| 4 | 866ad80c7ffffff | {"type":"Polygon","coordinates":[[[30.26964435... | 0.069827 | 6195003 | 169.925201 | 832.059998 | 751.781555 | 668.376648 | 389.565125 | 483.499878 | ... | 82.641251 | 48.474693 | 22.666998 | 11.139693 | 12.547176 | 5887.773416 | 6180.161799 | 12067.935215 | 12149.723 | 848.372900 |
5 rows × 45 columns
Convert to geopandas geodataframe
df["geometry"] = df["geometry"].apply(lambda geom: from_geojson(geom))
gdf = gpd.GeoDataFrame(df, geometry="geometry", crs="EPSG:4326")
gdf.head()
| hex_id | geometry | pop_flood_pct | ogc_fid | sum_pop_f_0_2020 | sum_pop_f_10_2020 | sum_pop_f_15_2020 | sum_pop_f_1_2020 | sum_pop_f_20_2020 | sum_pop_f_25_2020 | ... | sum_pop_m_60_2020 | sum_pop_m_65_2020 | sum_pop_m_70_2020 | sum_pop_m_75_2020 | sum_pop_m_80_2020 | sum_pop_m_2020 | sum_pop_f_2020 | sum_pop_2020 | pop | pop_flood | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 866ad8087ffffff | POLYGON ((30.10399 -2.3257, 30.11745 -2.29278,... | 0.050492 | 6194996 | 318.054565 | 1747.088623 | 1923.321899 | 1251.022461 | 990.854004 | 1019.086853 | ... | 206.265076 | 122.173080 | 48.124065 | 19.658331 | 24.460270 | 12000.383184 | 13345.491516 | 25345.874701 | 26163.988 | 1321.075323 |
| 1 | 866ad808fffffff | POLYGON ((30.16219 -2.34876, 30.17563 -2.31584... | 0.107920 | 6194997 | 146.099442 | 825.470764 | 822.941284 | 574.661377 | 440.004395 | 454.970947 | ... | 122.356934 | 70.530174 | 25.380550 | 14.198433 | 14.668016 | 5801.703337 | 6068.657375 | 11870.360712 | 12280.807 | 1325.350181 |
| 2 | 866ad8097ffffff | POLYGON ((30.05471 -2.36355, 30.06817 -2.33062... | 0.039941 | 6194998 | 276.175537 | 1284.143066 | 1353.945312 | 1069.982666 | 753.500305 | 765.620483 | ... | 140.353241 | 89.866302 | 38.299133 | 19.963312 | 19.103374 | 9014.848974 | 10019.483501 | 19034.332476 | 18063.918 | 721.491072 |
| 3 | 866ad809fffffff | POLYGON ((30.11291 -2.38662, 30.12636 -2.35369... | 0.035242 | 6194999 | 190.571701 | 1002.784302 | 1023.486633 | 749.586670 | 556.957886 | 575.871033 | ... | 142.405670 | 81.206604 | 30.017090 | 15.441221 | 16.896353 | 7089.448963 | 7611.102129 | 14700.551092 | 14131.070 | 498.011870 |
| 4 | 866ad80c7ffffff | POLYGON ((30.26964 -2.33398, 30.28307 -2.30106... | 0.069827 | 6195003 | 169.925201 | 832.059998 | 751.781555 | 668.376648 | 389.565125 | 483.499878 | ... | 82.641251 | 48.474693 | 22.666998 | 11.139693 | 12.547176 | 5887.773416 | 6180.161799 | 12067.935215 | 12149.723 | 848.372900 |
5 rows × 45 columns
m = gdf.explore(
column="sum_pop_2020",
tooltip="sum_pop_2020",
cmap="YlGnBu",
legend=True,
scheme="naturalbreaks",
legend_kwds=dict(colorbar=True, caption="Population", interval=False),
style_kwds=dict(weight=0, fillOpacity=0.8),
name="Population by Hexagon",
)
flm.LayerControl("topright", collapsed=False).add_to(m)
m
Make this Notebook Trusted to load map: File -> Trust Notebook
A faster interactive map can be created using lonboard
gdf.loc[:, "geometry"] = gdf.geometry.representative_point()
# Define custom breaks and corresponding RGBA colors
breaks = [
gdf["sum_pop_2020"].min(),
1,
1000,
10000,
50000,
100000,
200000,
gdf["sum_pop_2020"].max(),
]
colors = np.array(
[
[211, 211, 211, 255], # Light gray for 0
[255, 255, 0, 255], # Yellow for 1-1000
[255, 165, 0, 255], # Orange for 1000-10000
[255, 0, 0, 255], # Red for 10000-50000
[128, 0, 128, 255], # Purple for 50000-100000
[0, 0, 255, 255], # Blue for 100000-200000
[0, 0, 139, 255], # Dark blue for 200000+
]
)
# Function to assign colors based on custom bins
def assign_color(value, breaks, colors):
for i in range(len(breaks) - 1):
if breaks[i] <= value < breaks[i + 1]:
return colors[i]
return colors[-1] # In case value exceeds all breaks
# Map sum_pop_2020 values to colors using the custom function
gdf["color"] = gdf["sum_pop_2020"].apply(lambda x: assign_color(x, breaks, colors))
colors = np.uint8(gdf["color"].tolist())
# Create the scatterplot layer with the assigned colors
layer = ScatterplotLayer.from_geopandas(gdf, get_radius=2000, get_fill_color=colors)
m = Map(layer)
m