Lower Columbia River Rainbow Trout Spawning 2019
The suggested citation for this analytic appendix is:
Thorley, J.L. and Amies-Galonski, E. (2020) Lower Columbia River Rainbow Trout Spawning 2019. A Poisson Consulting Analysis Appendix. URL: https://www.poissonconsulting.ca/f/949693135.
Background
Each spring in the Lower Columbia River (LCR) below Hugh L. Keenleyside Dam (HLK) and in the Lower Kootenay River (LKR) below Brilliant Dam, thousands of Rainbow Trout spawn.
Since 1992, BC Hydro has stabilized the spring discharge releases from HLK to protect Rainbow Trout redds from dewatering.
The primary goal of the current analysis to estimate the abundance of spawners, egg survival and the stock-recruitment relationship. We also examine the survival of dewatered eggs and the temperatures in the gravels.
Data Preparation
The aerial spawner counts were conducted by Mountain Water Research. The age-1 recruitment estimates were provided by the Fish Population Indexing Program. The remaining data were collected by Mountain Water Research, Poisson Consulting and Nuupqu.
The study area was divided into seven sections: Norns Creek Fan (NCF), NCF to LKR, LKR, LKR to Genelle, Genelle. Redd and spawner counts upstream of Norns Creek Fan and downstream of Genelle were excluded from the section totals because they constitute less than 0.1% of the total count and were not surveyed in all years. The redd and spawner counts for the Right Upstream Bank above Robson Bridge were also excluded as they appear to be primarily driven by viewing conditions (and constitute less than 2.5% of the total). Viewing conditions were classified as Good or Poor. If information on the viewing conditions was not available a decline in the redd count of more than one third of the cumulative maximum count for a particular section was assumed to be caused by poor viewing conditions.
The mapped peak redd and fish counts are the peak counts for that site.
The data were prepared for analysis using R version 3.6.2 (R Core Team 2017).
Statistical Analysis
Model parameters were estimated using Bayesian methods. The estimates were produced using JAGS (Plummer 2015) and STAN (Carpenter et al. 2017). For additional information on Bayesian estimation the reader is referred to McElreath (2016).
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).
The parameters are summarised in terms of the point estimate, standard deviation (sd), the z-score, lower and upper 95% confidence/credible limits (CLs) and the p-value (Kery and Schaub 2011, 37, 42). The estimate is the median (50th percentile) of the MCMC samples, the z-score is \(\mathrm{mean}/\mathrm{sd}\) and the 95% CLs are the 2.5th and 97.5th percentiles. A p-value of 0.05 indicates that the lower or upper 95% CL is 0.
The results are displayed graphically by plotting the modeled relationships between particular variables and the response(s) 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 typical values (expected values of the underlying hyperdistributions) (Kery and Schaub 2011, 77–82). When informative the influence of particular variables is expressed in terms of the effect size (i.e., percent change in the response variable) with 95% confidence/credible intervals (CIs, Bradford, Korman, and Higgins 2005).
The analyses were implemented using R version 3.6.2
(R Core Team 2019) and the
mbr
family of packages.
Model Descriptions
Area-Under-The-Curve
The spawner abundance and spawn timings were estimated from the aerial fish and redd counts for the five sections (in three segments) using an Area-Under-The-Curve (AUC) model.
Key assumptions of the AUC model include:
- Spawner abundance varies by river section.
- Spawner abundance varies randomly by year and section within year.
- Spawner observer efficiency is between 0.8 and 1.0.
- Number of redds per spawner is between 1 and 2.
- Spawner residence time is between 14 and 21 days as determined in a previous year’s analysis.
- Redd residence time is between 30 and 40 days.
- Spawner arrival and departure times are normally distributed.
- Spawner arrival duration (SD of normal distribution) varies randomly by river segment.
- Peak spawner arrival timing varies randomly by year.
- The residual variations in the spawner and redd counts are described by separate Negative Binomial distributions.
Stock-Recruitment
The relationship between the number of spawners and the resultant number of age-1 fish was estimated using a Beverton-Holt stock-recruitment model (Walters and Martell 2004):
\[ R = \frac{\alpha \cdot S}{1 + \beta \cdot S} \quad,\]
where \(S\) is the spawners (stock), \(R\) is the recruits, \(\alpha\) is the recruits per spawner at low density and \(\beta\) determines the density-dependence.
Key assumptions of the stock-recruitment model include:
- The recruits per spawner at low density (\(\alpha\)) is normally distributed with a mean of 90 and a SD of 50.
- The recruits per spawner varies with the percent of redds dewatered.
- The residual variation in the number of recruits is log-normally distributed.
The mean of 90 for \(\alpha\) was based on an average of 2,900 eggs per female spawner, a 50:50 sex ratio, 50% egg survival, 50% post-emergence fall survival, 50% overwintering survival and 50% summer survival.
The carrying capacity is \(\alpha / \beta\).
Model Templates
Area-Under-The-Curve
data {
int<lower=0> nObs;
int<lower=0> nSection;
int<lower=0> nSegment;
int<lower=0> nYear;
int Year[nObs];
int Section[nObs];
int Segment[nObs];
real Doy[nObs];
int Fish[nObs];
int Redds[nObs];
parameters {
vector<lower=3,upper=9>[nSection] bFishAbundanceSection;
real<lower=0> sFishAbundanceYear;
vector[nYear] bFishAbundanceYear;
real<lower=0> sFishAbundanceSectionYear;
vector[nSection * nYear] bFishAbundanceSectionYear;
real<lower=0.8,upper=1.0> bFishObserverEfficiency;
real<lower=0.0,upper=2.0> bReddObserverEfficiency;
real<lower=1, upper=2> bReddPerFish;
real<lower=14, upper=21> bFishResidenceTime;
real<lower=30, upper=40> bReddResidenceTime;
real<lower=100, upper=150> bPeakFishArrivalTiming;
real<lower=0,upper=21> sPeakFishArrivalTimingYear;
vector[nYear] bPeakFishArrivalTimingYear;
real<lower=log(10), upper=log(50)> bFishArrivalDuration;
real<lower=0> sFishArrivalDurationSegmentYear;
vector[nSegment * nYear] bFishArrivalDurationSegmentYear;
real<lower=0, upper=2> bDispersionRedds;
real<lower=0, upper=2> bDispersionFish;
model {
vector[nObs] eFishAbundance;
vector[nObs] ePeakFishArrivalTiming;
vector[nObs] eFishArrivalDuration;
vector[nObs] eRedds;
vector[nObs] eFish;
sFishAbundanceYear ~ normal(0, 1);
bFishAbundanceYear ~ normal(0, sFishAbundanceYear);
sFishAbundanceSectionYear ~ normal(0, 1);
bFishAbundanceSectionYear ~ normal(0, sFishAbundanceSectionYear);
bPeakFishArrivalTimingYear ~ normal(0, sPeakFishArrivalTimingYear);
bFishArrivalDurationSegmentYear ~ normal(0, 1);
bFishArrivalDurationSegmentYear ~ normal(0, sFishArrivalDurationSegmentYear);
for (i in 1:nObs) {
eFishAbundance[i] = exp(bFishAbundanceSection[Section[i]] + bFishAbundanceYear[Year[i]] + bFishAbundanceSectionYear[((Section[i] - 1) * Year[i]) + Year[i]]);
ePeakFishArrivalTiming[i] = bPeakFishArrivalTiming + bPeakFishArrivalTimingYear[Year[i]];
eFishArrivalDuration[i] = exp(bFishArrivalDuration + bFishArrivalDurationSegmentYear[((Segment[i] - 1) * Year[i]) + Year[i]]);
eRedds[i] = eFishAbundance[i] * bReddPerFish * bReddObserverEfficiency * (normal_cdf(Doy[i], ePeakFishArrivalTiming[i], eFishArrivalDuration[i]) - normal_cdf(Doy[i], ePeakFishArrivalTiming[i] + bReddResidenceTime, eFishArrivalDuration[i]));
eFish[i] = eFishAbundance[i] * bFishObserverEfficiency * (normal_cdf(Doy[i], ePeakFishArrivalTiming[i], eFishArrivalDuration[i]) - normal_cdf(Doy[i], ePeakFishArrivalTiming[i] + bFishResidenceTime, eFishArrivalDuration[i]));
}
Redds ~ neg_binomial_2(eRedds, 1/bDispersionRedds);
Fish ~ neg_binomial_2(eFish, 1/bDispersionFish);
..
Block 1. Model description.
Stock-Recruitment
.model {
bAlpha ~ dnorm(90, 50^-2) T(1,)
bAlphaDewatered ~ dnorm(0, 2^-2)
bBeta ~ dnorm(-5, 5^-2)
sRecruits ~ dunif(0, 2)
for(i in 1:length(Stock)) {
log(eAlpha[i]) <- log(bAlpha) + bAlphaDewatered * Dewatered[i]
log(eBeta[i]) <- bBeta
eRecruits[i] <- eAlpha[i] * Stock[i] / (1 + eBeta[i] * Stock[i])
Recruits[i] ~ dlnorm(log(eRecruits[i]), sRecruits^-2)
}
..
Block 2. Model description.
Results
Tables
Area-Under-The-Curve
Table 1. Parameter descriptions.
Parameter | Description |
---|---|
bDispersionFish |
Clustering parameter for Fish |
bDispersionRedds |
Clustering parameter for Redds |
bFishAbundanceSection[i] |
Intercept for log(eFishAbundance) by
Section |
bFishAbundanceSectionYear[i] |
Effect of Section by Year on
bFishAbundanceSection |
bFishAbundanceYear[i] |
Effect of i th Year on
bFishAbundanceSection |
bFishArrivalDuration[i] |
Intercept for log(eFishArrivalDuration) |
bFishArrivalDurationSegmentYear[i] |
Effect of Segment by Year on
bFishArrivalDuration |
bFishObsEfficiency |
Fish observer efficiency |
bFishResidenceTime |
Fish residence time (days) |
bPeakFishArrivalTiming |
Intercept for ePeakFishArrivalTiming |
bPeakFishArrivalTimingYear[i] |
Effect of i th Year on
bPeakFishArrivalTiming |
bReddObserverEfficiency |
Redd observer efficiency |
bReddPerFish |
Number of Redds per Fish |
bReddResidenceTime |
Redd residence time (days) |
Doy[i] |
Day of the year of i th count |
eFishAbundance[i] |
Expected Fish abundance for i th count |
eFishArrivalDuration[i] |
Expected SD of Fish arrival timing for
i th count |
ePeakFishArrivalTiming[i] |
Expected Doy of peak Fish arrival for
i th count |
Fish[i] |
Observed number of Fish on i th count |
Redds[i] |
Observed number of Redds on i th count |
Section[i] |
Section of i th count |
Segment[i] |
Segment of i th count |
sFishAbundanceSectionYear |
SD of bFishAbundanceSectionYear |
sFishAbundanceYear |
SD of bFishAbundanceYear |
sFishArrivalDurationSegmentYear |
SD of bFishArrivalDurationSegmentYear |
sPeakFishArrivalTimingYear |
SD of bPeakFishArrivalTimingYear |
Year[i] |
Year of i th count |
Table 2. Model coefficients.
term | estimate | sd | zscore | lower | upper | pvalue |
---|---|---|---|---|---|---|
bDispersionFish | 0.4906261 | 0.0289764 | 16.970361 | 0.4352173 | 0.5505956 | 0.0006662 |
bDispersionRedds | 0.1962949 | 0.0135124 | 14.570889 | 0.1728795 | 0.2247813 | 0.0006662 |
bFishAbundanceSection[1] | 4.5256828 | 0.2291133 | 19.760193 | 4.0933065 | 4.9916457 | 0.0006662 |
bFishAbundanceSection[2] | 7.3017620 | 0.2253793 | 32.417061 | 6.8831049 | 7.7551923 | 0.0006662 |
bFishAbundanceSection[3] | 6.6582501 | 0.2320588 | 28.729528 | 6.2208770 | 7.1547150 | 0.0006662 |
bFishAbundanceSection[4] | 7.2087913 | 0.2321086 | 31.054027 | 6.7678231 | 7.6821928 | 0.0006662 |
bFishAbundanceSection[5] | 7.2121764 | 0.2354190 | 30.656207 | 6.7498749 | 7.7203229 | 0.0006662 |
bFishArrivalDuration | 3.1720548 | 0.0347185 | 91.350599 | 3.1029772 | 3.2391417 | 0.0006662 |
bFishObserverEfficiency | 0.9048826 | 0.0580127 | 15.569918 | 0.8047635 | 0.9944007 | 0.0006662 |
bFishResidenceTime | 19.8487449 | 1.2174736 | 16.059572 | 16.5137296 | 20.9577506 | 0.0006662 |
bPeakFishArrivalTiming | 118.2296004 | 2.3068480 | 51.265016 | 113.9813952 | 122.9932436 | 0.0006662 |
bReddObserverEfficiency | 0.6007467 | 0.1320291 | 4.640964 | 0.4076091 | 0.8855104 | 0.0006662 |
bReddPerFish | 1.4209319 | 0.2883805 | 5.028460 | 1.0180470 | 1.9639243 | 0.0006662 |
bReddResidenceTime | 31.4750506 | 1.7308721 | 18.468174 | 30.0672955 | 36.3407426 | 0.0006662 |
sFishAbundanceSectionYear | 0.3788588 | 0.0449997 | 8.497811 | 0.3064768 | 0.4851855 | 0.0006662 |
sFishAbundanceYear | 0.7994207 | 0.1513213 | 5.414558 | 0.5693903 | 1.1641755 | 0.0006662 |
sFishArrivalDurationSegmentYear | 0.1599298 | 0.0248510 | 6.536656 | 0.1194830 | 0.2165317 | 0.0006662 |
sPeakFishArrivalTimingYear | 8.3716831 | 1.7103795 | 5.038135 | 5.9050545 | 12.5903873 | 0.0006662 |
Table 3. Model summary.
n | K | nchains | niters | nthin | ess | rhat | converged |
---|---|---|---|---|---|---|---|
858 | 18 | 3 | 500 | 1 | 251 | 1.016 | TRUE |
Stock-Recruitment
Table 4. Parameter descriptions.
Parameter | Description |
---|---|
bAlpha |
Intercept for eAlpha |
bAlphaDewatered |
Effect of Dewatered on log(bAlpha) |
bBeta |
Intercept for log(eBeta) |
Dewatered[i] |
Proportional redd dewatering in i th spawn year |
eAlpha |
Expected number of recruits at low density |
eBeta |
Expected density-dependence |
eRecruits[i] |
Expected Recruits |
Recruits[i] |
Number of age-1 recruits from i th spawn year |
sRecruits |
SD of residual variation in Recruits |
Stock[i] |
Number of spawners in i th spawn year |
Table 5. Model coefficients.
term | estimate | sd | zscore | lower | upper | pvalue |
---|---|---|---|---|---|---|
bAlpha | 101.8179211 | 41.7508757 | 2.529679 | 32.7516372 | 195.0634631 | 0.0006662 |
bAlphaDewatered | 0.2096527 | 0.0669056 | 3.105676 | 0.0760181 | 0.3329382 | 0.0086609 |
bBeta | -5.4231079 | 0.4760855 | -11.508782 | -6.6290218 | -4.7315143 | 0.0006662 |
sRecruits | 0.2550770 | 0.0532394 | 4.945341 | 0.1853825 | 0.3949024 | 0.0006662 |
Table 6. Model summary.
n | K | nchains | niters | nthin | ess | rhat | converged |
---|---|---|---|---|---|---|---|
18 | 4 | 3 | 500 | 1000 | 1118 | 1 | TRUE |
Dewatering
Table 7. Reduction dates, magnitude of reduction, number and general location of dewatered redds in 2019.
Reduction Date | HLK Discharge Start (m3/s) | HLK Discharge End (m3/s) | BRD Discharge Start (m3/s) | BRD Discharge End (m3/s) | Location | Dewatered Redds |
---|---|---|---|---|---|---|
2019-02-22 | 1354 | 734 | 400 | 398 | NA | 0 |
2019-03-08 | 1289 | 494 | 425 | 459 | Genelle Ch. E | 7 |
2019-03-08 | 1289 | 494 | 425 | 459 | Norns Creek Fan | 9 |
2019-05-15 | 752 | 430 | 1289 | 1326 | Norns Creek Fan | 3 |
2019-05-18 | 426 | 297 | 1483 | 1479 | RUB, d/s Waldie Island | 3 |
2019-05-18 | 426 | 297 | 1483 | 1479 | Norns Creek Fan | 71 |
2019-06-03 | 649 | 341 | 1804 | 1842 | Norns Creek Fan | 8 |
Figures
Maps
Sensors
Spawning
Sensor Data
Real Time Stations
Gravel Temperature Stations
Relative Gravel Temperature Stations
Norns Creek Fan
The Oxbow
Absolute Gravel Temperature Stations
Norns Creek Fan
The Oxbow
Egg Mortality
Area-Under-The-Curve
Norns Creek Fan
NCF To LKR
Lower Kootenay River
LKR To Genelle
Genelle
Stock-Recruitment
Norns Creek
Acknowledgements
The organisations and individuals whose contributions have made this analysis report possible include:
- BC Hydro
- Philip Bradshaw
- James Baxter
- Guy Martel
- Margo Dennis
- Darin Nishi
- Nuupqu
- Mark Fjeld
- Natalie Morrison
- Mountain Water Research
- Jeremy Baxter
- Poisson Consulting
- Robyn Irvine
- Dam Helicopters
- Duncan Wassick
- Highland Helicopters
- Phil Hocking
- Mark Homis
- Golder Associates
- Dustin Ford
- Demitria Burgoon
- David Roscoe
- Additional Support
- Clint Tarala
- Crystal Lawrence
- Gary Pavan
- Gerry Nellestijn
- Jay Bowers
References
Bradford, Michael J, Josh Korman, and Paul S Higgins. 2005. “Using Confidence Intervals to Estimate the Response of Salmon Populations (Oncorhynchus Spp.) to Experimental Habitat Alterations.” Canadian Journal of Fisheries and Aquatic Sciences 62 (12): 2716–26. https://doi.org/10.1139/f05-179.
Brooks, Steve, Andrew Gelman, Galin L. Jones, and Xiao-Li Meng, eds. 2011. Handbook for Markov Chain Monte Carlo. Boca Raton: Taylor & Francis.
Carpenter, Bob, Andrew Gelman, Matthew D. Hoffman, Daniel Lee, Ben Goodrich, Michael Betancourt, Marcus Brubaker, Jiqiang Guo, Peter Li, and Allen Riddell. 2017. “Stan : A Probabilistic Programming Language.” Journal of Statistical Software 76 (1). https://doi.org/10.18637/jss.v076.i01.
Gelman, Andrew, Daniel Simpson, and Michael Betancourt. 2017. “The Prior Can Often Only Be Understood in the Context of the Likelihood.” Entropy 19 (10): 555. https://doi.org/10.3390/e19100555.
Kery, Marc, and Michael Schaub. 2011. Bayesian Population Analysis Using WinBUGS : A Hierarchical Perspective. Boston: Academic Press. http://www.vogelwarte.ch/bpa.html.
McElreath, Richard. 2016. Statistical Rethinking: A Bayesian Course with Examples in R and Stan. Chapman & Hall/CRC Texts in Statistical Science Series 122. Boca Raton: CRC Press/Taylor & Francis Group.
Plummer, Martyn. 2015. “JAGS Version 4.0.1 User Manual.” http://sourceforge.net/projects/mcmc-jags/files/Manuals/4.x/.
R Core Team. 2017. “R: A Language and Environment for Statistical Computing.” Vienna, Austria: R Foundation for Statistical Computing. https://www.R-project.org/.
———. 2019. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing. https://www.R-project.org/.
Walters, Carl J., and Steven J. D. Martell. 2004. Fisheries Ecology and Management. Princeton, N.J: Princeton University Press.