3. The model
3.2. Observation models
3.2.1 Electrofishing process
Electrofishing is a removal sampling method, based on the fishes’ sensitivity towards
electromagnetic fields and the paralysis caused be the electromagnetic current (Halsband 1967).
The effect of the electric current is also influenced by the conductivity of the water (Alabaster &
Hartley 1962).
The typical sampling procedure is carried out by two to three persons. One of the samplers is carrying the generator on their back and operating the anode staff, which is used to attract and stun the fish. The others follow the anode operator, using a handheld net to scoop the stunned fish in to a water filled container, where the fish are stored alive for further handling. (Taylor et al. 2002) In the sampling process parr are assumed to be caught individually with certain probability of capture, i.e. catchability (w). With this assumption the amount of parrs in a sampling area can be estimated with a binomial distribution, where the number of trials (N) is the true amount of fish in the sampled area (formula 18).
(18) In order to facilitate the MCMC-simulation the binomial distribution of catchy,a,r is approximated using a Poisson distribution (formula 19).
(19)
34
Catchability is typically estimated based on the declining amount of captured fish per removal, and it is assumed to be equal for all fish in the sampled population and not varying between successive removals (Bohlin 1983). As pointed out by Mäntyniemi et al. (2005b) this can result in an
underestimation of the population size, which can be corrected by using Bayesian models with unequal catchability estimates.
In this model the catchability is assumed to be correlated with the fish’s length. In theory the fishes affinity to swim towards the electrical current increases with fish size. Thus the longer the fish, the easier it is to catch it with electrofishing gear (Reynolds 1996; Reynolds & Simpson 1978; Vibert 1967). This has been found to hold true for brown trout, especially in streams with low conductivity (Borgstrom & Skaala 1993). While this relationship may be true in theory, in practice the larger fish sense the electrical current farther away, and tend to be spooked away from the sampling area. Very large individuals are also too big to be scooped up with the sampling net. Consequently their
existence is observed, but they are not recorded in the data.
In order to model the different effects that the fishing gear has on larger individuals, two logistic regression lines are used as priors for the catchability of individual fish. The first regression line was fitted to have increasing catchability values with increasing length, and the other one was fitted with decreasing values, with increasing length. Both models were then weighted using Bayesian model averaging (Gelman et al. 2014; Carlin & Chip 1995), resulting in final catchability estimate being the weighted average of the two models (formula 20, figure 13).
35
(20) Typically age has a strong positive correlation with the length of the fish, especially in young
juvenile fish. For the sake of simplicity age was used as an estimate for fish length, but in
populations exiting stunted growth, this not always true (Jobling & Reinsnes 1986). In the datasets used for this study there was no evidence of stunted growth.
In this model all yearly electrofishing surveys within a river are treated as one single-removal trial, with different catchabilities between age-groups. This parameterizing of the model is justified by the lower amount of model parameters, and the fact that the mean values used for the catchability estimate are already based on declining catches between successive removals (Bohlin 1983).
Figure 13. Two opposite priors for age-specific catchability (w) in electrofishing experiments.
36 3.2.2. Smolt trapping
Smolt trapping is typically done using a fyke net or a rotary trapping device, usually referred to as a
‘smoltscrew’. In the dataset used in this study both trap types were used (see section 2.3.). In the smolt trapping experiments the trap was placed beneath the lowest rapids, or in the river mouth, where there is sufficient current.
Salmonids tend to lose their territorial behavior after smoltification and gather into schools for their downstream migration (Hoar 1988). According to Svendsen et al. (2007) the smolts swim in the center and near the bottom, but other studies suggest that salmonid smolts swim near the surface (Moore et al. 1998; Davidsen et al. 2005). According to Jonsson & Jonsson (2011) orientation of smolts with the current depends on the velocity of the current. In fast currents the smolts move downstream tail first, head facing the current. In slower currents the smolts actively swim
downwards head facing downstream (Jonsson & Jonsson 2011). This would suggest that the smolts are unable to detect the trap in areas of fast current velocities.
In this model each fish is assumed to enter the trap individually and therefore the process can be modeled as a binomially distributed variable, where the true number of smolts leaving the river is the amount of independent trials (formula 21):
(21) In order to ease the MCMC simulation, the computationally difficult binomial distribution of the smolt-trapping processes is also approximated using a Poisson distribution (formula 22):
(22) The catchability of an independent smolt was given a beta-distribution, with a simple mean value derived from the trap’s width per river’s width ratio (formula 23):
37 The catchability in smolt trapping is typically estimated using the Petersen-Lincoln method
(Petersen 1896; Lincoln 1930), in which some or all of the fish caught in the trap are marked with a visual tag, and then transported upstream to be released back into the river. Different batches of fish tagged on different dates are differentiated from one another by using different colored tags, or reshaping the tag, when using a T-anchor tag. The catchability is then estimated to be the ratio between the number of re-captured individuals and the number of fish released. The findings of Moore et al. (1998) and Davidsen et al. (2005) suggest that the use of mandarins for estimating catchability in the River Ingarskila data is justified (see section 2.3.1.).
In this model the tagging experiment is defined as a binomially distributed process with the same chance of success as in the trapping procedure (formula 24):
(24) The size of the tagged population is assumed to be constant in the Petersen-Lincoln model, but in reality certain tagging methods can harm the fish causing some of them to perish before they can leave the river. This is taken into account by adding an extra variable Tsurv into the model. Tsurv is calculated as the proportion of tagged (Tagged) fish surviving the tagging procedure (surv).
Probability of survival is given a loose beta-distributed prior, with the expected value decided river specifically, depending on the methods used by researchers collecting the data. This way the binomial distribution describing the process takes the form:
38
(25)
3.2.3. Model for previously published data
The previously published datasets (appendix C.) of observed mean smolt ages (Χmean age) are fed into the model as observations of the weighted average of smolt ages (µ mean age) in every given river.
The confidence intervals reported by L’Abée-Lund et al. (1989) are also included into the analysis as variation in the observation processes. Thus Χ Smolt mean age (formula 26) can be seen as a normally-distributed variable with an expected value µ Smolt age (formula 27)and standard deviation σ Χ Smolt mean age (formula 28): For rives, where the 95 % confidence interval was not reported and the study rivers σ Χ mean age was given a flat prior between 0 and 1 (formula 29):
(29)
39
This parameterizing forms a hierarchical link from the previous datasets, into the population dynamics model. Since the µr parameter in the smoltification model is assumed to correlate with latitude, the model is able to adjust the peak age of smoltification based on the observed mean age in every given river. If mortality was assumed to be zero, then:
40 4. Results
4.1. Model checking with simulated data
In order to validate that the model works as intended it was run with an artificial 16 year dataset, created by running a fully deterministic version of the model. In the model version used for model checking the relationship between smoltification and latitude (see section 3.3) was ignored and the models for Petersen-Lincoln tag-recapture and historical observations were left out.
The priors for catchabilities in the observation models were assigned with wide beta distributions, with expected value of 0.2. The x-parameter controlling the amount of new-born individuals was given a wide log.normally distributed prior with mode values corresponding to the fixed values used for data simulation (40 and 30), and standard deviations corresponding to 24.7 on the original scale. The idea was to give these key-parameters uninformative priors, so that the model would have to learn them almost only based on the information gained from the observations.
The MCMC simulation using the simulated data was first run with two chains for 5.0x105 iterations per chains, while saving every 1 000th sample. This took approximately 66 minutes on Windows 7 computer, fitted with 3.4 GHz Intel i7-4770 processor and 32 GB of memory. The chains seemed to be somewhat mixed, but heavily autocorrelated (figure 14).
41
Figure 14. Convergence diagnostics catchability after 5.0x105 iterations.
The MCMC simulation was re-run with two chains, now for 2.0 x107 iterations per chain while saving every 20 000th sample, resulting in a sample size of 1000 iterations. This took approximately 45 hours on Windows 7 computer, fitted with 3.4 GHz Intel i7-4770 processor and 32 GB of
memory.
First 5.0 x106 iterations were discarded in order to eliminate the burn-in phase of the simulation.
After this the simulation was deemed as converged, even though the running means of the two chains were not completely converged in some parameters. The simulation seemed to struggle especially with the x - parameter controlling the 0+ parr densities in the river and the catchability (w ) parameters in the electrofishing model (figures 15 and 16).
42
Figure 15. Convergence diagnostics and posterior density of x - parameter, depicting slightly different running means (center right graph) in the two simulation chains.
Figure 16. Convergence diagnostics and posterior density of w - parameter, depicting slightly different running means (center right graph) in the two simulation chains.
43
The convergence diagnostics for these parameters also suggested that the chains were still
somewhat autocorrelated (top right graph in figures 15 and 16). This was interpreted as a result of inefficient sampling by the algorithms used by JAGS. Further simulation with higher thinning values would have significantly increased the time needed for simulation. Since the posterior distributions were already centered near the true values used for data creation and the differences between the running means were considered small, further simulation was thought to be
unnecessary.
Figures 17-20 show the prior and the posterior distributions of the target variables. Red lines indicate the true known value of the parameter, used to simulate the data. In the graphs the solid black lines are the mean, and the dashed lines the median values of the posterior distribution (invisible lines indicate that the values are equal).
Figure 17. Prior (blue) and posterior distributions (grey) of smolt trapping catchability. Vertical lines indicate posterior mean (black), posterior median (dashed) and the true value (red).
44
Figure 18. Prior (blue) and posterior distributions (grey) of electrofishing catchability. Vertical lines indicate posterior mean (black), posterior median (dashed) and the true value (red).
45
Figure 19 Prior (blue) and posterior distributions (grey) of 0+ parr density. Vertical lines indicate posterior mean (black), posterior median (dashed) and the true value (red).
46
Figure 20 Prior (blue) and posterior distributions (grey) of σ and τ. Vertical lines indicate posterior mean (black), posterior median (dashed) and the true value (red).
Pictures 21 and 22 show box-plot graphs for the posterior distributions of the total amount of smolts and parrs. The red line indicates the true amounts obtained from the deterministic model. Dots are simulated values grated than 3/2 times the upper quartile, interpreted as outliners by R’s boxplot
47 function. These values are not actual outliner.
Figure 21. Posterior distributions of total size of the parr stock at the beginning of each year. Red line indicates the true value obtained from the deterministic model. Dots represent simulated values 3/2 times higher than the upper quartile, interpreted as outliners by the plotting function. The whiskers represent lowest and highest values, excluding outliners.
48
Figure 22. Posterior distributions of total yearly smolt run. Red line indicates the true value obtained from the deterministic model. Dots represent simulated values 3/2 times higher than the upper quartile, interpreted as outliners by the plotting function. The whiskers represent lowest and highest values, excluding outliners.
4.2. Analysis of the River Pirita and River Ingarskila datasets
The trout stocks of the study rivers was simulated for a 16 year period (1999-2014) with
electrofishing data from 2005 to 2013 for River Pirita and 2006 to 2013 for River Ingarskila, and smolt trapping data from 2009 to 2013 for River Pirita, and 2012 to 2013 for River Ingarskila. The
49
trout stocks of the 41 rivers providing the previously published observations were simulated for a 16 year period, without electrofishing and smolt trapping observations.
The River Pirita fyke net smolt catch of 2014 was omitted from the analysis. In order to assess the predictive performance of the model the simulated observations for 2014 were compared with the observed true data.
In order to help the reader understand the connections between different variables and better interpret the results presented in the following chapters, a graphical summary of the model is provided in appendix A.
4.2.1. MCMC simulation
The MCMC simulation was run with two chains for 1.0 x106 iterations while saving every 1 000th sample, resulting in a sample size of 1 000 iterations per chain. This took approximately 96 hours on Windows 7 computer, fitted with 3.4 GHz Intel i7-4770 processor and 32 GB of memory. A burn-in phase was evident in some parameters (bottom graphs of figures 23 and 24), because of which the first 40 % of the iterations were discarded, resulting in a sample of 600 iterations per chains.
50
Figure 23. Simulation diagnostics of parameter alpha_pii_Pirita.
Figure 24. Simulation diagnostics of alpha_p.
51 Figure 25. Simulation diagnostics of peak_p_Pirita.
4.3.2. Posterior densities
Figures 26 and 27 show the posterior densities for the estimated yearly amounts of brown trout smolts migrating downwards from both study rivers. The orange lines mark the credible intervals of the estimates (Gelman et al. 2014). The smolt run estimates in the beginning of the series (1999-2000) were extremely high in River Pirita and extremely low in River Ingarskila, due to the high uncertainty. The highest simulated value in the time series was 119092for in both rivers.
In River Pirita between the years 2001 and 2011 the estimates for the smolt run are somewhat leveled between 300 and 400 individuals per year. Between 2011 and 2014 the numbers of downward migrating smolts to drop between 200 and 100 individuals.
In River Ingarskila the trend seems to be reversed at first, with initial smolt numbers staying between 200 and just under 100 individuals. Between 2001 and 2003 there is very high uncertainty in the smolt run estimates. Between 2004 and 2010 the estimates stabilize between circa 400 and 300 individuals per year. In 2011 and 2012 there is a sharp drop, with 95 % intervals ranging
52
between circa 150 and almost 0 individuals. In 2013 and 2014 the number smolts is estimated to be around 100 and 200 individuals.
Figure 26. 200 draws from the posterior density of yearly smolt run in River Pirita.
53
Figure 27. 200 draws from the posterior density of yearly smolt run in River Ingarskila.
Figure 28 show the posterior distribution for the density 0+ trout parrs per 100 m2 in both study rivers. The yearly mean values of 0+ parr density varied between 22.14and 48 individuals per 100 m2 in River Pirita and between 22 and 48 in River Ingarskila.
54
Figure 28. Posterior distribution of yearly 0+ trout parrs in the autumn period in River Pirita and Ingarskila.
Figure 29 shows the prior (red lines) and posterior (blue lines) densities of ρ, the variable describing trout’s affinity to migrate to sea at certain age in both study rivers, in the absence of mortality. In river Ingarskila the highest probability seems to be concentrated around age 3, and around age 4 in River Pirita. There is visibly more uncertainty remaining in the posterior distribution of River Ingarskila.
55
Figure 29.Prior and posterior distributions for age-specific probability of smoltification (ρ) in the absence of mortality in study rivers.
Figure 30 shows the prior (red lines) and posterior (blue lines) densities of π, the variable used for estimating age-specific yearly survival probabilities. In river Pirita the probabilities are concentrated around 0.2 for newborn parrs, with relatively small uncertainty. Yearly survival probabilities rise as age increases and are increasingly more uncertain after age 3. In river Ingarskila the probability estimates are more uncertain for age 0 parrs and the slope of the logistic regression lines seem to be steeper than in River Pirita. For fish aged 6 and higher the probability estimates seem to approach 1.
56
Figure 30.Prior and posterior distributions for age-specific yearly survival probability estimate (π) in the absence of smoltification in both study rivers.
Figures 31 and 32 show the resulting combinations of π and ρ parameters, used to calculate the actual estimates for survival and smoltification probabilities in formulas 5 and 7. In river Ingarskila the highest probability of smoltification seems to be concentrated around age 3, and around age 4 in River Pirita. There is visibly more uncertainty remaining in the posterior distribution of River Ingarskila. Both posteriors resemble the posterior distributions ρ.
57
Figure 31. Posterior distributions for age-specific probability of smoltification (ρπ1/4) in study rivers.
The age-specific survival probability P(survival) is the joint probability of survival and not-smoltifying. In other words it can be interpreted also as the probability of staying in the river another year. In River Pirita there is very little uncertainty for ages 1 to 3 and a very clear drop at age 3, after which probability increases. In River Ingarskila the posterior distribution of P(survival) is more leveled and contains more uncertainty at all ages. Posterior distributions in both rivers resemble an inverted version of P(smolt).
58
Figure 32. Posterior distributions for age-specific probability of survival ((1-ρ)π1/4) in study rivers.
Figures 33 – 35 represent the prior (red lines) and posterior (blue lines) densities of parameter µ, σ and τ in study rivers. Parameter µ, which determines the location of the probability peak in
parameter ρ, seems to have stabilized around age 3 and 3.5 in River Pirita and around age 4 and 4.5 in River Ingarskila. The posterior distribution is much narrower in River Pirita.
59
Figure 33. Prior (red line) and posterior (blue line) distributions for parameter µ study rivers.
Parameter σ determines the range of ages, where individual trout can become migratory. In River Pirita the posterior distribution σ is quite narrow and centered approximately around 1.5. In River Ingarskila the distribution is significantly wider, indicating much higher uncertainty. The peak of the posterior distribution is located near 3.
60
Figure 34. Prior (red line) and posterior (blue line) distributions for parameter σ study rivers.
Parameter τ controls maximum height of parameter ρ and on a population scale it can be interpreted as the proportion of population that migrates out to sea, in the absence of mortality. In River Pirita posterior distribution of τ is located close to 1 and centered around 0.95. In River Ingarskila the posterior distribution is much wider and contains more uncertainty. Posterior mean value is close to 0.45.
61
Figure 35 Prior (red line) and posterior (blue line) distributions for parameter τ study rivers.
Figures 36 and 37 show posterior (blue lines) and prior (dashed red lines) of parameter απ and βπ in study rivers respectively. Parameter απ determines the intercept term for the logistic regression curve used to calculate the values for yearly probability estimate for survival (π). Parameter βπ determines the slope of the curve. The posterior distribution απ is centered approximately around -1 in River Pirita and around 0 in River Ingarskila. Posterior distribution of απ contains more
uncertainty in River Ingarskila, than in River Pirita. Posterior distributions of parameter βπ contain similar amounts uncertainty in both rivers. River Pirita’s posterior is slightly higher in density and its mean value is smaller than the posterior mean in River Ingarskila.
62
Figure 36. Prior (red line) and posterior (blue line) distributions for parameter απ study rivers.
63
Figure 37. Prior (red line) and posterior (blue line) distributions for parameter βπ study rivers
Figures 38 and 39 show the posterior densities of the catchability parameters in both observation models. Catchability in smolt trapping experiments (figure 38) seems to vary between 0 and 0.5 in both rivers throughout the time series.
Electrofishing catchability estimates vary between fish ages in both rivers (figure 39). In River Pirita there is a clear notch in catchability in age 6.This is also seen in River Ingarskila at age 2 and 8. Total range of catchability estimates in both rivers is between almost 1 and just under 0.1.