Data systems

tabulai/

mlx-geospatial

Point-in-polygon, spatial joins and haversine distances on Apple silicon GPUs through MLX, with an sjoin for GeoDataFrames.

What’s new here

This library moves the point-against-polygon part onto the Mac GPU using MLX. Results can stay as MLX arrays, so zone labels and distance features can feed MLX operations directly.

What it does

It covers a narrow set of operations: point-in-polygon containment, spatial joins of points to polygons, and haversine distances (pairwise, all-pairs and k-nearest). Polygons with holes and multipolygons work. Polygon overlays and polygon-to-polygon predicates do not.

You can pass NumPy arrays, Shapely geometries or GeoDataFrames. For repeated queries you build a PolygonSet and a GridIndex once and query many times. sjoin supports inner and left joins, and all_matches=True returns every overlapping polygon instead of the first.

The GPU path uses float32 in local coordinates, so points near polygon edges can land on either side. With host inputs, exact_boundary=True adds a CPU check to match Shapely’s behavior, at the cost of some host work.

Who it’s for

Anyone on an Apple silicon Mac who assigns many GPS points to zones, finds nearby places, or builds spatial features for an MLX model. If you run on Linux or Windows, or need overlays and polygon-to-polygon predicates, this is not the tool. The README itself notes that a prebuilt BallTree can be faster for nearest-neighbor queries against large reference sets.

Try it

It needs an Apple silicon Mac, Python 3.10+ and MLX 0.32.0+. NumPy and Shapely install automatically.

pip install mlx-geospatial
# for GeoPandas support
pip install "mlx-geospatial[geo]"

Then find which polygon contains each point (coordinates are longitude, latitude):

import numpy as np
from shapely.geometry import box
import mlx_geospatial as ms
 
points = np.array([
    [-74.00, 40.72],
    [-73.98, 40.75],
    [-73.90, 40.80],
])
zones = [
    box(-74.02, 40.70, -73.99, 40.74),
    box(-74.00, 40.73, -73.96, 40.78),
]
 
labels = ms.point_in_polygon(points, zones)
print(labels)  # [ 0  1 -1]

Each label is a polygon’s position in zones, and -1 means no match.

How mature is it

Young and active. The repo was created on 2026-08-29, has 12 stars, 1 fork, 2 contributors and 2 releases, the latest being v0.1.1 on 2026-10-20. It had 60 commits in the last 90 days, the last push was 2026-10-20, and there are no open issues or pull requests. The license is MIT.

The README’s benchmarks are the author’s own and are synthetic. On an M1 Ultra and an M4 Max, end-to-end spatial joins ran 3.78 to 18.66 times faster than native GeoPandas on one machine and 2.55 to 16.51 times faster on the other, across six workloads. The README says they are not a guarantee for other datasets, and it reports 721 passing tests on each machine.