Choropleth Mapping in Python: A Practical Guide to Geographic Data Science

Q: Which classification scheme should I use for skewed data?

For heavy-tailed or power-law distributed data such as income, population, or wealth, Head-Tail Breaks is designed specifically for this case, recursively partitioning around the mean. Quantiles is a solid general-purpose alternative.

Elizabeth Sramek Published 12 October 2026 Updated 25 September 2026 4 MIN READ

A choropleth map is a geographic map that encodes statistical values as color across predefined regions — countries, states, ZIP codes — and in Python, the standard toolkit for building one is GeoPandas for the geometry and PySAL’s mapclassify package for deciding how to bucket your values into color classes. That second part — classification — is where most choropleth maps I review actually go wrong, and it’s almost never a geometry problem.

The One-Line Map That Hides a Hard Decision

Making a choropleth in GeoPandas looks deceptively simple:

import geopandas as gpd

ax = mx.plot(
    column="PCGDP1940",
    scheme="Quantiles",
    cmap="YlGn",
    legend=True,
    legend_kwds={"fmt": "{:.0f}"}
)
ax.set_axis_off()

That single scheme="Quantiles" argument is doing all the analytical heavy lifting, and it’s the part beginners skip past without realizing they’ve made a choice. Change it, and the same data can tell visibly different stories on the map.

Classification Schemes: The Decision That Actually Matters

I keep a mental shortlist of when to reach for each classifier, built from redoing the same regional analysis with three different schemes early in my career and watching the “hot spots” move around the map each time.

import mapclassify

# Equal-width bins — simple, but sparse if data has outliers
ei5 = mapclassify.EqualInterval(mx["PCGDP1940"], k=5)

# Roughly equal count per class — avoids sparse bins
q5 = mapclassify.Quantiles(mx["PCGDP1940"], k=5)

# Best for roughly normal distributions
msd = mapclassify.StdMean(mx["PCGDP1940"])

# Good for skewed, heavy-tailed data (income, population)
ht = mapclassify.HeadTailBreaks(mx["PCGDP1940"])

# Optimal partitioning for a fixed k — my default when I need one clean answer
fj5 = mapclassify.FisherJenks(mx["PCGDP1940"], k=5)

My rule of thumb: for economic data like GDP or income, which is almost always right-skewed, I start with HeadTailBreaks or Quantiles and avoid EqualInterval entirely — equal-width bins on skewed data tend to dump 80% of regions into a single class and make the map look far more uniform than the underlying data actually is.

How to Compare Classifiers Objectively

Rather than eyeballing which map “looks right,” I compare the Absolute Deviation around Class Medians (ADCM) across classifiers — lower is a tighter fit:

import numpy

class5 = q5, ei5, ht, fj5
fits = numpy.array([c.adcm for c in class5])
for name, fit in zip(["Quantiles", "EqualInterval", "HeadTailBreaks", "FisherJenks"], fits):
    print(f"{name}: {fit:.1f}")

This turns “which classification looks better” from a subjective call into a number I can put in a methods appendix — which matters a great deal when a client asks why the map draws attention to one region over another.

Choosing a Color Palette to Match Your Classification

The palette family has to match the data type, the same rule that applies to any Seaborn chart — see our guide to Seaborn color palettes for the general version of this rule.

# Sequential — one direction, low to high
ax = mx.plot(column="PCGDP1940", scheme="Quantiles", cmap="Blues",
             edgecolor="k", linewidth=0.1, legend=True)

# Diverging — meaningful midpoint (e.g. growth rate above/below zero)
ax = mx.plot(column="growth_rate", scheme="Quantiles", cmap="RdYlGn", legend=True)

# Qualitative — categories with no inherent order
ax = mx.plot(column="region_type", categorical=True, legend=True)

Comparing Maps Across Years With Pooled Classification

A mistake I made on an early client project: classifying each year’s map independently, which meant the color scale shifted from year to year and made a genuinely flat trend look dramatic. The fix is pooled classification — one shared set of bins across all time periods:

years = ["PCGDP1940", "PCGDP1960", "PCGDP1980", "PCGDP2000"]
pooled = mapclassify.Pooled(mx[years].values, classifier="Quantiles", k=5)

fig, axs = plt.subplots(1, 4, figsize=(20, 5))
for i, y in enumerate(years):
    mx.plot(
        y,
        scheme="UserDefined",
        classification_kwds={"bins": pooled.global_classifier.bins},
        legend=True,
        ax=axs[i]
    )

If you’re building any kind of time-series choropleth — year-over-year sales by region, changing demographics — pooled classification isn’t optional. Without it, you’re not comparing the data across time; you’re comparing however each year happened to bucket independently.

Frequently Asked Questions

What’s the difference between quantile and equal interval classification?

Equal interval splits the full value range into equal-width bins, which can leave some classes nearly empty if the data has outliers. Quantile classification instead puts roughly the same number of observations in each class, which avoids sparse bins but can produce very different interval widths between classes.

Which classification scheme should I use for skewed data?

For heavy-tailed or power-law distributed data — income, population, wealth — Head-Tail Breaks is designed specifically for this case, recursively partitioning around the mean. Quantiles is a solid general-purpose alternative.

Why does my choropleth map look misleading?

The most common cause is a classification scheme mismatched to your data’s distribution — usually Equal Interval applied to skewed data, which compresses most regions into one class. Try Quantiles or Fisher-Jenks and compare the ADCM fit statistic before finalizing.