⚙️ How-To Biodiversity & species data Geometry & topology Biodiversity & Species Data · step 13 Notebook By · Updated

How to draw a species range from occurrence points

On this page

Problem statement

A range map turns scattered points into a polygon, and the polygon is used for things the points cannot support: area of occupancy, IUCN threat assessment, overlap with protected areas, gap analysis. Every method of drawing it makes a different assumption about the space between the points.

A convex hull assumes the species occupies everything between its outermost records, including the ocean. A grid of occupied cells assumes it occupies nothing it was not recorded in. An alpha hull sits between them, controlled by a parameter that changes the area by a factor of several.

Quick answer

Compute the two IUCN measures and look at the gap between them:

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

lynx = search(taxonKey=2435263, hasCoordinate="true")
clean = lynx.dropna(subset=["decimalLongitude", "decimalLatitude"])
clean = clean[clean.coordinateUncertaintyInMeters.fillna(0) <= 10_000]
pts = gpd.GeoDataFrame(clean, crs=4326,
                       geometry=gpd.points_from_xy(clean.decimalLongitude,
                                                   clean.decimalLatitude)).to_crs(6933)   # equal area

eoo = MultiPoint(list(pts.geometry)).convex_hull            # extent of occurrence
cells = {(int(x // 2000), int(y // 2000))
         for x, y in zip(pts.geometry.x, pts.geometry.y)}   # IUCN's 2 km grid
aoo = len(cells) * 4                                        # km² per occupied cell

print(f"{len(pts):,} usable records")
print(f"extent of occurrence (convex hull): {eoo.area / 1e6:>12,.0f} km²")
print(f"area of occupancy (2 km cells):     {aoo:>12,.0f} km²")
print(f"ratio: {eoo.area / 1e6 / aoo:,.0f}:1")
4,423 usable records
extent of occurrence (convex hull):  120,427,581 km²
area of occupancy (2 km cells):            7,640 km²
ratio: 15,763:1

Those two numbers are both standard, and 15,763 times apart. A convex hull over uncleaned records covers a quarter of the planet, because a few outliers on the far side of the world stretch it; the occupancy measure cannot move by more than one cell per record. A range map has to say which of the two it is drawing, and a gap this large says the records need more cleaning first.

The same occurrence points enclosed by a convex hull, an alpha hull and occupied grid cells. three polygons from one point setconvex hullthe whole outer envelopealpha hullconcavities up to alphaoccupied cellsonly what was recordedExtent of occurrence 120 million km²; area of occupancy 7,640 km² spatialworkflow.io
Three polygons from one point set; the difference between them is an assumption.

Decide the purpose, clean and thin, choose alpha openly

1. Decide what the polygon is for

Extent of occurrence for a threat assessment, occupied cells for an atlas, a modelled surface for habitat. They are not interchangeable.

2. Clean and thin first

A single mislocated record can add a million square kilometres to a convex hull. Range polygons are maximally sensitive to outliers.

3. Work in an equal-area projection

Areas in degrees are meaningless, and a species spanning 30° of latitude is badly distorted in Web Mercator.

4. Choose alpha deliberately, and show the sensitivity

An alpha hull's shape is controlled by one number. Publishing a polygon without it is publishing a shape nobody can reproduce.

5. Clip to land, and to the species' habitat

A marine gap in a terrestrial range is an artefact of the hull, not a distribution.

6. Report the record count and the date range

A range polygon from 200 records over a century is a different object from one from 5,000 records in five years.

Alpha hull area at four values of alpha. one parameter, eleven-fold in areaalpha 500 km9.4 M km², 1 partalpha 200 km5.9 M km², 2 partsalpha 100 km2.6 M km², 24 partsalpha 50 km0.8 M km², 70 parts05M km²10M km²convex hull: 120 M km²At 50 km the polygon has stopped being a range and become the fieldwork spatialworkflow.io
One parameter, four areas: alpha is the range map's biggest single decision.

Alpha hulls and the effect of one outlier

The data used on this page

Occurrence records are pulled from the GBIF occurrence search API by the code on this page. gbif.search is the paging function from How to download GBIF occurrences in Python, saved as a module by the setup. GBIF reprocesses continuously, so the setup answers each request with the reply it got on 9 October 2026; a live pull differs in the last digits.

import hashlib
import io
import os
import zipfile

import geopandas as gpd
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", 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"]


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])
''')

if not os.path.exists("bio/ne_land.gpkg"):                 # Natural Earth 1:110m land, the land mask
    blob = requests.get("https://naciscdn.org/naturalearth/110m/physical/ne_110m_land.zip",
                        timeout=120).content
    zipfile.ZipFile(io.BytesIO(blob)).extractall("bio/ne_land")
    gpd.read_file("bio/ne_land/ne_110m_land.shp").to_file("bio/ne_land.gpkg", driver="GPKG")
print("gbif.py written")
gbif.py written

Example 1 — an alpha hull, and how much alpha matters

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

def alpha_shape(points_xy, alpha_m):
    """Union of Delaunay triangles whose circumradius is below alpha."""
    tri = Delaunay(points_xy)
    keep = []
    for simplex in tri.simplices:
        a, b, c = points_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
        circum = la * lb * lc / (4 * area)
        if circum < alpha_m:
            keep.append(box(0, 0, 1, 1))                    # placeholder, replaced below
            keep[-1] = MultiPoint([tuple(a), tuple(b), tuple(c)]).convex_hull
    return unary_union(keep)

xy = np.c_[pts.geometry.x.to_numpy(), pts.geometry.y.to_numpy()]
xy = np.unique(np.round(xy / 1000) * 1000, axis=0)          # 1 km dedupe for the triangulation
print(f"{len(xy):,} distinct 1 km positions")
print(f"convex hull area: {MultiPoint([tuple(p) for p in xy]).convex_hull.area / 1e6:,.0f} km²")
for alpha_km in (500, 200, 100, 50):
    shape = alpha_shape(xy, alpha_km * 1000)
    print(f"  alpha {alpha_km:>3} km: {shape.area / 1e6:>10,.0f} km², "
          f"{len(shape.geoms) if shape.geom_type == 'MultiPolygon' else 1} part(s)")
1,982 distinct 1 km positions
convex hull area: 120,432,268 km²
  alpha 500 km:  9,357,839 km², 1 part(s)
  alpha 200 km:  5,928,029 km², 2 part(s)
  alpha 100 km:  2,615,032 km², 24 part(s)
  alpha  50 km:    821,761 km², 70 part(s)

Alpha moves the area by a factor of eleven and the polygon from one piece to seventy. There is no value that is correct: 500 km draws a blanket over the whole boreal zone, 50 km draws the sampling. Somewhere in between is a range map, and which one gets published is a judgement that has to be stated as one.

Example 2 — what one outlier does

hull_all = MultiPoint([tuple(p) for p in xy]).convex_hull
centre = np.array([xy[:, 0].mean(), xy[:, 1].mean()])
far = np.argsort(-np.linalg.norm(xy - centre, axis=1))

print(f"all {len(xy):,} positions: {hull_all.area / 1e6:>12,.0f} km²")
for drop in (1, 5, 10, 25):
    kept = np.delete(xy, far[:drop], axis=0)
    shrunk = MultiPoint([tuple(p) for p in kept]).convex_hull
    print(f"  drop the {drop:>2} most distant: {shrunk.area / 1e6:>10,.0f} km² "
          f"({shrunk.area / hull_all.area - 1:>6.1%})")
all 1,982 positions:  120,432,268 km²
  drop the  1 most distant: 63,836,925 km² (-47.0%)
  drop the  5 most distant: 19,191,158 km² (-84.1%)
  drop the 10 most distant: 19,149,554 km² (-84.1%)
  drop the 25 most distant: 18,869,626 km² (-84.3%)

The convex hull is decided by a handful of records on its edge, which are exactly the records most likely to be vagrants, escapes or mislocations. That sensitivity is why IUCN assessments pair it with the area of occupancy, which no single record can move by more than one cell.

Example 3 — clipping to land

land = gpd.read_file("bio/ne_land.gpkg").to_crs(6933)
hull_frame = gpd.GeoSeries([eoo], crs=6933)
nearby = land[land.intersects(eoo)]
clipped = eoo.intersection(nearby.union_all())
print(f"land polygons touching the hull: {len(nearby)}")
print(f"convex hull:        {eoo.area / 1e6:>12,.0f} km²")
print(f"clipped to land:    {clipped.area / 1e6:>12,.0f} km²")
print(f"sea removed:        {(eoo.area - clipped.area) / 1e6:>12,.0f} km² "
      f"({1 - clipped.area / eoo.area:.0%})")
print(f"area of occupancy:  {aoo:>12,.0f} km² "
      f"({aoo / (clipped.area / 1e6):.2%} of the clipped hull)")
land polygons touching the hull: 48
convex hull:         120,427,581 km²
clipped to land:      72,406,028 km²
sea removed:          48,021,553 km² (40%)
area of occupancy:         7,640 km² (0.01% of the clipped hull)

Two fifths of the convex hull is ocean, and after clipping it still covers 72 million km² — because a handful of outlying records stretch it across the northern hemisphere. The area of occupancy is 0.01% of that. Neither number is wrong; the pair of them is a warning that this point set has not been cleaned enough to draw a range from, which is the most useful thing the calculation can tell you.

Extent and occupancy both

Why extent and occupancy are both needed

Extent of occurrence measures how widely a species is spread, which is what matters for its exposure to a spatially correlated threat. Area of occupancy measures how much it actually occupies, which is what matters for how much habitat exists. A species can have a continental extent and a tiny occupancy — a specialist on scattered outcrops — and the pair of numbers says that where either alone does not.

Why alpha hulls need their parameter published

An alpha hull is a convex hull with holes and concavities allowed up to a scale of alpha. As alpha falls the polygon tightens onto the points, and the area falls continuously; there is no natural value. Reporting a polygon without alpha is like reporting a smoothed curve without the bandwidth.

Why the grid size in occupancy is not arbitrary

IUCN specifies a 2 km grid precisely because the measure is otherwise unstable: occupancy computed at 1 km and at 10 km differ by an order of magnitude for the same records, and the criterion thresholds are calibrated to the standard cell. A locally sensible cell size is scientifically reasonable and not comparable with anything.

Convex hull area as the most distant records are removed. the hull is decided by a handful of recordsall 1,982 positions120 M km²−1 recorddrop 1 most distant64 M (−47%)−4 more÷ 3.3drop 5 most distant19 M (−84%)−20 moredrop 25 most distant19 MAfter the first five, removing twenty more changes almost nothing spatialworkflow.io
Dropping the five most distant records removes 84% of the extent.

Dateline ranges, seasons and introduced populations

  • Dateline-crossing ranges break a planar hull; use an appropriate projection or split the polygon.
  • Seasonal ranges — breeding, wintering, passage — are different polygons and should not be merged.
  • Introduced populations inflate the extent of a native range; filter by establishment means.
  • establishmentMeans in GBIF distinguishes native, introduced and vagrant records where publishers supply it.
  • Alpha shapes can produce holes; whether a hole is real habitat or a sampling gap is not decidable from the points.
  • IUCN criteria use EOO and AOO with specific thresholds; the calculation has a standard and should follow it.

FAQ

How do I draw a species range from points?

Choose between a convex hull, an alpha hull and occupied grid cells, and say which — they differ by orders of magnitude.

What is the difference between EOO and AOO?

Extent of occurrence is the area within the outer boundary of the records; area of occupancy is the area of the 2 km cells they fall in.

What alpha should I use?

There is no natural value. Show the area at several and publish the one you used.

Why is my hull covering the ocean?

Because a convex hull fills everything between the outer records. Clip it to land or to suitable habitat.

How much can one record change the range?

A convex hull can move by hundreds of thousands of square kilometres on a single distant record.

Which measure does IUCN use?

Both, with different thresholds, and with a specified 2 km grid for occupancy.