Create a choropleth map: a practical guide to classification, color, and geographic data

A choropleth map assigns colors to geographic areas according to an associated numeric or ordinal value. This guide shows how to build one with GeoPandas and Matplotlib, choose an appropriate classification scheme, handle missing values and coordinate reference systems, design a readable legend, and validate that the map communicates the data rather than exaggerating it.

Choropleth maps are useful when the geography itself is the unit of comparison: counties, census tracts, municipalities, watersheds, service areas, districts, or other polygons. Each polygon receives a color based on a value such as population density, percentage of households, incident rate, permit count, or average travel time. GeoPandas supports choropleths through the `column`, `cmap`, `legend`, `scheme`, and related plotting arguments, while Matplotlib controls color normalization and figure output. (geopandas.org) The important work happens before the final plotting call. You must decide whether the variable should be mapped as a count, rate, percentage, or index; confirm that the values and geometries refer to the same geographic units and time period; choose a classification method that fits the distribution; and make missing or suppressed values visible. A polished map is not automatically an accurate map.

How EDK fits

For: GIS analysts, Python developers, data journalists, researchers, public-sector data teams, and anyone who needs to communicate a geographic comparison from polygon data.

A choropleth map is one way to turn structured geographic observations into an interpretable view of the physical world. The map does not merely display geometry: it joins a measurement to places, establishes a spatial comparison, and lets a reader inspect geographic patterns such as where service requests are concentrated, which districts have higher construction activity, or how infrastructure conditions vary between neighborhoods. In the broader world-model perspective, choropleths are a presentation layer over geographic computation. The underlying workflow may require finding the relevant boundaries, matching observations to those boundaries, aggregating records, calculating rates, ranking regions, filtering by time, and comparing one geography with another. E2’s editorial connection is therefore genuine when the goal is to find, understand, or compute about real places. This guide does not claim that the external libraries described here provide those broader capabilities automatically; GeoPandas and Matplotlib handle the local data preparation and visualization workflow described below.

  • •Map the rate of sanitation, noise, tree-maintenance, parking, or building-related service requests by neighborhood or administrative area after aggregating point or address records into shared polygons.
  • •Compare physical-world features such as roads, streams, parks, schools, hydrants, charging stations, or construction permits by region using counts normalized by population, area, or another defensible denominator.
  • •Use a choropleth as the final view of a larger geographic computation: filter observations by date, join them to boundaries, calculate rates, rank regions, and expose the resulting pattern for inspection.
Primary library
GeoPandas, using Matplotlib for rendering
Best geometry type
Polygons or multipolygons representing meaningful comparison areas
Mapped value
Numeric or ordinal attribute associated with each polygon
Optional dependency
mapclassify for classification schemes such as quantiles and natural breaks
Common CRS warning
Area and distance calculations should use an appropriate projected CRS rather than geographic degrees

The EDK outcome

Before you start

Python 3 with GeoPandas and Matplotlib installed.

A polygon dataset containing the geographic units to compare.

A numeric or ordinal field, or a separate table that can be aggregated and joined to the polygons.

A stable geographic identifier for joins whenever one exists.

Knowledge of the value’s definition, observation period, denominator, and missing-data rules.

A known CRS for the geometry. If the data have no CRS, determine it from authoritative metadata before assigning one; `set_crs()` labels coordinates but does not transform them. (docs.geopandas.org)

Learning objectives

Distinguish a choropleth from a categorical polygon map and from a point-density map.

Prepare one value per polygon and identify when counts should be converted to rates or densities.

Choose between continuous color mapping, quantiles, equal intervals, natural breaks, and explicit user-defined classes.

Create a reproducible choropleth with GeoPandas and Matplotlib.

Validate CRS, joins, legends, missing values, and interpretation before publishing.

Verification

The polygon count is unchanged after joins unless multipart features were intentionally dissolved or exploded.

The mapped field is numeric or explicitly converted to numeric, and its null, zero, minimum, maximum, and median values have been reviewed.

The legend shows the same class boundaries used by the renderer.

A rate or density includes its denominator and unit in the title, caption, legend, or accompanying documentation.

The CRS is present, and all area or distance calculations were performed in an appropriate projected CRS.

The map has been compared with at least one alternate classification or an unclassified continuous version.

Missing and suppressed values are visually distinguishable from true low values and true zeros.

The output file opens correctly and retains readable labels, boundaries, source notes, and attribution where required.

Troubleshooting

Blank or misplaced geometries: inspect `polygons.crs`, geometry bounds, and whether the layer contains empty or invalid geometries. Do not assign a CRS by guesswork; `set_crs()` labels existing coordinates, while `to_crs()` transforms them. (docs.geopandas.org)

The map is one solid color: check whether the mapped field is numeric, whether all values are identical, and whether a previous conversion turned values into strings. Print `polygons["rate_per_1000"].dtype` and `describe()`.

The legend does not match expectations: inspect the selected `scheme`, `k`, `vmin`, `vmax`, and any `classification_kwds`. A classified map and a continuous colorbar communicate different things.

The map exaggerates large rural areas: replace raw counts with a rate or density where the analytical question calls for comparison rather than total volume.

Areas disappear after a join: check the join key, identifier types, duplicate keys, and whether the table uses a different boundary vintage. Use `validate="one_to_one"` or another explicit merge validation where appropriate.

Rates look unstable: inspect denominator size, small-number effects, reporting coverage, and suppression rules. Consider a minimum denominator, confidence interval, smoothing method, or an uncertainty layer instead of presenting the rate without context.

Colors are difficult to distinguish: reduce the number of classes, use a palette designed for ordered data, and test print, grayscale, and color-vision-deficiency conditions. ColorBrewer provides sequential and diverging scheme guidance and accessibility-oriented selection options. (colorbrewer2.org)

Area values are implausible: check whether area was calculated while the geometries were still in a geographic CRS. GeoPandas documents that area should be calculated after transforming to a planar CRS with suitable units. (docs.geopandas.org)

Next steps

Add a second layer for context, such as roads, water, boundaries, or labeled places, while keeping the choropleth visually dominant.

Publish the computed polygon values as GeoJSON or another documented geographic format for interactive inspection.

Add time controls or produce a consistent series of maps using fixed class breaks so changes reflect the data rather than changing scales.

Move aggregation and spatial joins into PostGIS or another spatial data system when the source data becomes too large for an in-memory notebook workflow.

Document the geographic unit, denominator, boundary vintage, and classification method as part of the dataset rather than treating the map image as the only deliverable.

What you need before creating the map

You need a polygon layer and a value associated with each polygon. The value may already be present as an attribute, or it may be produced by an aggregation. For example, a boundary layer might contain one row per neighborhood while a separate table contains service requests. In that case, count or summarize the requests by neighborhood, then join the result to the polygons using a stable identifier. Do not join on display names alone when a formal geographic identifier is available.

Decide what the value means before choosing colors. A count answers “how many records occurred here?” A rate answers “how often did this occur relative to a population or exposure?” A density answers “how much occurred per unit of area?” These are different questions. Mapping raw construction-permit counts, for example, may emphasize large regions, while permits per square mile or per thousand residents may better support comparison.

  • •Polygon layer: GeoJSON, GeoPackage, Shapefile, PostGIS query, or another format supported by GeoPandas.
  • •Numeric or ordinal field: counts, rates, percentages, densities, scores, or categories with an explicit ordering.
  • •A known CRS for the geometry. GeoPandas describes the CRS as the information that relates coordinate values to locations on Earth. (docs.geopandas.org)
  • •Python with GeoPandas, Matplotlib, and optionally `mapclassify` for classified choropleths.

Choose the visual encoding and classification

For a continuous variable, GeoPandas can map values directly through the `column` argument. Its default color mapping is continuous and linear between the observed minimum and maximum. You can change the range with `vmin` and `vmax`, or use a classification scheme with `scheme` when discrete classes are more appropriate. GeoPandas documents schemes including quantiles, equal intervals, natural breaks, Fisher–Jenks, percentiles, standard-deviation methods, and user-defined breaks. (geopandas.org)

Use quantiles when you want each class to contain roughly the same number of geometries. This can make regional differences easy to see, but the numeric ranges may be uneven. Equal intervals divide the numeric range into equal-width bins, which makes the legend straightforward but can leave some classes nearly empty when the data are skewed. Natural-breaks methods can reflect clusters in the observed values, but the resulting thresholds may be less intuitive and can change substantially when the dataset changes.

  • •Sequential: counts, rates, percentages, density, elevation, or other low-to-high variables.
  • •Diverging: change around zero, difference from a target, or values above and below a meaningful benchmark.
  • •Qualitative: nominal categories with no numeric order; this is usually not a choropleth of a continuous variable.
  • •User-defined breaks: policy thresholds, risk bands, service-level limits, or other externally meaningful categories.

Ordered workflow with GeoPandas

The following example assumes a polygon file with a numeric field called `rate_per_1000`. Replace the path and field name with your own data. The example uses a classified map with a legend and a restrained base style.

For a basic map, `GeoDataFrame.plot()` is sufficient. For more controlled color behavior, Matplotlib normalization objects can map data values to a color scale. `BoundaryNorm` is useful when explicit boundaries should determine discrete color intervals. (matplotlib.org)

  • •Load the polygons and inspect the schema.
  • •Confirm the CRS and geometry type.
  • •Validate the mapped field and decide how to handle missing values.
  • •Select a palette and classification method.
  • •Render the map with an explicit legend, title, and source note.
  • •Save the figure at a deliberate resolution and review it at its intended display size.

Validation, interpretation, and publication

A choropleth should be validated as both a data product and a visual product. First, compare the number of polygons before and after any joins. Unexpected row multiplication often means that the join key was not unique. Next, compare summary statistics for the source field and the plotted field. Check the minimum, maximum, median, number of nulls, number of zeros, and number of polygons outside any manually selected display range.

Inspect the map alongside a table of the highest and lowest values. A visually dark polygon should correspond to a high value according to the legend, not simply to a large polygon. If the map displays a rate, include the denominator and rate definition in the title, subtitle, caption, or metadata. If some areas have very small denominators, consider suppression, uncertainty flags, or a minimum-count rule rather than presenting unstable rates as precise comparisons.

  • •Check that every polygon has at most one mapped value after joining.
  • •Verify that the legend labels match the actual class boundaries.
  • •Use a distinct treatment for missing values instead of assigning them the lowest numeric class.
  • •Test the palette in grayscale and with a color-vision-deficiency simulator when accessibility matters.
  • •Compare the final map with an alternate classification before making substantive claims about geographic patterns.

Workflow

01

1. Install the plotting and classification dependencies

Install GeoPandas, Matplotlib, and mapclassify. The last package is optional for unclassified continuous maps but is needed for the classification schemes exposed through GeoPandas.

python -m pip install geopandas matplotlib mapclassify
02

2. Load the polygon layer and inspect the schema

Start by checking geometry, CRS, column names, null counts, and basic statistics. Do not select a color scheme until you understand the mapped field.

import geopandas as gpd

polygons = gpd.read_file("data/areas.gpkg", layer="areas")

print(polygons.head())
print(polygons.dtypes)
print("CRS:", polygons.crs)
print("Geometry types:", polygons.geometry.geom_type.value_counts())
print(polygons["rate_per_1000"].describe())
print("Missing:", polygons["rate_per_1000"].isna().sum())
03

3. Prepare one defensible value per polygon

If the value is in another table, aggregate first and join on a stable identifier. Verify that the join did not multiply polygon rows. Keep the original count and denominator when deriving a rate so the result can be audited.

events_by_area = (
    events.groupby("area_id", as_index=False)
    .agg(event_count=("event_id", "count"), population=("population", "first"))
)
events_by_area["rate_per_1000"] = (
    events_by_area["event_count"] / events_by_area["population"] * 1000
)

polygons = polygons.merge(events_by_area, on="area_id", how="left", validate="one_to_one")
assert len(polygons) == polygons["area_id"].nunique()
04

4. Check the CRS before calculating area or distance

For display, retain a CRS appropriate to the output. For area or distance-derived denominators, transform to a suitable projected CRS first. GeoPandas warns that area in a geographic CRS can be invalid because the units are degrees rather than linear units. (docs.geopandas.org)

if polygons.crs is None:
    raise ValueError("The polygon layer has no CRS; determine it from source metadata first.")

# Example only: choose a projection appropriate to the study area.
projected = polygons.to_crs("EPSG:5070")
projected["area_sq_km"] = projected.geometry.area / 1_000_000
05

5. Render a first continuous choropleth

Create a baseline before experimenting with classification. This makes it easier to tell whether a later problem comes from the data or from the classing method.

import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(10, 8))
polygons.plot(
    column="rate_per_1000",
    cmap="viridis",
    legend=True,
    edgecolor="white",
    linewidth=0.3,
    missing_kwds={"color": "lightgrey", "label": "No data"},
    ax=ax,
)
ax.set_title("Service-request rate by area")
ax.set_axis_off()
plt.tight_layout()
plt.savefig("choropleth-continuous.png", dpi=200, bbox_inches="tight")
06

6. Compare classification schemes

Use the same field, palette, map extent, and title while comparing classing methods. This isolates the effect of classification. GeoPandas supports schemes through mapclassify and allows the number of classes to be controlled with `k`. (geopandas.org)

fig, axes = plt.subplots(2, 2, figsize=(14, 10))

for ax, scheme in zip(
    axes.flat,
    ["quantiles", "equal_interval", "natural_breaks", "fisher_jenks"],
):
    polygons.plot(
        column="rate_per_1000",
        scheme=scheme,
        k=5,
        cmap="YlOrRd",
        legend=True,
        edgecolor="white",
        linewidth=0.25,
        missing_kwds={"color": "lightgrey", "label": "No data"},
        ax=ax,
    )
    ax.set_title(scheme)
    ax.set_axis_off()

plt.tight_layout()
plt.savefig("choropleth-classification-comparison.png", dpi=200)
07

7. Use explicit breaks when the thresholds matter

If the map must reflect defined operational or policy bands, provide the breaks directly instead of relying on a data-driven classifier.

polygons.plot(
    column="rate_per_1000",
    scheme="UserDefined",
    classification_kwds={"bins": [5, 10, 25, 50]},
    cmap="YlOrRd",
    legend=True,
    edgecolor="white",
    linewidth=0.3,
    missing_kwds={"color": "lightgrey", "label": "No data"},
)
plt.title("Rate bands defined by reporting thresholds")
plt.axis("off")
plt.tight_layout()
plt.savefig("choropleth-user-defined.png", dpi=200)
08

8. Add metadata and perform a final review

Before publication, record the source, date range, boundary version, field definition, denominator, classification method, number of classes, palette, missing-data treatment, and CRS. Review the image at the size and medium where readers will encounter it.

metadata = {
    "value_field": "rate_per_1000",
    "definition": "events per 1,000 residents",
    "classification": "quantiles",
    "classes": 5,
    "palette": "YlOrRd",
    "missing_value_style": "lightgrey",
    "crs": str(polygons.crs),
}
print(metadata)

Questions teams ask

What data is best for a choropleth map?

Use a polygon layer with one defensible value per polygon. Numeric or ordinal variables work best: rates, percentages, densities, counts, indexes, or measured change. If your source data are points or events, aggregate them to the polygon layer first and document the aggregation. When polygon populations or areas differ, a normalized measure is usually more comparable than a raw count.

Should I use quantiles or equal intervals?

Neither is universally correct. Quantiles distribute observations across classes and can reveal relative ranking, but class ranges may be uneven. Equal intervals preserve simple numeric ranges but may produce empty or crowded classes for skewed data. Compare both against a histogram and the map’s intended question. Use user-defined breaks when an external standard or operational threshold is more important than statistical convenience. GeoPandas exposes these choices through the `scheme` argument when `mapclassify` is available. (geopandas.org)

Do I need to reproject the data before making a choropleth?

Not always for visual display, but you should understand the CRS. GeoPandas can plot geographic coordinates, and it can transform geometries with `to_crs()`. For area or distance calculations, use a suitable projected CRS because geographic coordinates are expressed in degrees and can produce invalid measurements for those operations. (docs.geopandas.org)

How should missing values appear?

Do not silently treat missing values as zero or as the lowest class. Keep missing, suppressed, invalid, and true-zero states distinct when the source supports that distinction. Render missing polygons with a neutral color or hatch, and explain the treatment in the legend or caption.

Can I create a choropleth from service requests, crimes, permits, or infrastructure records?

Yes, if the records can be reliably assigned to the target polygons and the resulting statistic is defined clearly. You might map requests per thousand residents, permits per square mile, or infrastructure features per district. The important steps are geocoding or spatially joining the records, choosing a denominator when appropriate, checking coverage and time period, and documenting uncertainty.

Keep exploring

Sources

More in Guides