Classifying Mixed-Use Overlaps with GeoPandas Overlay

Parcels often straddle more than one land-use designation, and deciding whether a lot counts as single-use or genuinely mixed-use depends on how its area splits across those designations. This guide intersects parcels with a land-use zone layer using a GeoPandas overlay, measures the area share each designation claims within every parcel, then classifies each lot as single-use or mixed-use against a configurable dominance threshold. The output is a per-parcel table that planning staff can sort, filter, and defend during entitlement review.

Prerequisites #

• Parcel polygons with a unique identifier for every lot
• A land-use zone layer carrying a use-category attribute
• A dominance threshold defining when a single category makes a lot single-use
• GeoPandas 1.0+, Shapely 2.0+, and PyProj 3.4+ installed
• A projected metric CRS selected for the area of interest

Step-by-step #

The walkthrough uses a Philadelphia-region dataset and projects to UTM Zone 18N (EPSG:26918), whose unit is the meter, so every area is reported in square meters. Substitute the CRS appropriate to your jurisdiction before running any area math.

Step 1: Project both layers to a metric CRS #

Area calculations are only valid in a projected system, so reproject the parcels and the land-use layer to a shared metric frame and repair geometry before intersecting.

import geopandas as gpd

METRIC_CRS = "EPSG:26918"  # UTM Zone 18N, units = meters

parcels = gpd.read_file("parcels.gpkg").to_crs(METRIC_CRS)
landuse = gpd.read_file("landuse_zones.gpkg").to_crs(METRIC_CRS)

# Repair invalid rings so the overlay does not raise a topology error
parcels["geometry"] = parcels.geometry.make_valid()
landuse["geometry"] = landuse.geometry.make_valid()

# Record each parcel's full area for later share calculations
parcels["parcel_area"] = parcels.geometry.area

Step 2: Intersect parcels with land-use zones #

The overlay splits each parcel wherever a land-use boundary crosses it, producing one row per parcel-and-category fragment. keep_geom_type discards stray points or lines that can appear along shared edges.

# Geometric intersection yields one fragment per parcel/land-use pair
fragments = gpd.overlay(
    parcels[["parcel_id", "parcel_area", "geometry"]],
    landuse[["use_category", "geometry"]],
    how="intersection",
    keep_geom_type=True,
)

# Area of each fragment, still in square meters thanks to the metric CRS
fragments["frag_area"] = fragments.geometry.area

Step 3: Compute each category’s area share per parcel #

Aggregate the fragment areas by parcel and category, then divide by the parcel area to express every designation as a percentage of the lot.

import pandas as pd

# Sum fragment areas for each parcel/category combination
shares = (
    fragments.groupby(["parcel_id", "use_category"])
    .agg(cat_area=("frag_area", "sum"), parcel_area=("parcel_area", "first"))
    .reset_index()
)

# Share of the parcel occupied by each land-use category
shares["share_pct"] = (shares["cat_area"] / shares["parcel_area"] * 100).round(2)

Step 4: Classify each parcel as single-use or mixed-use #

A parcel is single-use when one category clears the dominance threshold; otherwise the mix is meaningful enough to label it mixed-use. Capturing the dominant category and its share keeps the decision transparent.

DOMINANCE_PCT = 85.0  # a lot is single-use if one category exceeds this

# For each parcel, find the largest category share
top = shares.sort_values("share_pct", ascending=False).drop_duplicates("parcel_id")

top["classification"] = (
    top["share_pct"].ge(DOMINANCE_PCT).map({True: "SINGLE_USE", False: "MIXED_USE"})
)
result = top[["parcel_id", "use_category", "share_pct", "classification"]]
result = result.rename(columns={"use_category": "dominant_category"})

Verification #

Confirm that the shares are internally consistent before publishing the classification. Every parcel’s category shares should sum to roughly 100 percent, allowing small slivers of rounding, and no share should exceed 100.

totals = shares.groupby("parcel_id")["share_pct"].sum()
assert (totals.between(98, 102)).mean() > 0.99  # allow minor sliver loss

# Confirm the working layer is metric so areas are square meters
assert fragments.crs.axis_info[0].unit_name == "metre"
print(result["classification"].value_counts())

Spot-check a handful of parcels flagged as mixed-use against an aerial or the source zoning map. If shares consistently fall short of 100 percent, the land-use layer likely has gaps, which is worth flagging separately.

Common Pitfalls #

• Measuring area in a geographic CRS. Computing geometry.area on EPSG:4326 returns square degrees, making every share meaningless. Project to a metric CRS first and assert the unit as the verification step does.
• Sliver fragments inflating category counts. Tiny overlay slivers along shared boundaries can register spurious categories. Drop fragments below a minimum area, for example frag_area < 1.0, before aggregating shares.
• Gaps in the land-use layer. If zones do not fully tile the study area, parcel shares will not sum to 100 and lots may be misclassified. Validate coverage, or treat the unclassified remainder as its own category.

Frequently Asked Questions #

Why use an overlay instead of a spatial join for this task? #

A join tells you which zones a parcel touches but not how much of the parcel each one covers. An overlay physically splits the parcel along zone boundaries so you can measure area shares, which is exactly what a mixed-use determination requires. This is the core technique behind land use intersection mapping.

How do I choose the dominance threshold? #

The threshold encodes local policy, so set it from the ordinance rather than by intuition. A common convention treats a lot as single-use when one category covers 85 to 90 percent of it, but jurisdictions with fine-grained mixed-use goals may lower that bar. Keep the value in a config file so planners can tune it without editing code.

Can I report every category rather than only the dominant one? #

Yes. The per-parcel shares table already holds every category, so you can pivot it to a wide format or emit the full breakdown alongside the classification. Retaining the complete distribution strengthens the audit trail when a determination is challenged.

What happens when a parcel splits evenly across three uses? #

No category clears the threshold, so the lot is labeled mixed-use, which is the intended outcome. If you need to distinguish two-way from three-way mixes, count the categories above a minimum share and add that count as a secondary attribute.

Choosing the Overlay Operation #

overlay supports several set operations and they answer different questions, so picking the wrong one produces a plausible frame with the wrong rows in it.

intersection keeps only the overlapping parts and is what a classification needs: one row per parcel-district pair, carrying attributes from both. union keeps everything from both layers including the non-overlapping remainder, which is useful for finding gaps and rarely what a classification wants. identity keeps all of the left layer plus the attributes of the right where they overlap, which is the right choice when every parcel must appear in the output even if it matches no district — often exactly what you want, because a parcel matching no district is a finding rather than a row to drop. difference and symmetric_difference are the tools for auditing coverage rather than for classifying.

The practical recommendation is identity for classification work, precisely because it preserves the unmatched parcels that intersection silently removes. A run whose parcel count drops between input and output has lost something, and losing it quietly is how an unzoned parcel ends up with no verdicts and no explanation.

Attribute Collisions and Where They Bite #

Overlaying two layers that both have a name column, or both have area, produces suffixed columns whose meaning is no longer obvious three functions later.

Rename before overlaying rather than untangling afterwards. Prefixing each layer’s attributes with its role — parcel_id, parcel_area, district_code, district_name — costs one line and makes every downstream reference self-documenting. It also removes the specific bug where an area column carried through from the input is mistaken for the area of the intersected part, which is a different number and the one the share calculation actually needs.

Recompute areas after the overlay rather than trusting any that came in with the inputs. The intersected geometry is new, and its area is the only one that reflects it; an input area attribute describes the whole parcel and will produce shares that exceed one when summed.

Validating the Classification Output #

Three checks catch nearly everything that goes wrong in an overlay-based classification, and all three are cheap.

Parcel conservation: every input parcel appears in the output at least once, and none appears with a total share exceeding one. Share coverage: the per-parcel shares sum to approximately one, with any shortfall attributable to a recorded gap in the district layer. And geometry type stability: the output contains only polygons, with any line or point fragments — produced where two boundaries are exactly coincident — filtered explicitly and counted rather than left to surprise a later area calculation.

def check_overlay(parcels, parts, id_col="parcel_id", tol=0.01):
    """Conservation, coverage and geometry-type checks in one pass."""
    shares = parts.groupby(id_col)["share"].sum()
    missing = set(parcels[id_col]) - set(shares.index)
    over = shares[shares > 1 + tol]
    under = shares[shares < 1 - tol]
    non_poly = (~parts.geom_type.isin(["Polygon", "MultiPolygon"])).sum()
    return {"missing_parcels": sorted(missing), "over_covered": over.index.tolist(),
            "under_covered": under.index.tolist(), "non_polygon_fragments": int(non_poly)}

Run it every time. The output is usually empty, and the run where it is not is the one that would otherwise have shipped a quietly wrong classification into a compliance report.

Related #

Part of: Land use intersection mapping

• Apportioning split-zoned parcels by area — using these shares to apportion a standard.
• Deciding which parcels a rule applies to — turning shares into an applicability decision.
• Geometry validation and topology repair — removing the defects that create sliver fragments.
• Optimizing spatial joins for 100k parcel datasets — making the overlay fast at county scale.

Conclusion #

Intersecting parcels with land-use zones, measuring each category’s area share, and comparing the dominant share to a configurable threshold produces a defensible mixed-use classification for an entire jurisdiction in one pass. Because the shares table preserves the full distribution, the same output supports downstream products such as density heatmaps from parcel centroids. For deeper coverage of overlay semantics, consult the GeoPandas set-operations guide and the Shapely geometry manual.