How to use Geopandas to count intersecting features where there is a matching ID

Question

Is it possible to use sjoin() or another method to get the counts of interesecting features between two geodataframes where the two geodataframes have a common id? For example, i would like to add a column to a polygon geodataframe showing the counts of points within each polygon where the point has a matching identifier to the polygon.

Answer 1

Spatial join them, drop non matching on the common field, count and join/merge the result back to the input polygon df:

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
import geopandas as gpd

poly = gpd.read_file(r"/home/bera/GIS/data/Indexrutor/index_100km.shp")
point = gpd.read_file(r"/home/bera/Desktop/GIStest/pnts.geojson")

#Both dataframes have the field BK which is the field match them by
matchfield = "BK"

#Spatial join them and keep the fields BK_left, BK_right
sj = gpd.sjoin(left_df=poly, right_df=point, how="left", predicate="intersects")[[f"{matchfield}_left", f"{matchfield}_right"]]

 #Select the rows where BK match
sj = sj.loc[sj[f"{matchfield}_left"]==sj[f"{matchfield}_right"]]

#Polygons with multiple points in them are duplicated because of left join, so we can count unique BK:s.
pointcount = sj.groupby(f"{matchfield}_left").size().rename("pointcount")
# BK_left
# 63_1    1
# 64_2    1
# 69_3    2
# ...

#Join the pointcount Series back to the polygon dataframe
poly = gpd.pd.merge(left=poly, right=pointcount, left_on=matchfield, right_index=True, how="left")
poly["pointcount"] = poly["pointcount"].fillna(0) #Polygons with no point should be 0 and not NaN
poly.to_file(r"/home/bera/Desktop/GIStest/poly_with_pointcount.geojson")