<p>Great you are using our methods. I have extracted the code I used for making the p0 predictions in Tables 2 and 4 for our paper. I hope this helps. Please let me know how you get on with this.</p>
<p>Cheers, Cornelia</p>
<pre><code>### MRDS visual vs acoustic
# code written by Cornelia Oedekoven, including modified functions from mrds package

##########################################################
library(Distance)
library(mrds)

##########################################################
# Predicting p(0) in the trial configuration
# p(0) estimates per beaufort state
newdata.trial&lt;-data.frame(distance=rep(0,6),beaufort=c(0:5))
# predict will give you the average detection probability p
# not p(0)
# doesn't matter which distance you use
g0.trial&lt;-round(predict(object = best.trial,newdata = newdata.trial,integrate = F)$fitted, 3)
g0.trial

# average p from relative det fct g(y)
summary(best.trial)$ds.summary$average.p

# p(0) for the different beaufort states 
# =&gt; to obtain p(0) for the different beaufort states we take the p predicitons per beaufort 
# and divide them by the p from the relative detection function
p0.beaufort&lt;-g0.trial/summary(best.trial)$ds.summary$average.p
# Table 2 (Rankin et al. 2020)
p0.beaufort

###################################################################
#### Predicting p0 in the independent observer configuration
#' @param object ddf.io object
#' @param newdata New data for making predictions
my.predict.io &lt;-
function(object,newdata=NULL,compute=FALSE,int.range=NULL,...){

  model &lt;- object
  if(is.null(newdata)){
    xmat &lt;- model$mr$mr$data
  }else{
    compute &lt;- TRUE
    xmat &lt;- newdata
  }
  
  xmat$distance &lt;- 0
  ddfobj &lt;- model$ds$ds$aux$ddfobj
  
  # for gamma models need to find where p(x)=1 (apex), set that as distance
  if(ddfobj$type==&quot;gamma&quot;){
    xmat$distance &lt;- rep(apex.gamma(ddfobj),2)
  }
  
  # calculate ps for each part of the model
  xmat$offsetvalue &lt;- 0
  p.0 &lt;- predict(model$mr,newdata=xmat,integrate=FALSE,compute=compute)
  if(is.null(newdata)){
    pdot &lt;- predict(model$ds,esw=FALSE,compute=compute,
                    int.range=int.range)$fitted
  }else{
    pdot &lt;- predict(model$ds,newdata=newdata[newdata$observer==1,],
                    esw=FALSE,compute=compute,int.range=int.range)$fitted
  }
  
  return(list(p.0 = p.0, pdot = pdot))
}

# example with Beaufort (Table 4 from Rankin et al. 2020) 
bft=c(0:5)
newdata.io&lt;-data.frame(distance=rep(0,12),beaufort=rep(bft,each=2),observer=rep(c(1,2),times=6))
io.pred&lt;-data.frame(matrix(NA,3,(length(bft)+1)))
colnames(io.pred)&lt;-c(&quot;obs&quot;,bft)
io.pred[,1]&lt;-c(&quot;obs1&quot;,&quot;obs2&quot;,&quot;both&quot;)
io.pred[1,2:(length(bft)+1)]&lt;-my.predict.io(best.io,newdata.io)$p.0$p1
io.pred[2,2:(length(bft)+1)]&lt;-my.predict.io(best.io,newdata.io)$p.0$p2
io.pred[3,2:(length(bft)+1)]&lt;-my.predict.io(best.io,newdata.io)$p.0$fitted
print(io.pred, row.names = F)
</code></pre>

CC BY-SA 4.0
2022-10-06T07:43:33.963Z
As it’s currently written, your answer is unclear. Please [edit] to add additional details that will help others understand how this addresses the question asked. You can find more information on how to write good answers [in the help center](/help/how-to-answer).
Community
CC BY-SA 4.0
2022-10-05T14:57:45.887Z
2022-10-05T14:57:45.887Z
Cornelia Oedekoven
<p>Shipboard line-transect distance sampling methods are commonly used to estimate populations of cetaceans. One of the assumptions of these methods is that all animals on the trackline are detected, but this is not always the case (especially for cryptic species).  Mark Recapture Distance Sampling (MRDS) is a method that uses a secondary observation team to identify the fraction of animals detected on the trackline when it is suspected that animals have been missed. During typical visual line transect methods, MRDS use a secondary visual observation team to understand the fraction of animals on the trackline that are missed by the primary team. In this paper (<a href="https://link.springer.com/article/10.1007/s10651-020-00443-7" rel="nofollow noreferrer">https://link.springer.com/article/10.1007/s10651-020-00443-7</a>), acoustic detections from a towed hydrophone array are used as the secondary 'observation' platform.</p>
<p>One component of the methods presented in this paper estimates the trackline detection probability p(0) for different levels of an additional variable such as Beaufort Sea States (as presented in Table 2). This was not part of the MRDS output in R. How can these methods be replicated? Specifically, if I have completed MRDS analysis for my dataset, how can I estimate trackline detection probability for individual Beaufort sea states (or other variables with multiple levels)?</p>

CC BY-SA 4.0
2022-10-04T01:57:16.640Z
Welcome, Melissa! Could you add some information about what package you are using in R?
selene
CC BY-SA 4.0
2022-10-05T07:07:00.507Z
Thank you Selene. I am using the mrds package in R (Laake et al. 2018).
Melissa Costagliola-Ray
<analysis><density-estimation>
Mark Recapture Distance Sampling, estimating trackline detection probability p(0) for additional variables such as Beaufort Sea States
# Mark Recapture Distance Sampling, estimating trackline detection probability p(0) for additional variables such as Beaufort Sea States

(Asked by: username_0 on 2022-10-03)

Shipboard line-transect distance sampling methods are commonly used to estimate populations of cetaceans. One of the assumptions of these methods is that all animals on the trackline are detected, but this is not always the case (especially for cryptic species). Mark Recapture Distance Sampling (MRDS) is a method that uses a secondary observation team to identify the fraction of animals detected on the trackline when it is suspected that animals have been missed. During typical visual line transect methods, MRDS use a secondary visual observation team to understand the fraction of animals on the trackline that are missed by the primary team. In this paper (<https://link.springer.com/article/10.1007/s10651-020-00443-7>), acoustic detections from a towed hydrophone array are used as the secondary 'observation' platform.

One component of the methods presented in this paper estimates the trackline detection probability p(0) for different levels of an additional variable such as Beaufort Sea States (as presented in Table 2). This was not part of the MRDS output in R. How can these methods be replicated? Specifically, if I have completed MRDS analysis for my dataset, how can I estimate trackline detection probability for individual Beaufort sea states (or other variables with multiple levels)?

**Tags:** analysis, density-estimation

> **Comments:**
> **username_1**: Welcome, username_0! Could you add some information about what package you are using in R?  
> **username_0**: Thank you username_1. I am using the mrds package in R (Laake et al. 2018).  

---

## 1 answers

### Answer by username_2 on 2022-10-05. Score: 3

Great you are using our methods. I have extracted the code I used for making the p0 predictions in Tables 2 and 4 for our paper. I hope this helps. Please let me know how you get on with this.

Cheers, username_2

```
### MRDS visual vs acoustic
# code written by username_2 username_2, including modified functions from mrds package

##########################################################
library(Distance)
library(mrds)

##########################################################
# Predicting p(0) in the trial configuration
# p(0) estimates per beaufort state
newdata.trial<-data.frame(distance=rep(0,6),beaufort=c(0:5))
# predict will give you the average detection probability p
# not p(0)
# doesn't matter which distance you use
g0.trial<-round(predict(object = best.trial,newdata = newdata.trial,integrate = F)$fitted, 3)
g0.trial

# average p from relative det fct g(y)
summary(best.trial)$ds.summary$average.p

# p(0) for the different beaufort states 
# => to obtain p(0) for the different beaufort states we take the p predicitons per beaufort 
# and divide them by the p from the relative detection function
p0.beaufort<-g0.trial/summary(best.trial)$ds.summary$average.p
# Table 2 (Rankin et al. 2020)
p0.beaufort

###################################################################
#### Predicting p0 in the independent observer configuration
#' @param object ddf.io object
#' @param newdata New data for making predictions
my.predict.io <-
function(object,newdata=NULL,compute=FALSE,int.range=NULL,...){

  model <- object
  if(is.null(newdata)){
    xmat <- model$mr$mr$data
  }else{
    compute <- TRUE
    xmat <- newdata
  }
  
  xmat$distance <- 0
  ddfobj <- model$ds$ds$aux$ddfobj
  
  # for gamma models need to find where p(x)=1 (apex), set that as distance
  if(ddfobj$type=="gamma"){
    xmat$distance <- rep(apex.gamma(ddfobj),2)
  }
  
  # calculate ps for each part of the model
  xmat$offsetvalue <- 0
  p.0 <- predict(model$mr,newdata=xmat,integrate=FALSE,compute=compute)
  if(is.null(newdata)){
    pdot <- predict(model$ds,esw=FALSE,compute=compute,
                    int.range=int.range)$fitted
  }else{
    pdot <- predict(model$ds,newdata=newdata[newdata$observer==1,],
                    esw=FALSE,compute=compute,int.range=int.range)$fitted
  }
  
  return(list(p.0 = p.0, pdot = pdot))
}

# example with Beaufort (Table 4 from Rankin et al. 2020) 
bft=c(0:5)
newdata.io<-data.frame(distance=rep(0,12),beaufort=rep(bft,each=2),observer=rep(c(1,2),times=6))
io.pred<-data.frame(matrix(NA,3,(length(bft)+1)))
colnames(io.pred)<-c("obs",bft)
io.pred[,1]<-c("obs1","obs2","both")
io.pred[1,2:(length(bft)+1)]<-my.predict.io(best.io,newdata.io)$p.0$p1
io.pred[2,2:(length(bft)+1)]<-my.predict.io(best.io,newdata.io)$p.0$p2
io.pred[3,2:(length(bft)+1)]<-my.predict.io(best.io,newdata.io)$p.0$fitted
print(io.pred, row.names = F)

```

> **Comments:**
> **username_3**: As it’s currently written, your answer is unclear. Please [edit] to add additional details that will help others understand how this addresses the question asked. You can find more information on how to write good answers [in the help center](/help/how-to-answer).  

---

