Internals: plotting

This document illustrates how plotting is performed in mrds. Simply plotting a detection function is relatively straight forward, however when adding covariates and double observer components, this can get complicated.

It’s assumed you have some familiarity with mrds and a lot of familiarity with R. This is not supposed to be a replacement for the mrds documentation, but rather further explanation for those hacking on mrds.

The main point covered here at the moment is how the detection function is evaluated to plot the lines, not any of the other content in the plot.

Terms of reference

  • $y$ refers to distance: this could be perpendicular or radial.
  • $g(y, \mathbf{z})$ a general detection function (missing the $\mathbf{z}$ if there are no covariates), we assume that the detection function has some parameters ($\mathbf{\theta}$, say) but they aren’t important here.
  • $p_i = \int_0^w g(y,\mathbf{z}_i) \pi(y) \text{d}y$ (i.e. the fitted values of the detection function), where for line transects $\pi(y) = \frac{1}{w}$ and for point transects $\pi(y)=\frac{2y}{w^2}$.

Single observer – plot.ds

detection function only (dsmodel), no covariates

When there are no covariates, the detection function is simply evaluated over the grid (between left and right truncations).

detection function only (dsmodel), with covariates

When covariates are present, the detection function is evaluated at a series of distances between left and right truncations (say seq(left,width,100)). Then for each distance ($y_t$), the detection function is evaluated for each observed covariate combination ($\mathbf{z}_i$). Then the detection function values are averaged at each distance, weighted by the fitted probabilities of detection ($p_i$). For a given $y_t$, we evaluate:

\[\frac{\left(\sum_{i=1}^n\frac{ g(y_t, \mathbf{z}_i,\mathbf{\theta})}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\]

Double observer plots

All of plot.io, plot.io.fi, plot.rem, plot.rem.fi, plot.trial and plot.trial.fi call plot_cond when making plots of the conditional detection function, which in turn calls average.line.cond. They also each call plot_uncond for plots of the unconditional detection function, this in turn calls average.line.

Note that for gamma detection function models (and presumably other models where $g(y,\mathbf{z})=1$ for $x \neq 0$), one must set the distance in the predictions in average.line to be the apex rather than zero (since this is where detection is assumed to be certain).

Buckland et al (2004, Chapter 6) is the key reference. Pages 130-140 are most useful.

average.line.cond

Calculates the conditional detection function lines.

Calls predict per distance (over a grid of distances, as above), then takes the mean per distance. This means that for each distance, the detection function is evaluated at all observed covariate combinations.

For io, trial and rem methods, the predictions use the $mr part of the model. It produces a detection function line for each observer for io, io.fi, rem and rem.fi methods, but only for the first observer for the trial and trial.fi methods.

average.line

Calculates the unconditional detection function lines.

As with average.line.cond, predictions are made over a grid of distances for each covariate combinations. However, for the unconditional plots, the averages are weighted.

average.line – Point independence models

For point independence models (io, rem and trial), need to use the delta function given in Buckland et al (2004) p. 131, since we assume independence at a point (y=0).

The g0 part of the delta depends on the model. (Note that g0 in the code is $\hat{p}^c_.(0,\mathbf{z})$ in equation (6.28), Buckland et al (2004), Chapter 6.)

Method g0 prediction model
io model$mr
rem model$mr
trial model$mr$mr

So in the code, the lines:

detfct.pooled.values <- detfct(newdat$distance[newdat$observer==1],
                                 ddfobj,width=model$meta.data$width-
                                              model$meta.data$left)
deltax <- detfct.pooled.values/(cond.det$fitted/g0)

translate to the equations (6.27) and substituting in (6.28) in ADS:

\(\delta(y,\mathbf{z}) = \frac{\hat{p}^c_.(y,\mathbf{z})}{\hat{p}_.(y,\mathbf{z})} \qquad \text{(6.27)}\) \(\hat{p}_.(y,\mathbf{z}) = \hat{p}^c_.(0,\mathbf{z})\hat{g}_.(y,\mathbf{z}) \qquad \text{(6.28)}\) \(\delta(y,\mathbf{z})= \frac{\hat{g}_.(y,\mathbf{z}) \hat{p}^c_.(0,\mathbf{z})}{\hat{p}^c_.(y,\mathbf{z})}\)

So deltax in the code is $1/\delta(y,\mathbf{z})$ in maths. We need $\delta$ since we can only estimate $\hat{p}^c_.$ (Buckland et al 2004, p 131).

Depending on observer, computes average values detection function values are calculated for the given distance using Buckland et al (2004) equation (6.50). Observer options are observer 1 or 2 (obs==1 or obs==2), duplicates (obs==3) or “pooled observation team” (obs==4).

  • obs==1: sum(p1*deltax/prob.det)/sum(1/prob.det) \(\frac{\left(\sum_{i=1}^n\frac{ p_{1|2}(y) / \delta(x_i,\mathbf{z})}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)
  • obs==2: sum(p2*deltax/prob.det)/sum(1/prob.det) (not produced for trial) \(\frac{\left(\sum_{i=1}^n\frac{ p_{2|1}(y) / \delta(x_i,\mathbf{z})}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)
  • obs==3 (duplicates): sum(g0*detfct.pooled.values/prob.det)/sum(1/prob.det) \(\frac{\left(\sum_{i=1}^n\frac{ p_{1|2}(y) + p_{2|1}(y) - p_{1|2}(y) p_{2|1}(y)}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)
  • otherwise (i.e. obs==4, “pooled observer team”): sum(p1*p2*deltax/prob.det)/sum(1/prob.det) \(\frac{\left(\sum_{i=1}^n\frac{ (p_{1|2}(y) p_{2|1}(y))/ \delta(x_i,\mathbf{z})}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)

average.line – Full independence models

Things are simpler for the full independence model, since there is no $\delta$ to consider.

Again, depending on observer, computes (simple) average values based on observed covariates for a given distance. Observer options are observer 1 or 2 (obs==1 or obs==2), duplicates (obs==3) or “pooled observation team” (obs==4).

  • obs==1: sum(p1/prob.det)/sum(1/prob.det) \(\frac{\left(\sum_{i=1}^n\frac{p_{1|2}(y)}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)
  • obs==2: sum(p2/prob.det)/sum(1/prob.det) (not produced for trial.fi) \(\frac{\left(\sum_{i=1}^n\frac{p_{2|1}(y)}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)
  • obs==3 (duplicates): sum(g0*detfct.pooled.values/prob.det)/sum(1/prob.det) \(\frac{\left(\sum_{i=1}^n\frac{ p_{1|2}(y) + p_{2|1}(y) - p_{1|2}(y) p_{2|1}(y)}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)
  • otherwise (i.e. obs==4, “pooled observer team”): sum(p1*p2/prob.det)/sum(1/prob.det) \(\frac{\left(\sum_{i=1}^n\frac{p_{1|2}(y) p_{2|1}(y)}{p_i} \right) }{ \left( \sum_{i=1}^n\frac{1}{p_i} \right)}\)

References

Buckland, ST, DR Anderson, KP Burnham, JL Laake, DL Borchers, and L Thomas. Advanced Distance Sampling, Oxford University Press, 2004.