install.packages("Rcapture")Lab 12 – Estimating abundance with closed-population capture-mark-recapture data
Due by 5:00pm on Friday
The purpose of this lab is to learn how to estimate abundance using mark-recapture data from surveys of closed populations. A closed population does not experience recruitment, mortality, immigration or emigration. The closure assumption is usually only valid over very short time periods. We will learn how to work with data from open populations later.
Work through the problems below, put your code and your answers in a Word file, and then upload it to ELC. Feel free to use R Markdown or Quarto if you wish.
Assignment
Exercise I: Lincoln-Peterson estimation
Suppose you capture, mark, and release 100 largemouth bass (Micropterus salmoides) at Lake Herrick. The next day, you return and capture 50 individuals, 25 of which were marked on the first occasion. What is the Lincoln-Peterson estimate of abundance (\(N\))? Do the calculations in R and show your work.
Exercise II: Closed-population models in R
- Use the R package “Rcapture” to fit mark-recapture models, including those shown the table below, to the stinkpot (Sternotherus odoratus) data found in the file
CH-SO-Andy07.inp. These data were collected using baited “hoop traps” at Andy’s Pond in 2007. Software instructions are in the next section. - The model with the lowest AIC is considered the best in the set of models. Which model has the lowest AIC? Do your results suggest that there was a behavioral (i.e., a trap-happy or trap-shy) response? Explain.
- Why do you think the best model was superior to \(M_0\)? Look at the capture histories when answering this question.
- Report the parameter estimates from the best model. Include the standard error for the abundance estimate.
| Model | Description |
|---|---|
| \(M_0\) | The most basic model in which \(p\) and \(c\) are constant |
| \(M_t\) | \(p\) differs among sampling occasions and \(p_t=c_t\). |
| \(M_b\) | Behavioral response model in which \(p\) and \(c\) differ. Can describe trap happiness or trap shiness. |
Parameter definitions
- \(p\) – capture probability. The probability of capturing an individual on a single occasion
- \(p_t\) – capture probability on occasion \(t\)
- \(c\) – recapture probability. The probability of capturing an individual that has been captured previously.
- \(n\) – the number of individuals captured
- \(N\) – abundance. The number of individuals in the population.
Example analysis in R
Open R and install the “Rcapture” package using the following command:
Now load the package like so:
library(Rcapture)The example data are formatted for program MARK, but we can import them using read.table. The only trick is to tell R that the capture histories should be treated as a character string, rather than as a numeric variable. The colClasses arguments let’s us do that.
capture.histories <- read.table("CH-example.inp", sep=" ",
colClasses=c("character", "character"),
col.names=c("ch", "freq"))It’s a small dataset, so we can display it in full:
capture.histories ch freq
1 101001 1;
2 100001 1;
3 100001 1;
4 010000 1;
5 101000 1;
6 010001 1;
7 001000 1;
8 001011 1;
9 000100 1;
10 001001 1;
11 001000 1;
12 000010 1;
13 000010 1;
14 000101 1;We need to convert the capture histories from a character vector to a matrix. The following code does the trick.
ch.mat <- t(sapply(capture.histories$ch,
function(x) as.integer(strsplit(x, "")[[1]])))
dimnames(ch.mat) <- list(paste0("Animal", 1:nrow(capture.histories)),
paste0("Time", 1:6))Now the capture histories look like this:
ch.mat Time1 Time2 Time3 Time4 Time5 Time6
Animal1 1 0 1 0 0 1
Animal2 1 0 0 0 0 1
Animal3 1 0 0 0 0 1
Animal4 0 1 0 0 0 0
Animal5 1 0 1 0 0 0
Animal6 0 1 0 0 0 1
Animal7 0 0 1 0 0 0
Animal8 0 0 1 0 1 1
Animal9 0 0 0 1 0 0
Animal10 0 0 1 0 0 1
Animal11 0 0 1 0 0 0
Animal12 0 0 0 0 1 0
Animal13 0 0 0 0 1 0
Animal14 0 0 0 1 0 1We will use the closedp function to fit the closed population models. This returns abundance estimates, standard errors (stderr), and AIC values. You can ignore the models that we didn’t cover in class. Focus primarily on M0, Mt, and Mb.
mods <- closedp(ch.mat)
mods
Number of captured units: 14
Abundance estimations and model fits:
abundance stderr deviance df AIC BIC infoFit
M0 18.0 3.1 29.771 61 57.612 58.890 OK
Mt 17.4 2.8 22.692 56 60.533 65.007 OK
Mh Chao (LB) 18.0 3.1 29.771 61 57.612 58.890 OK
Mh Poisson2 15.6 2.3 28.704 60 58.545 60.462 OK
Mh Darroch 15.0 2.0 28.836 60 58.677 60.594 OK
Mh Gamma3.5 14.7 1.6 28.896 60 58.737 60.654 OK
Mth Chao (LB) 17.4 2.8 22.692 56 60.533 65.007 OK
Mth Poisson2 15.5 2.2 21.839 55 61.681 66.793 OK
Mth Darroch 15.0 2.0 21.962 55 61.803 66.915 OK
Mth Gamma3.5 14.7 1.7 22.015 55 61.856 66.969 OK
Mb 17.2 4.4 29.730 60 59.572 61.489 OK
Mbh 16.7 6.6 29.698 59 61.539 64.095 OKYou can extract parameter estimates for a model like this:
round(mods$parameters$M0, digits=3) N p
estimate: 17.962 0.223For model \(M_0\), capture probability is the same for all individuals during all time periods. This won’t be the case for the other models. Here are the estimates for model \(M_t\):
round(mods$parameters$Mt, digits=3) N p1 p2 p3 p4 p5 p6
estimate: 17.404 0.23 0.115 0.345 0.115 0.172 0.402And model \(M_b\):
round(mods$parameters$Mb, digits=3) N p c
estimate: 17.151 0.246 0.217We can get standard errors (stderr) for abundance estimates using the function closedpCI.t:
Mt <- closedpCI.t(ch.mat, m="Mt")
Mt$results # Abundance and standard error abundance stderr deviance df AIC BIC infoFit
Mt 17.40413 2.846915 22.69221 56 60.53333 65.00673 0