The solution you found will work, but it is indeed very inefficient. The first statement creates one huge, super complex multipolygon from your big list of (simple) polygons. The second statement will compare each single point with this huge polygon to see if it is inside, which is very costly.
A more efficient way to do this is to use geopandas.sindex.query. Under the hood, it will first use a spatial index to search (very quickly) which points are within the minimum bounding box of which of your original, so relatively simple and small, polygons. Only for those point-polygon combinations it will go on to do the actual check if the point is actually intersecting with the (relatively simple) polygon. A similar number of intersections might have to be calculated if most points are near polygons, but it is a comparison of a point with a small simple polygon instead of with a huge super complex one, which is a lot faster.
I use the "intersects" predicate instead of "within" because you stated your needs like this: which does not overlap/within/intersect with. A point lying on a border of a polygon will be treated as "intersects", but won't be "within".
The call looks like this:
points_in_poly_idx = polys.sindex.query(points.geometry, predicate="intersects")[0]
Once you have a list of points that intersect, you can use pandas indexing to find the points that are not in this list:
points_outside_polys = points.loc[~points.index.isin(points_in_poly_idx)]
Full sample script:
import geopandas as gpd
from matplotlib import pyplot as plt
from shapely import box, Point
# Create some test data
polys = gpd.GeoDataFrame(
geometry=[
box(xmin=xmin, ymin=ymin, xmax=xmin + 10, ymax=ymin + 10)
for xmin in range(0, 100, 10)
for ymin in range(0, 100, 10)
],
crs=31370,
)
points = gpd.GeoDataFrame(
geometry=[Point(50, 50), Point(55, 55), Point(200, 50)], crs=31370
)
# Find all intersections between points and polygons, using a spatial index for speed.
# The first result of `query` contains all integer indexes of all points where at least
# one intersection was found.
points_in_poly_idx = polys.sindex.query(points.geometry, predicate="intersects")[0]
# Exclude all points where an intersection was found.
points_outside_polys = points.loc[~points.index.isin(points_in_poly_idx)]
print(len(points_outside_polys))
# One point outside of the polygons
fig, ax = plt.subplots()
polys.plot(ax=ax, edgecolor="blue", facecolor="none")
points.plot(ax=ax, color="red")
points_outside_polys.plot(ax=ax, color="green")
plt.show()
Result:
