Geodata DareData Blog

count_hospitals_in_1km.py

# Public Hospitals in Lisbon
hospitals_url = "https://opendata.arcgis.com/datasets/172678f193144512860a397fde991361_4.geojson" # GeoJSON
hospitals_gdf = gpd.read_file(hospitals_url).to_crs(epsg=3857)
hospitals_gdf.head()

# Buffer the house locations by 1km
house_data_gdf_buffer = (
   house_data_gdf
   .copy()
   .assign(geometry_buffer = lambda d: d.buffer(1000))
   .set_geometry("geometry_buffer")
)
 
# To get the a statistic of a geometry within another geometry, 
# we can use a spatial join and then aggregate the values
house_data_gdf["hospitals_in_1km"] = (
    gpd.sjoin(
        house_data_gdf_buffer,
        hospitals_gdf,
        how="left",
        op="contains"
   )
   .groupby("house_id", as_index=False)
   .agg({"OBJECTID": "count"})
   .rename(columns={"OBJECTID":"Hospitals"})
   .loc[:, "Hospitals"]
)