🔧 Fix Biodiversity & species data Geometry & topology Biodiversity & Species Data · step 20 Notebook By · Updated

A range polygon has the wrong shape

On this page

Problem statement

The range map covers the ocean, includes a country the species has never been recorded in, or breaks into seventy fragments that are obviously an artefact of where people walk. All three come from the same place: a polygon is being asked to say something about the space between the records, and the records do not contain that information.

On uncleaned Canada lynx records a convex hull covers 120 million km² — a quarter of the planet's surface — because five records are on the wrong side of the world. Removing them removes 84% of the area.

Quick answer

Find the records that are doing the damage before changing the method:

import numpy as np, geopandas as gpd
from shapely.geometry import MultiPoint
from gbif import search

lynx = search(taxonKey=2435263, hasCoordinate="true").dropna(
    subset=["decimalLongitude", "decimalLatitude"])
pts = gpd.GeoDataFrame(lynx, crs=4326,
                       geometry=gpd.points_from_xy(lynx.decimalLongitude,
                                                   lynx.decimalLatitude)).to_crs(6933)

hull = MultiPoint(list(pts.geometry)).convex_hull
centre = np.array([pts.geometry.x.median(), pts.geometry.y.median()])
far = np.linalg.norm(np.c_[pts.geometry.x, pts.geometry.y] - centre, axis=1) / 1000

print(f"{len(pts):,} records, convex hull {hull.area / 1e6:,.0f} km²")
print(f"distance from the median position: median {np.median(far):,.0f} km, "
      f"99th percentile {np.percentile(far, 99):,.0f} km, max {far.max():,.0f} km")
order = np.argsort(-far)
print("\nthe five most distant records:")
print(pts.iloc[order[:5]][["decimalLatitude", "decimalLongitude", "year",
                           "basisOfRecord", "countryCode"]].to_string(index=False))
6,353 records, convex hull 120,427,581 km²
distance from the median position: median 796 km, 99th percentile 8,081 km, max 30,867 km

the five most distant records:
 decimalLatitude  decimalLongitude   year      basisOfRecord countryCode
       52.121667        177.558611 2014.0 PRESERVED_SPECIMEN          US
       50.282200         18.993900 2006.0 PRESERVED_SPECIMEN          PL
        0.000000          0.000000    NaN PRESERVED_SPECIMEN         NaN
        0.000000          0.000000    NaN PRESERVED_SPECIMEN         NaN
        0.000000          0.000000    NaN PRESERVED_SPECIMEN         NaN

A hull is decided entirely by its outermost points, so the first question is always which records those are and whether they should be in the dataset at all.

Range polygon problems with the symptom and the fix. four wrong shapesWhat the polygon doesWhat is actually wrongit covers a quarter of the planetOne record across the dateline;check the outermost fiveit includes a country with onerecordA specimen catalogued at itsholding institutionit covers ocean between islandsA convex hull cannot express acoastline — clip or use alphait breaks into seventy fragmentsAlpha is tighter than thesampling; the map is the fieldwork spatialworkflow.io
Four wrong shapes, four different causes.

Outliers as records, filter, then use occupancy

1. Look at the outliers as records, not as geometry

Check their country, date, basis and dataset. A "lynx" in Europe is usually a zoo animal, a museum specimen catalogued at its holding institution, or a different species under an old name.

2. Filter on establishment and country before drawing

establishmentMeans, countryCode and the species' known native range remove most of the damage before any geometry runs.

3. Use occupancy for anything quantitative

Occupied 2 km cells cannot be moved by one record by more than one cell; a hull can be moved by thousands of kilometres.

4. Choose a hull that can have concavities

An alpha hull follows the records; a convex hull cannot express a coastline, a mountain barrier or a gap.

5. Clip to habitat, not just to land

A terrestrial range clipped to land still covers the sea of grass a forest species does not use. Clipping to modelled habitat is better where a model exists.

6. Show the sensitivity

Publish the area with and without the outermost records, and at more than one alpha. The spread is the honest uncertainty of the polygon.

Hull area as records are filtered on non-geometric criteria. three filters, none of them geometricall 6,353 records120 M km²Canada and the US only50 M — four records removed+ since 199040 M+ uncertainty <= 10 km34 M, 1,719 records062.5M km²125M km²Removing four records takes 58% of the area; the geometry never changed spatialworkflow.io
The first five records removed take 84% of the area with them.

What the outliers are, and filtering before drawing

The data used on this page

Records come from the GBIF API via gbif.search, the paging function from the download guide plus eventDate and countryCode. Climate is WorldClim 2.1 bioclim at 10′ in bio/wc/; pa.protected_areas() is British Columbia's protected-area layer (EPSG:3005). GBIF changes daily, so the setup answers each request with the reply it got on 9 October 2026.

import hashlib
import io
import os
import zipfile

import requests

# GBIF adds and revises records daily. So that the numbers below reproduce, every GBIF request on
# this page is answered from the reply it got on 9 October 2026 (one file per request, named by a
# SHA-1 of its URL); any other request, or all of them once these lines are deleted, goes live.
SAVED = zipfile.ZipFile(io.BytesIO(requests.get(
    "https://data.spatialworkflow.io/gbif-api-2026-10-09.zip", timeout=300).content))
live_get = requests.get


def saved_get(url, params=None, **kwargs):
    prepared = requests.Request("GET", url, params=params).prepare()
    try:
        body = SAVED.read(hashlib.sha1(prepared.url.encode()).hexdigest()[:16] + ".json")
    except KeyError:
        return live_get(url, params=params, **kwargs)
    reply = requests.models.Response()
    reply.status_code, reply._content, reply.url = 200, body, prepared.url
    reply.headers["Content-Type"] = "application/json"
    return reply


requests.get = saved_get

os.makedirs("bio/wc", exist_ok=True)
with open("gbif.py", "w") as f:
    f.write('''
import pandas as pd, requests

API = "https://api.gbif.org/v1"
MAX_OFFSET = 100_000
FIELDS = ["key", "species", "speciesKey", "decimalLongitude", "decimalLatitude",
          "coordinateUncertaintyInMeters", "year", "month", "basisOfRecord",
          "datasetKey", "institutionCode", "issues", "eventDate", "countryCode"]


def count(**params):
    return requests.get(f"{API}/occurrence/search", params={**params, "limit": 0}, timeout=120).json()["count"]


def search(limit=None, page=300, **params):
    total = count(**params)
    if limit is None and total > MAX_OFFSET:
        raise ValueError(f"{total:,} records exceeds the {MAX_OFFSET:,} offset limit")
    rows, offset = [], 0
    target = total if limit is None else min(limit, total)
    while offset < target:
        body = requests.get(f"{API}/occurrence/search",
                            params={**params, "limit": min(page, target - offset), "offset": offset},
                            timeout=180).json()
        rows += body["results"]
        offset += len(body["results"])
        if body["endOfRecords"] or not body["results"]:
            break
    return pd.DataFrame([{f: r.get(f) for f in FIELDS} for r in rows])
''')
with open("pa.py", "w") as f:
    f.write('''
import os
import geopandas as gpd

URL = ("https://services6.arcgis.com/ubm4tcTYICKBpist/arcgis/rest/services/"
       "British_Columbia_Parks_Ecological_Reserves_and_Protected_Areas/FeatureServer/0/query"
       "?where=1%3D1&outFields=PROTECTED_LANDS_NAME&outSR=3005&f=geojson")


def protected_areas():
    if not os.path.exists("bio/bc_protected.gpkg"):
        gpd.read_file(URL).to_file("bio/bc_protected.gpkg")
    return gpd.read_file("bio/bc_protected.gpkg")
''')
if not os.path.exists("bio/wc/wc2.1_10m_bio_1.tif"):
    blob = requests.get("https://geodata.ucdavis.edu/climate/worldclim/2_1/base/wc2.1_10m_bio.zip", timeout=600).content
    zipfile.ZipFile(io.BytesIO(blob)).extractall("bio/wc")
print(len([n for n in os.listdir("bio/wc") if n.endswith(".tif")]), "WorldClim layers")
19 WorldClim layers

Example 1 — what the outliers are

suspects = pts.iloc[order[:20]]
print(f"{'country':>8}{'year':>7}{'basis':>22}{'distance km':>13}")
for r, d in zip(suspects.itertuples(), far[order[:20]]):
    print(f"{str(r.countryCode):>8}{r.year if r.year == r.year else 0:>7.0f}"
          f"{r.basisOfRecord:>22}{d:>13,.0f}")

by_country = pts.countryCode.value_counts()
print(f"\nrecords by country: {by_country.head(6).to_dict()}")
print(f"countries with fewer than 5 records: "
      f"{(by_country < 5).sum()} of {len(by_country)}")
 country   year                 basis  distance km
      US   2014    PRESERVED_SPECIMEN       30,867
      PL   2006    PRESERVED_SPECIMEN       15,581
     nan      0    PRESERVED_SPECIMEN       15,163
     nan      0    PRESERVED_SPECIMEN       15,163
     nan      0    PRESERVED_SPECIMEN       15,163
      CA   2024     HUMAN_OBSERVATION        8,692
      CA   1915    PRESERVED_SPECIMEN        8,690
      CA   1915    PRESERVED_SPECIMEN        8,690
      CA   1956    PRESERVED_SPECIMEN        8,676
      CA   1956    PRESERVED_SPECIMEN        8,676
      CA   1956    PRESERVED_SPECIMEN        8,676
      CA   2025     HUMAN_OBSERVATION        8,660
      CA   2025     HUMAN_OBSERVATION        8,636
      CA   1956    PRESERVED_SPECIMEN        8,630
      CA   1956    PRESERVED_SPECIMEN        8,628
      CA   1956    PRESERVED_SPECIMEN        8,577
      CA   1956    PRESERVED_SPECIMEN        8,577
      CA   2026     HUMAN_OBSERVATION        8,577
      CA   2025     HUMAN_OBSERVATION        8,572
      CA   2022     HUMAN_OBSERVATION        8,507

records by country: {'US': 3816, 'CA': 2533, 'PL': 1}
countries with fewer than 5 records: 1 of 3

The list names its own causes. The most distant record is a preserved specimen in the Aleutians at 177°E — on the far side of the dateline, which is what makes the hull wrap the planet. The next is a single specimen in Poland, almost certainly catalogued at its holding institution rather than where it was collected. Then three records at exactly 0, 0. One country with one record, four records that are not observations of a wild animal in its range, and a quarter of the planet's surface.

Example 2 — filtering before drawing

def hull_km2(frame):
    return MultiPoint(list(frame.geometry)).convex_hull.area / 1e6

native = pts[pts.countryCode.isin(["CA", "US"])]
recent = native[native.year >= 1990]
precise = recent[recent.coordinateUncertaintyInMeters.fillna(0) <= 10_000]

for label, frame in [("all records", pts), ("Canada and the US", native),
                     ("+ since 1990", recent), ("+ uncertainty <= 10 km", precise)]:
    print(f"{label:<26}{len(frame):>7,} records{hull_km2(frame):>15,.0f} km²")
all records                 6,353 records    120,427,581 km²
Canada and the US           6,349 records     50,032,761 km²
+ since 1990                2,683 records     39,963,664 km²
+ uncertainty <= 10 km      1,719 records     33,894,034 km²

Three filters, each defensible on its own terms, take the hull from 120 million km² to 34 million — a factor of 3.6, and most of it in the first step, which removed four records. None of the three is a geometric operation, and together they do far more for the polygon than any choice of hull algorithm.

Example 3 — how many fragments an alpha hull should have

from scipy.spatial import Delaunay
from shapely.ops import unary_union

def alpha_shape(xy, alpha_m):
    tri = Delaunay(xy)
    keep = []
    for simplex in tri.simplices:
        a, b, c = xy[simplex]
        la, lb, lc = (np.linalg.norm(b - c), np.linalg.norm(a - c), np.linalg.norm(a - b))
        s = (la + lb + lc) / 2
        area = max(s * (s - la) * (s - lb) * (s - lc), 1e-12) ** 0.5
        if la * lb * lc / (4 * area) < alpha_m:
            keep.append(MultiPoint([tuple(a), tuple(b), tuple(c)]).convex_hull)
    return unary_union(keep)

xy = np.unique(np.round(np.c_[precise.geometry.x, precise.geometry.y] / 1000) * 1000,
               axis=0)
print(f"{len(xy):,} distinct 1 km positions after filtering")
for alpha_km in (400, 200, 100, 50, 25):
    shape = alpha_shape(xy, alpha_km * 1000)
    parts = len(shape.geoms) if shape.geom_type == "MultiPolygon" else 1
    largest = (max(g.area for g in shape.geoms) if parts > 1 else shape.area) / 1e6
    print(f"  alpha {alpha_km:>3} km: {shape.area / 1e6:>9,.0f} km², {parts:>3} part(s), "
          f"largest {largest / (shape.area / 1e6):.0%} of the area")
1,192 distinct 1 km positions after filtering
  alpha 400 km: 5,981,266 km²,   2 part(s), largest 100% of the area
  alpha 200 km: 3,716,572 km²,   6 part(s), largest 58% of the area
  alpha 100 km: 1,597,756 km²,  27 part(s), largest 29% of the area
  alpha  50 km:   502,370 km²,  58 part(s), largest 16% of the area
  alpha  25 km:   108,472 km²,  98 part(s), largest 23% of the area

As alpha tightens the polygon fragments, and the fragments are the sampling pattern rather than the species' range. The value worth publishing is usually the largest one whose shape still looks like a range, chosen by eye and declared as such.

One record can dominate a hull

Why one record can dominate

A convex hull is the smallest convex set containing the points, so its area is set by the extremes in every direction. Adding one record 5,000 km away adds a wedge whose area scales with the square of that distance. No amount of averaging, smoothing or cleaning of the other records changes it; only removing that one does.

Why fragments are not gaps

An alpha hull leaves a hole wherever the records are sparser than alpha. That happens both where the species is absent and where nobody looked, and the polygon cannot tell those apart. A fragmented range map is therefore a statement about the record distribution that will be read as a statement about the species.

Why a modelled surface is often the better answer

If the question is "where does this species live", a distribution model with environmental covariates uses information the points alone do not have, and its output degrades gracefully where records are sparse. Hulls remain the right tool for IUCN's extent of occurrence, which is defined as a hull and is supposed to be crude.

The same records under a wide and a tight alpha. alpha decides what the map is aboutalpha 400 km2 parts — a rangealpha 50 km70 parts — the fieldworkThe holes are where records are sparser than alpha: absence, or nobody looking spatialworkflow.io
Alpha decides whether the map shows a range or shows the fieldwork.

Zoo records, synonyms and migratory ranges

  • Dateline: a species spanning 180° gets a hull wrapping the wrong way in most projections; split it or use a polar projection.
  • Zoo and botanic garden records are a common outlier source; establishmentMeans and institutionCode catch many.
  • Old synonyms pull in records of another species; check the taxonomy before the geometry.
  • Migratory species need seasonal ranges; one polygon over all months is meaningless.
  • Convex hull area is not IUCN EOO unless computed on the right projection and the right records; the criteria have a specification.
  • Holes that are real — a lake, a mountain block — are worth keeping; holes that are sampling gaps are not.

FAQ

Why does my range polygon cover the ocean?

A convex hull fills the space between the outer records. Clip it to land, or use an alpha hull.

One record is making the range enormous — what do I do?

Check it as a record: country, date, basis, dataset. It is usually a zoo animal, a synonym or a georeferencing error.

Should I use a convex or an alpha hull?

A convex hull for IUCN extent of occurrence; an alpha hull, with its parameter published, for anything meant to look like a range.

Why does the alpha hull break into fragments?

Because the records are sparser than alpha in those places — which may be absence or may be nobody looking.

How do I report the uncertainty of a range polygon?

Give the area with and without the outermost records, and at more than one alpha.

Is a modelled range better?

For "where does it live", usually yes. For a criterion defined as a hull, no.