Checking a distance sampling model

distance sampling

model checking

abundance

ecology tutorial

The goodness-of-fit test accepts surveys whose density is 13 per cent wrong, and the assumption that matters most leaves no trace in the likelihood at all.

Author

Tidy Ecology

Published

2026-08-26

Three posts in this series have shown a distance sampling analysis returning a confident wrong number: a pooled detection function reversing which habitat is denser, a half-normal at a point inflating density by a sixth, a hierarchical model halving a population when the animals hold still.

The reasonable response is to ask what you would check. Distance sampling has a natural answer, because unlike most ecological models it makes a direct, testable prediction about data you can see: the histogram of distances should look like the fitted curve. So plot it, run a goodness-of-fit test, and if it passes, you are fine.

This post takes that programme seriously and measures how much protection it buys. The honest answer is: some, and a good deal less than the ritual implies.

library(ggplot2)

g_hn  <- function(x, s) exp(-x^2/(2*s^2))
mu_hn <- function(s, w) sqrt(2*pi)*s*(pnorm(w, 0, s) - 0.5)
ld_hn <- function(x, s, w) -x^2/(2*s^2) - log(mu_hn(s, w))
g_hr  <- function(x, s, b) 1 - exp(-(x/s)^(-b))
mu_hr <- function(s, b, w) integrate(g_hr, 0, w, s = s, b = b)$value
fit_hn <- function(x, w) exp(optimise(function(p) -sum(ld_hn(x, exp(p), w)),
                                      c(-8, 1), tol = 1e-10)$minimum)
w <- 0.1

The check that works

Start with the good news, because there is some. Bin the distances, fit the key function, compare observed counts against expected counts with a chi-square statistic. Before trusting a test to reject a bad model, check that it accepts a good one at the advertised rate: a test with the wrong size is not evidence, it is a coin with an opinion.

cut <- seq(0, w, length.out = 6)
pi_hn <- function(s) { P <- sqrt(2*pi)*s*(pnorm(cut[-1], 0, s) - pnorm(cut[-6], 0, s))/w
                       P/sum(P) }
gof_hn <- function(x) {
  s  <- fit_hn(x, w)
  ob <- tabulate(findInterval(x, cut, rightmost.closed = TRUE), 5)
  ex <- length(x)*pi_hn(s)
  c(p = pchisq(sum((ob - ex)^2/ex), df = 5 - 1 - 1, lower.tail = FALSE), sigma = s)
}
## a world where the half-normal is exactly right
s_ok <- 0.045; D0 <- 100; L_ok <- 266/(2*D0*mu_hn(s_ok, w))
cal <- t(sapply(1:400, function(i) {
  set.seed(61000 + i)
  N <- rpois(1, D0*2*w*L_ok); x <- runif(N, 0, w)
  gof_hn(x[runif(N) < g_hn(x, s_ok)])
}))
sprintf("model-true data, n approx 266: test rejects at 0.05 in %.1f%% of surveys | median p %.3f",
        100*mean(cal[, "p"] < 0.05), median(cal[, "p"]))

[1] "model-true data, n approx 266: test rejects at 0.05 in 6.0% of surveys | median p 0.493"

Well behaved: it rejects the truth at about the nominal rate, so it is a real test rather than a formality. And it does catch gross problems. A spike of animals recorded at zero distance, a histogram with a hole in it, distances heaped on round numbers: all of that shows up.

What it misses

Now give it the job it is actually asked to do. The true detection function has a shoulder, generated by a hazard-rate key. You fit a half-normal, which cannot make a shoulder. From the point transect post, we know this costs about 11 per cent on a line. Does the test tell you?

s_hr <- 0.05; b_hr <- 3
mu_t <- mu_hr(s_hr, b_hr, w); Ln <- 266/(2*D0*mu_t)
res <- t(sapply(1:400, function(i) {
  set.seed(63000 + i)
  N <- rpois(1, D0*2*w*Ln); x <- runif(N, 0, w)
  x <- x[runif(N) < g_hr(x, s_hr, b_hr)]
  g <- gof_hn(x)
  c(p = g["p"], bias = length(x)/(2*Ln*mu_hn(g["sigma"], w))/D0 - 1, n = length(x))
}))
acc <- res[res[, 1] > 0.05, ]
sprintf("half-normal fitted to hazard-rate data (n approx %.0f, 5 bins, df 3):", mean(res[, 3]))

[1] "half-normal fitted to hazard-rate data (n approx 268, 5 bins, df 3):"

sprintf("  rejects at 0.05: %.1f%% of surveys | median p-value %.3f | mean density bias %+.1f%%",
        100*mean(res[, 1] < 0.05), median(res[, 1]), 100*mean(res[, 2]))

[1] "  rejects at 0.05: 17.2% of surveys | median p-value 0.312 | mean density bias +12.4%"

sprintf("  among the surveys the test ACCEPTS: mean bias %+.1f%%, worst %+.1f%%",
        100*mean(acc[, 2]), 100*max(acc[, 2]))

[1] "  among the surveys the test ACCEPTS: mean bias +12.2%, worst +37.8%"

sprintf("  correlation between the p-value and |bias|: %.3f", cor(res[, 1], abs(res[, 2])))

[1] "  correlation between the p-value and |bias|: 0.024"

Read those three lines slowly.

The test rejects the wrong model in 17 per cent of surveys, so 83 per cent of the time it hands the wrong model a clean bill of health. That is a power problem, and power problems are familiar and forgivable.

The second line is the one that should worry you. Among the surveys that pass, the mean density bias is +12.2 per cent, which is the same as the mean bias across all surveys. Passing the test tells you nothing about whether your density is right. The worst passing survey is out by 38 per cent, and it passed.

And the third line says it exactly: the correlation between the p-value and the size of the error is 0.024. Not weak. Zero, to within the noise of 400 surveys. The statistic the test computes and the quantity you care about are unrelated.

A goodness-of-fit test asks whether the fitted curve describes the shape of the histogram. Your density estimate depends on the height of that curve at zero. Those are different questions, and a model can answer the first one well while getting the second one wrong, because five bins of a smooth curve leave plenty of room to be wrong in the corner where nobody is looking.

The assumption with no fingerprints

That was the diagnostic being weak. This one is different in kind.

Every density estimate in this series divides by an effective width computed on the assumption that $g(0) = 1$: everything on the transect line is detected. It is the first assumption in every textbook, and it is the one that fails most often in practice. Divers miss animals below them. Aerial observers miss animals under the plane. A bird singing directly overhead in the canopy is not obviously easier to detect than one at 20 m.

Suppose it is false, and the true detection function is $c \cdot g(x)$ for some $c < 1$. Fit your model to the distances and see what happens.

set.seed(4283)
s_t <- 0.05; L <- 30; Dtrue <- 60
N <- rpois(1, Dtrue*2*w*L); xa <- runif(N, 0, w)
x <- xa[runif(N) < g_hn(xa, s_t)]
base <- optimise(function(p) -sum(ld_hn(x, exp(p), w)), c(-8, 1), tol = 1e-12)
out <- t(sapply(c(1, 0.8, 0.6, 0.4), function(g0) {
  ## the world where only a fraction g0 of the animals on the line are seen
  nll <- function(p) { s <- exp(p)
    -sum(log(g0*g_hn(x, s)) - log(g0*mu_hn(s, w))) }
  o <- optimise(nll, c(-8, 1), tol = 1e-12)
  c(g0 = g0, sigma = exp(o$minimum), d_sigma = exp(o$minimum) - exp(base$minimum),
    loglik = -o$objective, d_loglik = -o$objective + base$objective,
    D = length(x)/(2*L*mu_hn(exp(o$minimum), w)*g0))
}))
sprintf("one survey, n = %d detections", length(x))

[1] "one survey, n = 223 detections"

print(round(out, 8))

      g0      sigma d_sigma   loglik d_loglik         D
[1,] 1.0 0.05705162       0 531.7441        0  56.47631
[2,] 0.8 0.05705162       0 531.7441        0  70.59539
[3,] 0.6 0.05705162       0 531.7441        0  94.12719
[4,] 0.4 0.05705162       0 531.7441        0 141.19078

Look at the middle columns. The fitted scale is the same in all four rows, and so is the log-likelihood: those differences are not small, they are zero to the last bit an optimiser can carry. Meanwhile the density runs from 56.5 to 141.2, a factor of 2.5.

This is not a weak diagnostic. It is an exact invariance, and it takes one line of algebra to see why. The likelihood is built from the probability density of the observed distances,

\[f(x) = \frac{c\,g(x)}{\int_0^w c\,g(t)\,\mathrm{d}t} = \frac{g(x)}{\mu}\]

and the constant cancels. Your distances are conditional on detection, so they know the shape of the detection function and nothing whatsoever about its height. Every estimate in the series scales as $\hat{D} = D/c$: +25 per cent at $c = 0.8$, +67 per cent at $c = 0.6$, +150 per cent at $c = 0.4$.

No goodness-of-fit test can see this, because the fit is identical. No residual can see it. No AIC comparison can see it, because both models have the same likelihood and the same parameter count. There is no diagnostic, no plot and no statistic in R that will tell you, because the information is not in the data. The assumption you most need is the one your survey cannot check.

Left panel shows four detection curves at different heights, all producing one identical density of observed distances. Right panel shows the fitted sigma and log-likelihood flat across g(0) while density rises steeply. — Figure 2: Left: four detection functions differing only in their height at zero, and the single distribution of observed distances that all four produce. Right: what changes and what does not as g(0) falls.

What to check, and what to decide instead

Sorting the assumptions by whether the data can speak to them produces a short and lopsided table.

Assumption	Can the distances check it	What settles it
Key function has the right shape	Weakly	Goodness of fit, AIC, shouldered key by default
Distances measured accurately	Partly	Repeat measurement of a subset
No heaping on round numbers	Yes	The histogram, and a rangefinder
Detection depends on the recorded covariates	Partly	AIC, once the covariate is recorded
Animals uniform with respect to the line	No	Random placement of the design
g(0) = 1	No	Double observer or mark-recapture distance sampling
Animals do not move before detection	No	Snapshot protocol, or a faster observer

Two of the seven rows are about analysis. The other five are decisions made before anyone leaves the office, and the pattern is the one this blog keeps running into: the strongest thing you can do for an estimate is to make an assumption true by construction rather than to test it afterwards.

For $g(0)$ specifically, the fix is a second look. Two observers working independently, or a double platform, gives you a mark-recapture sample within the distance sampling survey: animals seen by one observer and missed by the other carry the information about what both of them missed. That information is not in a single observer’s distances at any sample size, and it never will be. Buckland et al. (2015) work through the designs; the point here is only that it is a field decision, not an analysis one.

The honest close

Distance sampling earns its reputation. It needs one visit rather than three. It gives an absolute density rather than an index. Its central property is real: pooling unmodelled heterogeneity leaves the study-area total close to unbiased (Rexstad et al. 2023), which is a genuinely unusual thing for an estimator to do.

The failures in this series are consistent with all of that, and they share a shape. The estimator reads one number off a curve at the origin, and everything that goes wrong goes wrong by moving that number: pooling moves it the same way for every stratum and destroys the comparison while sparing the total; the wrong key moves it more at a point than on a line; a missing $g(0)$ scales it by a constant nobody can estimate from distances. In every case the fit stays clean, the standard error stays narrow, and the output looks exactly like an output that is right.

Which is the sentence this cluster has been building towards, and it is not really about distance sampling: the diagnostics available to you are the ones that operate on the data you collected, and the assumptions that hurt you most are the ones that determined what data you were able to collect. Check what you can. Decide the rest with a protocol, write down which is which, and report the number with the assumption attached rather than with an interval that has quietly forgotten it.

References

Buckland, S. T. (2006). Point transect surveys for songbirds: robust methodologies. The Auk 123(2), 345-357. DOI: 10.1093/auk/123.2.345

Buckland, S. T., Rexstad, E. A., Marques, T. A. and Oedekoven, C. S. (2015). Distance Sampling: Methods and Applications. Springer. ISBN 978-3-319-19218-5

Marques, F. F. C. and Buckland, S. T. (2003). Incorporating covariates into standard line transect analyses. Biometrics 59(4), 924-935. DOI: 10.1111/j.0006-341X.2003.00107.x

Rexstad, E., Buckland, S. T., Marshall, L. and Borchers, D. (2023). Pooling robustness in distance sampling: avoiding bias when there is unmodelled heterogeneity. Ecology and Evolution 13(1), e9684. DOI: 10.1002/ece3.9684

--- title: "Checking a distance sampling model" description: "The goodness-of-fit test accepts surveys whose density is 13 per cent wrong, and the assumption that matters most leaves no trace in the likelihood at all." date: 2026-08-26 12:00 categories: [R, distance sampling, model checking, abundance, ecology tutorial] image: thumbnail.png image-alt: "A scatter of goodness-of-fit p-values against density bias showing no relationship" --- Three posts in this series have shown a distance sampling analysis returning a confident wrong number: [a pooled detection function](../covariates-in-the-detection-function/) reversing which habitat is denser, [a half-normal at a point](../point-transect-distance-sampling/) inflating density by a sixth, [a hierarchical model](../hierarchical-distance-sampling/) halving a population when the animals hold still. The reasonable response is to ask what you would check. Distance sampling has a natural answer, because unlike most ecological models it makes a direct, testable prediction about data you can see: the histogram of distances should look like the fitted curve. So plot it, run a goodness-of-fit test, and if it passes, you are fine. This post takes that programme seriously and measures how much protection it buys. The honest answer is: some, and a good deal less than the ritual implies. ```{r setup} #| message: false #| warning: false library(ggplot2) g_hn <- function(x, s) exp(-x^2/(2*s^2)) mu_hn <- function(s, w) sqrt(2*pi)*s*(pnorm(w, 0, s) - 0.5) ld_hn <- function(x, s, w) -x^2/(2*s^2) - log(mu_hn(s, w)) g_hr <- function(x, s, b) 1 - exp(-(x/s)^(-b)) mu_hr <- function(s, b, w) integrate(g_hr, 0, w, s = s, b = b)$value fit_hn <- function(x, w) exp(optimise(function(p) -sum(ld_hn(x, exp(p), w)), c(-8, 1), tol = 1e-10)$minimum) w <- 0.1 ``` ## The check that works Start with the good news, because there is some. Bin the distances, fit the key function, compare observed counts against expected counts with a chi-square statistic. Before trusting a test to reject a bad model, check that it accepts a good one at the advertised rate: a test with the wrong size is not evidence, it is a coin with an opinion. ```{r calibrate} cut <- seq(0, w, length.out = 6) pi_hn <- function(s) { P <- sqrt(2*pi)*s*(pnorm(cut[-1], 0, s) - pnorm(cut[-6], 0, s))/w P/sum(P) } gof_hn <- function(x) { s <- fit_hn(x, w) ob <- tabulate(findInterval(x, cut, rightmost.closed = TRUE), 5) ex <- length(x)*pi_hn(s) c(p = pchisq(sum((ob - ex)^2/ex), df = 5 - 1 - 1, lower.tail = FALSE), sigma = s) } ## a world where the half-normal is exactly right s_ok <- 0.045; D0 <- 100; L_ok <- 266/(2*D0*mu_hn(s_ok, w)) cal <- t(sapply(1:400, function(i) { set.seed(61000 + i) N <- rpois(1, D0*2*w*L_ok); x <- runif(N, 0, w) gof_hn(x[runif(N) < g_hn(x, s_ok)]) })) sprintf("model-true data, n approx 266: test rejects at 0.05 in %.1f%% of surveys | median p %.3f", 100*mean(cal[, "p"] < 0.05), median(cal[, "p"])) ``` Well behaved: it rejects the truth at about the nominal rate, so it is a real test rather than a formality. And it does catch gross problems. A spike of animals recorded at zero distance, a histogram with a hole in it, distances heaped on round numbers: all of that shows up. ## What it misses Now give it the job it is actually asked to do. The true detection function has a shoulder, generated by a hazard-rate key. You fit a half-normal, which cannot make a shoulder. From the [point transect post](../point-transect-distance-sampling/), we know this costs about 11 per cent on a line. Does the test tell you? ```{r blind} s_hr <- 0.05; b_hr <- 3 mu_t <- mu_hr(s_hr, b_hr, w); Ln <- 266/(2*D0*mu_t) res <- t(sapply(1:400, function(i) { set.seed(63000 + i) N <- rpois(1, D0*2*w*Ln); x <- runif(N, 0, w) x <- x[runif(N) < g_hr(x, s_hr, b_hr)] g <- gof_hn(x) c(p = g["p"], bias = length(x)/(2*Ln*mu_hn(g["sigma"], w))/D0 - 1, n = length(x)) })) acc <- res[res[, 1] > 0.05, ] sprintf("half-normal fitted to hazard-rate data (n approx %.0f, 5 bins, df 3):", mean(res[, 3])) sprintf(" rejects at 0.05: %.1f%% of surveys | median p-value %.3f | mean density bias %+.1f%%", 100*mean(res[, 1] < 0.05), median(res[, 1]), 100*mean(res[, 2])) sprintf(" among the surveys the test ACCEPTS: mean bias %+.1f%%, worst %+.1f%%", 100*mean(acc[, 2]), 100*max(acc[, 2])) sprintf(" correlation between the p-value and |bias|: %.3f", cor(res[, 1], abs(res[, 2]))) ``` Read those three lines slowly. The test rejects the wrong model in `r sprintf("%.0f", 100*mean(res[, 1] < 0.05))` per cent of surveys, so `r sprintf("%.0f", 100*mean(res[, 1] > 0.05))` per cent of the time it hands the wrong model a clean bill of health. That is a power problem, and power problems are familiar and forgivable. The second line is the one that should worry you. Among the surveys that *pass*, the mean density bias is `r sprintf("%+.1f", 100*mean(acc[, 2]))` per cent, which is the same as the mean bias across all surveys. Passing the test tells you nothing about whether your density is right. The worst passing survey is out by `r sprintf("%.0f", 100*max(acc[, 2]))` per cent, and it passed. And the third line says it exactly: the correlation between the p-value and the size of the error is `r sprintf("%.3f", cor(res[, 1], abs(res[, 2])))`. Not weak. Zero, to within the noise of 400 surveys. The statistic the test computes and the quantity you care about are unrelated. > A goodness-of-fit test asks whether the fitted curve describes the shape of the histogram. Your density estimate depends on the height of that curve at zero. Those are different questions, and a model can answer the first one well while getting the second one wrong, because five bins of a smooth curve leave plenty of room to be wrong in the corner where nobody is looking. ```{r fig-blind} #| echo: false #| fig-width: 8 #| fig-height: 3.6 #| fig-cap: "Four hundred surveys of a shouldered detection function analysed with a half-normal. Left: the p-value against the resulting density error, with the accepted surveys in green. Right: the distribution of density error among surveys the test accepts, against the truth." #| fig-alt: "Left panel is a scatter of goodness-of-fit p-values against density bias showing a formless cloud with no trend. Right panel shows a histogram of density bias among accepted surveys, sitting entirely above the true value." ink <- "#16241d"; forest <- "#275139"; gold <- "#c9b458"; rust <- "#b5534e" dd <- data.frame(p = res[, 1], bias = 100*res[, 2], k = ifelse(res[, 1] > 0.05, "test accepts", "test rejects")) p1 <- ggplot(dd, aes(p, bias, colour = k)) + geom_hline(yintercept = 0, colour = ink, linetype = 2) + geom_vline(xintercept = 0.05, colour = rust, linetype = 3) + geom_point(size = 1, alpha = 0.7) + scale_colour_manual(values = c("test accepts" = forest, "test rejects" = rust)) + labs(x = "goodness-of-fit p-value", y = "density error (%)", colour = NULL, title = "No relationship at all", subtitle = sprintf("correlation of p with |error|: %.3f", cor(res[, 1], abs(res[, 2])))) + theme_minimal(base_size = 10) + theme(legend.position = "bottom", plot.title = element_text(colour = ink, face = "bold"), plot.subtitle = element_text(size = 8)) p2 <- ggplot(dd[dd$k == "test accepts", ], aes(bias)) + geom_histogram(bins = 30, fill = forest, alpha = 0.6) + geom_vline(xintercept = 0, colour = ink, linetype = 2, linewidth = 0.7) + labs(x = "density error (%)", y = "surveys", title = "The surveys that passed", subtitle = sprintf("mean error %+.1f%%, none of them near zero", 100*mean(acc[, 2]))) + theme_minimal(base_size = 10) + theme(plot.title = element_text(colour = ink, face = "bold"), plot.subtitle = element_text(size = 8)) grid::grid.newpage() grid::pushViewport(grid::viewport(layout = grid::grid.layout(1, 2))) print(p1, vp = grid::viewport(layout.pos.row = 1, layout.pos.col = 1)) print(p2, vp = grid::viewport(layout.pos.row = 1, layout.pos.col = 2)) ``` ## The assumption with no fingerprints That was the diagnostic being weak. This one is different in kind. Every density estimate in this series divides by an effective width computed on the assumption that $g(0) = 1$: everything on the transect line is detected. It is the first assumption in every textbook, and it is the one that fails most often in practice. Divers miss animals below them. Aerial observers miss animals under the plane. A bird singing directly overhead in the canopy is not obviously easier to detect than one at 20 m. Suppose it is false, and the true detection function is $c \cdot g(x)$ for some $c < 1$. Fit your model to the distances and see what happens. ```{r g0} set.seed(4283) s_t <- 0.05; L <- 30; Dtrue <- 60 N <- rpois(1, Dtrue*2*w*L); xa <- runif(N, 0, w) x <- xa[runif(N) < g_hn(xa, s_t)] base <- optimise(function(p) -sum(ld_hn(x, exp(p), w)), c(-8, 1), tol = 1e-12) out <- t(sapply(c(1, 0.8, 0.6, 0.4), function(g0) { ## the world where only a fraction g0 of the animals on the line are seen nll <- function(p) { s <- exp(p) -sum(log(g0*g_hn(x, s)) - log(g0*mu_hn(s, w))) } o <- optimise(nll, c(-8, 1), tol = 1e-12) c(g0 = g0, sigma = exp(o$minimum), d_sigma = exp(o$minimum) - exp(base$minimum), loglik = -o$objective, d_loglik = -o$objective + base$objective, D = length(x)/(2*L*mu_hn(exp(o$minimum), w)*g0)) })) sprintf("one survey, n = %d detections", length(x)) print(round(out, 8)) ``` Look at the middle columns. The fitted scale is the same in all four rows, and so is the log-likelihood: those differences are not small, they are zero to the last bit an optimiser can carry. Meanwhile the density runs from `r sprintf("%.1f", out[1, "D"])` to `r sprintf("%.1f", out[4, "D"])`, a factor of `r sprintf("%.1f", out[4, "D"]/out[1, "D"])`. This is not a weak diagnostic. It is an exact invariance, and it takes one line of algebra to see why. The likelihood is built from the probability density of the *observed* distances, $$f(x) = \frac{c\,g(x)}{\int_0^w c\,g(t)\,\mathrm{d}t} = \frac{g(x)}{\mu}$$ and the constant cancels. Your distances are conditional on detection, so they know the *shape* of the detection function and nothing whatsoever about its height. Every estimate in the series scales as $\hat{D} = D/c$: `r sprintf("%+.0f", 100*(1/0.8 - 1))` per cent at $c = 0.8$, `r sprintf("%+.0f", 100*(1/0.6 - 1))` per cent at $c = 0.6$, `r sprintf("%+.0f", 100*(1/0.4 - 1))` per cent at $c = 0.4$. No goodness-of-fit test can see this, because the fit is identical. No residual can see it. No AIC comparison can see it, because both models have the same likelihood and the same parameter count. There is no diagnostic, no plot and no statistic in R that will tell you, because the information is not in the data. The assumption you most need is the one your survey cannot check. ```{r fig-g0} #| echo: false #| fig-width: 8 #| fig-height: 3.6 #| fig-cap: "Left: four detection functions differing only in their height at zero, and the single distribution of observed distances that all four produce. Right: what changes and what does not as g(0) falls." #| fig-alt: "Left panel shows four detection curves at different heights, all producing one identical density of observed distances. Right panel shows the fitted sigma and log-likelihood flat across g(0) while density rises steeply." xs <- seq(0, w, length.out = 300) cur <- do.call(rbind, lapply(c(1, 0.8, 0.6, 0.4), function(g0) data.frame(x = xs, g = g0*g_hn(xs, s_t), k = sprintf("g(0) = %.1f", g0)))) p1 <- ggplot(cur, aes(x*1000, g, colour = k)) + geom_line(linewidth = 0.9) + geom_line(data = data.frame(x = xs, g = g_hn(xs, s_t)/mu_hn(s_t, w)/20.2, k = "observed f(x)"), colour = ink, linetype = 2, linewidth = 1) + scale_colour_manual(values = c("g(0) = 1.0" = forest, "g(0) = 0.8" = gold, "g(0) = 0.6" = rust, "g(0) = 0.4" = "#93a87f")) + labs(x = "perpendicular distance (m)", y = "detection probability", colour = NULL, title = "Four worlds, one dataset", subtitle = "dashed: the distances you record, scaled; identical in all four") + theme_minimal(base_size = 10) + theme(legend.position = "bottom", plot.title = element_text(colour = ink, face = "bold"), plot.subtitle = element_text(size = 8), legend.text = element_text(size = 7)) sc <- data.frame(g0 = rep(out[, "g0"], 3), v = c(out[, "sigma"]/out[1, "sigma"], out[, "loglik"]/out[1, "loglik"], out[, "D"]/out[1, "D"]), k = rep(c("fitted sigma", "log-likelihood", "estimated density"), each = 4)) p2 <- ggplot(sc, aes(g0, v, colour = k)) + geom_hline(yintercept = 1, colour = ink, linetype = 3) + geom_line(linewidth = 0.9) + geom_point(size = 2) + scale_colour_manual(values = c("fitted sigma" = forest, "log-likelihood" = gold, "estimated density" = rust)) + scale_x_reverse() + labs(x = "true g(0)", y = "relative to the g(0) = 1 analysis", colour = NULL, title = "Everything you can check is flat", subtitle = "the one thing you report is not") + theme_minimal(base_size = 10) + theme(legend.position = "bottom", plot.title = element_text(colour = ink, face = "bold"), plot.subtitle = element_text(size = 8)) grid::grid.newpage() grid::pushViewport(grid::viewport(layout = grid::grid.layout(1, 2))) print(p1, vp = grid::viewport(layout.pos.row = 1, layout.pos.col = 1)) print(p2, vp = grid::viewport(layout.pos.row = 1, layout.pos.col = 2)) ``` ## What to check, and what to decide instead Sorting the assumptions by whether the data can speak to them produces a short and lopsided table. ```{r checklist} #| echo: false knitr::kable(data.frame( Assumption = c("Key function has the right shape", "Distances measured accurately", "No heaping on round numbers", "Detection depends on the recorded covariates", "Animals uniform with respect to the line", "g(0) = 1", "Animals do not move before detection"), `Can the distances check it` = c("Weakly", "Partly", "Yes", "Partly", "No", "No", "No"), `What settles it` = c("Goodness of fit, AIC, shouldered key by default", "Repeat measurement of a subset", "The histogram, and a rangefinder", "AIC, once the covariate is recorded", "Random placement of the design", "Double observer or mark-recapture distance sampling", "Snapshot protocol, or a faster observer"), check.names = FALSE)) ``` Two of the seven rows are about analysis. The other five are decisions made before anyone leaves the office, and the pattern is the one this blog keeps running into: the strongest thing you can do for an estimate is to make an assumption true by construction rather than to test it afterwards. For $g(0)$ specifically, the fix is a second look. Two observers working independently, or a double platform, gives you a mark-recapture sample within the distance sampling survey: animals seen by one observer and missed by the other carry the information about what both of them missed. That information is not in a single observer's distances at any sample size, and it never will be. Buckland et al. (2015) work through the designs; the point here is only that it is a field decision, not an analysis one. ## The honest close Distance sampling earns its reputation. It needs one visit rather than three. It gives an absolute density rather than an index. Its central property is real: pooling unmodelled heterogeneity leaves the study-area total close to unbiased (Rexstad et al. 2023), which is a genuinely unusual thing for an estimator to do. The failures in this series are consistent with all of that, and they share a shape. The estimator reads one number off a curve at the origin, and everything that goes wrong goes wrong by moving that number: pooling moves it the same way for every stratum and destroys the comparison while sparing the total; the wrong key moves it more at a point than on a line; a missing $g(0)$ scales it by a constant nobody can estimate from distances. In every case the fit stays clean, the standard error stays narrow, and the output looks exactly like an output that is right. Which is the sentence this cluster has been building towards, and it is not really about distance sampling: the diagnostics available to you are the ones that operate on the data you collected, and the assumptions that hurt you most are the ones that determined what data you were able to collect. Check what you can. Decide the rest with a protocol, write down which is which, and report the number with the assumption attached rather than with an interval that has quietly forgotten it. ## Related tutorials - [Distance sampling on line transects](../distance-sampling-line-transect/) - [Covariates in the detection function](../covariates-in-the-detection-function/) - [Point transect distance sampling](../point-transect-distance-sampling/) - [Hierarchical distance sampling](../hierarchical-distance-sampling/) ## References Buckland, S. T. (2006). Point transect surveys for songbirds: robust methodologies. *The Auk* 123(2), 345-357. DOI: 10.1093/auk/123.2.345 Buckland, S. T., Rexstad, E. A., Marques, T. A. and Oedekoven, C. S. (2015). *Distance Sampling: Methods and Applications*. Springer. ISBN 978-3-319-19218-5 Marques, F. F. C. and Buckland, S. T. (2003). Incorporating covariates into standard line transect analyses. *Biometrics* 59(4), 924-935. DOI: 10.1111/j.0006-341X.2003.00107.x Rexstad, E., Buckland, S. T., Marshall, L. and Borchers, D. (2023). Pooling robustness in distance sampling: avoiding bias when there is unmodelled heterogeneity. *Ecology and Evolution* 13(1), e9684. DOI: 10.1002/ece3.9684