Species Distribution Models Explained

On this page

Problem statement

A species distribution model takes occurrence records and environmental layers and produces a map of suitability. The map looks like a prediction of where the species lives. It is a statement about a comparison: how the environment at the records differs from the environment at a set of points the analyst declared to be available.

Everything that makes these models useful and everything that makes them misleading follows from that. There are no absences, so "available" is asserted. The records are biased, so the comparison inherits the bias. And the output is on an arbitrary scale, so calling it a probability is a category error.

Quick answer

The parts, and which of them you choose rather than measure:

parts = [
    ("presences", "occurrence records, cleaned and thinned", "measured, with bias"),
    ("background", "points defining what was available", "asserted by you"),
    ("covariates", "climate, terrain, land cover at both", "measured, at a chosen grain"),
    ("feature shape", "linear, quadratic, hinge, tree splits", "chosen, and decisive"),
    ("algorithm", "logistic, MaxEnt, boosting, ensemble", "chosen, least important"),
    ("evaluation", "how folds are split", "chosen, and usually wrong"),
]
width = max(len(p[0]) for p in parts)
for name, what, status in parts:
    print(f"{name:<{width}}  {what:<44}  {status}")

print("\nmeasured: 2 of 6")
print("on the same lynx data, changing only the feature shape moved AUC "
      "from 0.54 to 0.93;\nchanging only the background extent moved it "
      "from 0.77 to 0.52")
presences      occurrence records, cleaned and thinned       measured, with bias
background     points defining what was available            asserted by you
covariates     climate, terrain, land cover at both          measured, at a chosen grain
feature shape  linear, quadratic, hinge, tree splits         chosen, and decisive
algorithm      logistic, MaxEnt, boosting, ensemble          chosen, least important
evaluation     how folds are split                           chosen, and usually wrong

measured: 2 of 6
on the same lynx data, changing only the feature shape moved AUC from 0.54 to 0.93;
changing only the background extent moved it from 0.77 to 0.52

Two of the six inputs are data. The other four are decisions, and the published map rarely records them.

The inputs to a species distribution model, separating measured data from decisions. two inputs are data; four are decisionsinputwhat it ismeasured or chosenpresencesoccurrence records, cleaned andthinnedmeasured, with biasbackgroundpoints defining what was availablechosen β€” decisivecovariatesclimate and terrain at bothmeasured, at a chosen grainfeature shapelinear, quadratic, hinge, treeschosen β€” decisivealgorithmlogistic, MaxEnt, boostingchosen β€” least importantevaluationhow the folds are splitchosen β€” usually wrongsame lynx data: feature shape moved AUC 0.54 to 0.93; background extent 0.77 to 0.52The map cannot tell a reader which of its inputs were decisions spatialworkflow.io
Two inputs are measured; four are chosen, and the map cannot tell you which is which.

The question, the accessible area, mechanistic covariates

1. Say what question the model answers

Where the species could live (potential), where it does live (realised), or where it might live later (projection). The same code answers all three and means different things.

2. Define the accessible area explicitly

Background points represent "available". Available to what, over what time, across which barriers? This is the largest single decision in the workflow.

3. Choose covariates with a mechanism in mind

A variable is defensible if there is a reason the species should respond to it. Correlated climate layers picked by a stepwise procedure produce a model that predicts well and explains nothing.

4. Let the model express an optimum

Species have environmental optima. Linear terms cannot represent one; quadratics, hinges and trees can.

5. Evaluate on held-out space, not held-out rows

Spatially blocked cross-validation is the only evaluation that answers "will this work somewhere I have not sampled".

6. Report the relative suitability as relative

Prevalence is not identifiable from presence-only data. Any statement of the form "70% probability of occurrence" is unsupported.

Three questions a distribution model can answer and what each needs. one method, three questionswhere could it live?where does it live?where might it livelater?backgroundthe reachable regionthe occupied regionthe reachable regionpresencesany periodrecent, at equilibriumrecentsupportspotential, relativerealised, relative!needs anextrapolation mapNone of the three supports a probability of occurrence spatialworkflow.io
One method, three questions; they need different backgrounds and support different claims.

No absences, two backgrounds, one variable

Example 1 β€” why there are no absences

import numpy as np

# A species present in 30% of an area, detected with probability p per visit.
rng = np.random.default_rng(0)
truth = rng.random(10_000) < 0.30
for p, visits in [(0.2, 1), (0.2, 3), (0.5, 1), (0.5, 3), (0.9, 1)]:
    detected = truth & (rng.random(10_000) < 1 - (1 - p) ** visits)
    seen_rate = detected.mean()
    print(f"detection {p:>4.1f} Γ— {visits} visit(s): "
          f"{seen_rate:.1%} of cells record a presence, "
          f"{(truth & ~detected).mean():.1%} are present and unrecorded")
print("\na cell with no record is a cell with no record")
detection  0.2 Γ— 1 visit(s): 5.8% of cells record a presence, 24.3% are present and unrecorded
detection  0.2 Γ— 3 visit(s): 14.6% of cells record a presence, 15.5% are present and unrecorded
detection  0.5 Γ— 1 visit(s): 14.6% of cells record a presence, 15.4% are present and unrecorded
detection  0.5 Γ— 3 visit(s): 26.4% of cells record a presence, 3.7% are present and unrecorded
detection  0.9 Γ— 1 visit(s): 26.7% of cells record a presence, 3.4% are present and unrecorded

a cell with no record is a cell with no record

Example 2 β€” what the background does to the comparison

# Two backgrounds over the same presences, in one dimension for clarity.
presence_temp = rng.normal(-2.0, 3.0, 2_000)              # the species' realised niche
continent = rng.normal(8.0, 12.0, 10_000)                 # a continental background
region = rng.normal(-1.0, 5.0, 10_000)                    # its accessible region

def separability(a, b):
    """The AUC of a one-variable classifier, computed directly."""
    order = np.argsort(np.r_[a, b])
    labels = np.r_[np.ones(len(a)), np.zeros(len(b))][order]
    ranks = np.arange(1, len(labels) + 1)
    r_pos = ranks[labels == 1].sum()
    auc = (r_pos - len(a) * (len(a) + 1) / 2) / (len(a) * len(b))
    return max(auc, 1 - auc)          # the discriminating direction

print(f"against a continental background: AUC {separability(presence_temp, continent):.3f}")
print(f"against a regional background:    AUC {separability(presence_temp, region):.3f}")
print("the species did not change; the question did")
against a continental background: AUC 0.790
against a regional background:    AUC 0.566
the species did not change; the question did

One variable, one set of presences, two backgrounds: 0.79 against a continental sample and 0.57 against the region the species could actually reach. The first number would be reported as a good model. It measures the fact that the species is not in the tropics.

Example 3 β€” the scale of the output

suitability = rng.beta(2, 5, 10_000)                       # a typical model surface
for prevalence in (0.05, 0.2, 0.5):
    scaled = suitability * prevalence / suitability.mean()
    print(f"if true prevalence is {prevalence:>4.0%}: "
          f"mean predicted occurrence {scaled.mean():.3f}, "
          f"max {min(scaled.max(), 1):.3f}")
print("\nthe model output is identical in all three cases β€”")
print("presence-only data contains no information about which is right")
if true prevalence is   5%: mean predicted occurrence 0.050, max 0.153
if true prevalence is  20%: mean predicted occurrence 0.200, max 0.614
if true prevalence is  50%: mean predicted occurrence 0.500, max 1.000

the model output is identical in all three cases β€”
presence-only data contains no information about which is right

The same suitability surface is consistent with a species occupying 5% or 50% of the landscape. Thresholding it into "suitable" and "unsuitable" imports a prevalence assumption without stating it.

Niches, autocorrelation and transferability

Why "niche" needs qualifying

The fundamental niche is where the species could persist physiologically. The realised niche is where it does persist given competitors, predators and history. A model fitted on occurrence records estimates something narrower than both: the environments where the species was recorded, which adds detectability and observer effort to the list. Calling the output a niche model is a claim about which of those three it estimates.

Why spatial autocorrelation breaks the evaluation

Environmental layers are smooth and records are clustered, so a record's neighbours share both its covariates and its label. Random cross-validation therefore tests the model on data it has effectively seen, and returns a score that a spatially blocked split will not reproduce. The gap between the two is a measure of how much of the model is memorisation.

Why transferability is the hard part

A model fitted in one region and applied to another, or fitted on today's climate and projected to 2070, is extrapolating: the covariate combinations it is asked about were not in the training data. All the flexible methods that fit well within their training range β€” boosting especially β€” extrapolate badly and confidently. That is why projections should be accompanied by a map of where the model is extrapolating at all.

Fundamental, realised and recorded niches as nested sets. occurrence records sample the smallest setfundamental nichewhere it could persist physiologicallyrealised nicheminus competitors, predators, historyrecorded nicheminus detectability and efforta model fitted on records estimates the inner set and is read as the outer one spatialworkflow.io
Three nested sets; occurrence records sample the smallest one.

Ensembles, sample size and logistic output

  • Ensembles average over algorithms, which averages over the least important decision.
  • Sample size: below about 30 presences, models are unstable and simple envelopes are competitive.
  • Non-equilibrium species β€” invasives, recovering populations β€” violate the assumption that the records show where the species can live.
  • Biotic interactions are absent from the covariates and often dominate the real distribution.
  • MaxEnt's logistic output was reinterpreted as "cloglog" in later versions; both are relative, not probabilities.
  • Presence-absence data, where it exists, supports far stronger claims β€” use it instead when you have it.

FAQ

What is a species distribution model?

A model that compares the environment where a species was recorded with the environment available to it, and maps the resulting contrast.

Why is there no absence data?

Occurrence records document presences. "Not recorded" conflates absence with nobody looking.

What does the output mean?

Relative suitability. Presence-only data cannot identify the probability of occurrence.

Which algorithm is best?

The question matters less than the background and the feature shape. Most comparisons of algorithms are comparisons of flexibility.

How should I evaluate one?

With spatially blocked folds, and by reporting the drop from the random-fold score.

Can I project one into the future?

With care, and with a map of where the projection is extrapolating beyond the training data.