The suggested citation for this analytic appendix is:
Hill, N.E., Pearson, A.D., & Thorley, J.L. (2023) Elk River WCT Abundance Analysis 2022. A Poisson Consulting Analytic Appendix. URL: https://www.poissonconsulting.ca/f/876409566.
The Elk River is one of eight upper Kootenay River tributaries (Bull, Michel, Skookumchuck, St. Mary, upper Kootenay, White and Wigwam) that are designated as Class II due to the importance of their recreational fisheries for Westslope Cutthroat Trout (Oncorhynchus clarkii lewisi). The population abundance of Westslope Cutthroat Trout in the Elk River is uncertain. To improve population estimates a mark-recapture analysis using Passive Integrated Transponder (PIT) tags deployed by both boat electrofishing and guides was recommended (J. L. Thorley 2021).
After a successful pilot study in 2021, data collection via boat electrofishing and guides continued in 2022. The goal of the present analysis is to use the additional data to get an improved estimate of Westslope Cutthroat Trout population abundance in the Elk River.
The data were provided by Nupqu Limited Partnership. The data were prepared for analysis using R version 4.2.2 (R Core Team 2022) and organized in a SQLite database. All river distances are based on the BC Freshwater Atlas layer.
Key assumptions of the data preparation included:
Model parameters were estimated using Bayesian methods. The estimates were produced using JAGS (Plummer 2003) and STAN (Carpenter et al. 2017). For additional information on Bayesian estimation the reader is referred to McElreath (2020).
Unless stated otherwise, the Bayesian analyses used weakly informative normal and half-normal prior distributions (Gelman, Simpson, and Betancourt 2017). The posterior distributions were estimated from 1500 Markov Chain Monte Carlo (MCMC) samples thinned from the second halves of 3 chains (Kery and Schaub 2011, 38–40). Model convergence was confirmed by ensuring that the potential scale reduction factor \(\hat{R} \leq 1.05\) (Kery and Schaub 2011, 40) and the effective sample size (Brooks et al. 2011) \(\textrm{ESS} \geq 150\) for each of the monitored parameters (Kery and Schaub 2011, 61).
Model adequacy was assessed via posterior predictive checks (Kery and Schaub 2011). More specifically, the proportion of zeros in the data and the first four central moments (mean, variance, skewness and kurtosis) in the deviance residuals were compared to the expected values by simulating new data based on the posterior distribution and assumed sampling distribution and calculating the deviance residuals.
Where computationally practical, the sensitivity of the posteriors to the choice of prior distributions was evaluated by increasing the standard deviations of all normal, half-normal and log-normal priors by an order of magnitude and then using \(\hat{R}\) to evaluate whether the samples were drawn from the same posterior distribution (Joseph L. Thorley and Andrusak 2017).
The parameters are summarised in terms of the point estimate, lower and upper 95% compatibility limits (Rafi and Greenland 2020) and the surprisal s-value (Greenland 2019). The estimate is the median (50th percentile) of the MCMC samples while the 95% CLs are the 2.5th and 97.5th percentiles. The s-value indicates how surprising it would be to discover that the true value of the parameter is in the opposite direction to the estimate (Greenland 2019). An s-value of \(>\) 4.32 bits, which is equivalent to a p-value \(<\) 0.05 (Kery and Schaub 2011; Greenland and Poole 2013), indicates that the surprise would be equivalent to throwing at least 4.3 heads in a row.
Model selection was based the heuristic of directional certainty (Kery and Schaub 2011). Fixed effects were included if their s-value was \(>\) 4.32 bits (Kery and Schaub 2011). Based on a similar argument, random effects were included if their standard deviation had lower 95% CLs \(>\) 5% of the estimate.
The results are displayed graphically by plotting the modeled relationships between individual variables and the response with the remaining variables held constant. In general, continuous and discrete fixed variables are held constant at their mean and first level values, respectively, while random variables are held constant at their average values (expected values of the underlying hyperdistributions) (Kery and Schaub 2011, 77–82). Primary explanatory variables are evaluated based on their estimated effect sizes with 95% CLs (Bradford, Korman, and Higgins 2005)
The analyses were implemented using R version 4.2.2
(R Core Team 2022) and the
mbr family of packages.
The data were analysed using a hierarchical Bayesian mark-recapture abundance model. Each zone was split into two kilometres sites except sites at the end of the zone were added to the previous site if the length was less than one kilometre. Captures from sites that had less then 10% of their length surveyed were not included in the analysis.
Key assumptions of the mark-recapture model include:
.model{
bDensity ~ dnorm(6, 2^-2)
bEfficiency ~ dnorm(-5, 2^-2)
bEfficiencySecondsPerMetre ~ dnorm(0, 1^-2)
bTheta ~ dexp(10)
bDensityAnnual[1] <- 0
for (i in 2:nannual) {
bDensityAnnual[i] ~ dnorm(0, 2^-2)
}
sDensityAnnualZoneID ~ dnorm(0, 1^-2) T(0,)
for (i in 1:nannual) {
for (j in 1:nzone_id) {
bDensityAnnualZoneID[i, j] ~ dnorm(0, sDensityAnnualZoneID^-2)
}
}
bEfficiencyAnnual[1] <- 0
for (i in 2:nannual) {
bEfficiencyAnnual[i] ~ dnorm(0, 2^-2)
}
for (i in 1:nObs) {
log(eDensity[i]) <- bDensity + bDensityAnnual[annual[i]] + bDensityAnnualZoneID[annual[i], zone_id[i]]
logit(eProb[i]) <- bEfficiency + bEfficiencySecondsPerMetre * seconds_per_metre[i] + bEfficiencyAnnual[annual[i]]
eEfficiency[i] ~ dbeta((2 * eProb[i]) / bTheta, (2 * (1 - eProb[i])) / bTheta)
recaptures[i] ~ dbin(eEfficiency[i], marked[i])
bAbundance[i] ~ dpois(eDensity[i] * site_length[i] * survey_proportion[i])
count[i] ~ dbin(eEfficiency[i], bAbundance[i])
}
Block 1. Model description.