Re-Analysis of 2018 Lookout Point Reservoir Data

Field data from Lookout Point Reservoir collected in 2018 were reanalyzed to determine whether mortality could be statistically attributed to two different sources. These data are described in detail in Kock and others (2019b). The two sources of mortality considered were (1) an index of sampling occasion (time) as a surrogate measure of reduced predation risk with increasing fish size and (2) adult copepod incidence rate in sampled fish. We did not attempt to model field data collected in 2017 because information on copepod infection was not collected. Adult copepod infection rate was the proportion of captured juvenile Chinook salmon with visible adult, female copepods attached averaged across release groups. The average across groups was used to reduce the number of parameters required in the model and thus increase power to detect an effect. Infection rate increased through time in 2018 from relatively low levels in June and July to very high levels (exceeding 90 percent) in September (fig. 2).

We hypothesized that the two sources of mortality should have countervailing effects. That is, absent any effect of copepod infection, monthly survival probability should increase over time because of decreased predation as juvenile Chinook salmon grew to larger sizes. Survival estimates from the 2017 field study supported this hypothesis as Kock and others (2019a) found support in the data for a positive trend in survival with increasing time (fig. 1A). Contrastingly, absent any effect of size-dependent predation, we hypothesized that survival should decrease with increasing copepod infection rates. Emerging evidence has demonstrated the potential for copepod infections to substantially decrease survival of juvenile Chinook salmon (Beeman and others, 2015; Herron and others, 2018, Neal and others, 2021). In particular, Neal and others (2021) demonstrated through laboratory-controlled infection of juvenile Chinook salmon that fish exposed to greater levels of copepod density and higher water temperatures experienced elevated mortality. The declining trend in monthly survival estimates from 2018 (Kock and others, 2019b; fig. 1B) coincident with rising copepod infection rates suggest that copepod induced mortality is a plausible explanation for the overserved trend.

We used the PBT N-mixture model previously described in Kock and others (2019a, b) to assess the support in the data for countervailing effects from these two covariates. We briefly review the structure of that model here. The model consists of two related equations:

The likelihood of the data $y_{ijk}$ given the survival and collection probability parameters $S_k$ and p_jk is the product of equations (1) and (2) over all PBT family groups $i\in \left\{1,\dots ,R\right\},$ and all primary sampling occasions $k\in \left\{1,\dots ,K\right\},$ where the relation between unconditional per-sample capture probability p_jk and ${\pi }_{ijk}$ is governed by the recursive equation

This model assumes the following:

The capture probability parameter, p, was set to vary among sampling gear used (electrofishing, gill net, and box trap), and whether fish were nearshore or offshore. We assumed that fry were nearshore during the first two primary sampling occasions and offshore thereafter (Monzyk and others, 2015a). Capture probability was modeled as:

For survival, we fit three models to the data. The first model corresponded to a “predation-only” model which assumed that monthly survival increased or decreased at a fixed rate over time. Here, time is treated as surrogate for fish growing larger and changing behavior to better evade predators. The second model corresponded to “copepod-only” model which assumed that monthly survival depended only on the proportion of captured fish with visible copepod infection. The third model included both effects and is represented as:

We used standard normal prior distributions for all slope and intercept parameters. Models were compared based on the whether 90 percent credible intervals (CI) for parameter estimates of the covariate effects contained only positive or negative values. If the CI contained both positive and negative values, then we concluded no effect of the covariate could be discerned. We also estimated and compared the six period-specific survival estimates for the three models: (1) mid-April to mid-May (S₁), (2) mid-May to mid-June (S₂), (3) mid-June to mid-July (S₃), (4) mid-July to mid-August (S₄), (5) mid-August to mid-September (S₅), and (6) mid-September to mid-October (S₆).

Multi-Year Data Simulation with Subsampling

We conducted a simulation study to assess the potential for a study design to separately attribute mortality to predation and copepod infection. There were two major motivating factors to the design of this simulation study. First, we designed the study to use an assessment for copepod infection based on a subsample of captured fish. Second, we designed the study to be conducted over the course of multiple field seasons with interannual variation in copepod infection rates. In the simulation study we varied the subsampling rate, the number of years and the degree of interannual variation.

We simulated a subsampling scheme out of concern that using adult copepod presence as a measure of infection rates may not be an accurate measure of copepod infection dynamics. Emerging evidence suggests that the presence of adult copepods may be a lagging indicator of copepod infection rates. For example, Neal and others (2021) demonstrated that substantial mortality occurred before adult, female copepods were visible in juvenile Chinook salmon held in a controlled laboratory setting and attributed this in part to damage to gills from feeding by microscopic juvenile copepods. Sanders and others (2021) suggested that “prevalence surveys that rely on counting grossly visible adult female copepods only greatly underestimate the impact on juvenile Chinook salmon survival.” Both Neal and others (2021) and Sanders and others (2021) characterized gill damage in juvenile Chinook salmon through histopathological examination of excised gills. However, as Sanders and others (2021) notes, histopathology is a low-throughput and labor-intensive method which requires special processing of each sampled fish. Because of the nature of this method, we considered that this method could only realistically be applied to a small subsample of fish collected in the field. Thus, in our simulation we based copepod infection prevalence on only a subset of fish collected.

To better evaluate the feasibility of identifying two competing sources of mortality, we structured our simulation study to mimic a field effort conducted over the course of 2 or more years. A major challenge to estimating covariate relations in observational studies is the phenomenon of multicollinearity, defined as a correlation between two or more predictor variables. Data from the 2018 field study demonstrated such collinearity between copepod infection and time as adult copepod prevalence increased monotonically over the course of the study. The consequence of multicollinearity is that a statistical model will produce unreliable estimates of the coefficients for correlated variables. Overall predictions from a model under multicollinearity will tend to be accurate (in this case, time-interval specific survival estimates), but estimates for separate covariate effects (that is, time and copepod infection) may not be. For example, estimates for the coefficients will tend to have large standard errors. The only remedy for multicollinearity is to design a study in such a way that the correlation between two variables can be broken. In observational studies, where it is usually not possible to experimentally manipulate covariates, one way to break this collinearity is to collect data over a broad range of conditions such that the correlation between variables is lessened. In the case of the PBT N-mixture model for Lookout Point this means sampling over 2 or more years.

The extent to which copepod prevalence varies from year to year in Lookout Point is unknown. There is only limited evidence to suggest that the relation between copepods and time may differ in different years. Data on adult copepod prevalence were not collected by Kock and others (2019a) in 2017 but was collected in 2018 (Kock and others, 2019b). Additionally, Sanders and others (2021) used similar methods to collect juvenile Chinook salmon in 2020 to evaluate copepod prevalence on sampled fish over a similar time frame. Kock and others (2019b) collected adult copepod prevalence data on age-0 hatchery and natural-origin Chinook salmon from June to October in 2018 and Sanders and others (2021) between July and December 2020. Copepod prevalence for natural-origin Chinook salmon collected by Kock and others (2019a) was less than 15 percent in mid-July and greater than 50 percent in the remaining months, whereas Sanders and others (2021) observed copepod prevalence of approximately 20 percent in mid-July but between 40 and 50 percent in the remaining months. Whereas hatchery-origin and age-1 Chinook salmon had higher copepod prevalence in Kock and others (2019b), a comparison cannot be made to Sanders and others (2021) for these groups. Thus, we have evidence for 2 years of qualitatively similar copepod prevalence dynamics, but with some slight differences in timing and absolute levels.

Given the unknowns of copepod infection dynamics, we structured our simulation to evaluate multiple scenarios of interannual variation in copepod infection rates. We generated three families of copepod infection curves with respect to time: “Low,” “Moderate,” and “High.” These curves were drawn from the Richard distribution (Richards, 1959), a generalized logistic function obtained to model the growth of organisms over time and which has been used in epidemiological models (Lee, 2020). For each year in a simulation a new growth curve was drawn from the family specified by the simulation scenario to determine the probability that an uninfected fish at sampling occasion t will be infected before sampling occasion t + 1. The curve families are depicted in figure 3 and were obtained by drawing parameters of the Richard curve from different sets of bounded uniform distributions. In contrast to existing field datasets, we envisioned the infection rates detailed here to produce the proportion of the juvenile Chinook salmon population infected by juvenile copepods, rather than the proportion showing grossly visible adult copepods.

We produced a series of simulated datasets with common characteristics varying only the number of years in the study, the type of copepod infection curves in each year, and the number of fish subsampled during each sampling occasion. For simplicity, all fish in the simulation were released as a single group at the start of the simulated study. The common characteristics included the number of distinct PBT family groups $\left(R = 70\right)$, the total number of juveniles released $\left(N = 175,000\right)$, the number of primary $\left(K=6\right)$ and secondary $\left(J= 4\right)$ sampling occasion, the capture probability during each secondary sampling occasion $\left(p= 0.00075\right)$, and the parameters determing survival. The survival parameters were logit-scaled and included an intercept (${\beta }_0=1.1)$, the slope of the time effect (${\beta }_1=0.22)$, and the slope of the copepod effect (${\beta }_2=-2.49)$. These correspond to a baseline survival of 75 percent at $t = 0$ and ${\pi }_t=0$ (where ${\pi }_t$ is the prevalence of copepod infection in the population), a survival of 95 percent at $t = 5$ and ${\pi }_t=0$, and a survival of 20 percent at $t = 0$ and ${\pi }_t=1$ (that is, 100 percent of population infected). Figure 4 depicts the combined effect of both covariates as well as the marginal effect of each (time in blue, copepod prevalence in red) for an example copepod prevalence curve. We simulated 96 datasets for analysis using a multi-year version of the PBT N-mixture model. These datasets consisted of three replicates of each specific combination of number of years, year-types, and subsampling rates. We considered two subsampling rates: 10 fish per sample period and 30 fish per sample period. We simulated datasets covering 2 and 3 years and considered each potential combination of year-types. For example, for the 2-year models the year types consisted of “Low-Low,” “Low-Moderate,” “Low-High,” “Moderate-Moderate,” “Moderate-High,” and “High-High.”

Each simulated dataset was constructed by iteratively stepping through the data generating process implied by the PBT N-mixture model. The simulations began with a release of $N$ fish divided among $R$ release groups at time $t=0$. Because we assume that copepod infection rate is zero for the first time step, these fish survive to time$t = 1$ with probability equal to logit(${\beta }_0$) = 0.75. Prior to sampling, surviving uninfected fish (all fish at time $t=0$) are “infected” with copepods with the probability of infection being determined by the copepod infection curve drawn for that simulation. Surviving fish were then “captured” through removal sampling over $J$ sampling occasions with the probability of being captured equal to $p$. Captured fish were then subject to random sampling to be assigned to the subsample group for “histopathology” to determine infection status. For the remaining time periods, the probability of survival from time $t$ to $t+1$ was calculated based on the population level infection rate at time $t$ (${\pi }_t$) according to the relation logit$\left(S_t\right)={\beta }_0+{\beta }_{1}t+ {\beta }_2{\pi }_t$ and the process of “infection,” capture, and subsampling was repeated. If fewer fish were captured than the assigned subsampling rate, then all captured fish were used for simulated histopathological examination.

Each simulated dataset was analyzed using a modified version of the PBT N-mixture model described above. The model was modified to apply to multiple years by looping over likelihoods for all years. Capture probability was also simplified because we ignored the complexities of multiple gear-types and off-shore/on-shore effects for the purposes of simulation. The final modifications accounted for uncertainty in the true copepod infection prevalence because of the subsampling effect and related survival to the latent true infection prevalence. That is for each timestep, we assumed $c_t ~ Binomal\left(Z_t, {\pi }_t\right)$, where $Z_t$ is the number of fish subsampled and $c_t$ is the number of these fish with signs of copepod infection. We then estimated survival according to the relation $\text{logit}\left(S_t\right)={\beta }_0+{\beta }_{1}t+ {\beta }_2{\pi }_t$. We evaluated each simulation in terms of Type S and Type M errors (Gelman and Carlin, 2014). Type S errors are errors of the sign of the estimate, if the 90 percent posterior uncertainty interval for ${\beta }_1$ and ${\beta }_2$ contained both positive and negative values or contained only negative values for ${\beta }_1$ and only positive values for ${\beta }_2$ then this would be an instance of Type S error. Type M errors are errors of the magnitude of the estimate, which would be the case where the true parameter value was not in 90 percent posterior uncertainty value. In evaluating the model’s ability to detect an effect of both time (predation) and copepod infection, Type S error takes precedence over Type M error. That is, if credible intervals are of the correct sign, but wrong magnitude (Type S error negative, type M error positive) then we can conclude that time and copepods have a positive and negative effect respectively; even if the absolute value of that effect is incorrect, we can use this as a basis for making management decisions. On the other hand, if the credible interval contains the true parameter value but also contains values of the opposite sign, then we cannot conclude there is any effect and thus cannot improve the information on which we base decisions. Because of the limited number of replications (a consequence of limited computational efficiency), our analysis of Type M and Type S error rates are only qualitative.

The results of the simulation study must be interpreted with caution because it makes several assumptions which may not hold in the field. Among these, we assumed that the survival effects would be constant from year to year, which is unlikely to be the case in nature. We also assumed that there was no interaction between time (a surrogate of fish size) and copepod infection which may not hold in nature. Examples of interactions between the two covariates include the following: (1) if larger fish are more likely to die from copepod infection than smaller fish (given equal infection rates); or (2) if larger fish are more likely to be infected but equally likely to die from infection. The true field implementation is also likely to be more complicated than the relatively simple simulation model obtained here; in particular, capture probability was set to be a constant. If it is necessary to account for additional parameters (for example, capture probabilities from different gear types) this may reduce the ability of the model to detect an effect. There is no empirical evidence for the effect of copepod infection prevalence as defined here (that is, infection damage prior to the appearance of adult copepods) and an effect of lesser magnitude is less likely to be detected. Finally, it is unclear whether such a subsampling approach is feasible in a field study; that is, we do not know if gill tissue recovered in the field will be of sufficient quality for histopathological examination.

Re-Analysis of 2018 Lookout Point Data

Re-analysis of the 2018 Lookout Point data was unable to identify separate effects of both time and copepod adult prevalence on survival. Whereas both the “predation only” and “copepod only” models detected a negative effect of the respective covariate, the model with both covariates estimated a negative effect for time that was larger than the effect of the “predation only” model and a positive effect with high uncertainty for the copepod model (fig. 5A). Additionally, the 90 percent posterior uncertainty model for the copepod effect in the combined model overlapped zero which can be interpreted as a non-significant effect. That is, there is some non-negligible probability that the effect of copepods could be either negative or positive. Survival estimates between the three models were similar, with the posterior interval overlapping for all time periods except the third and last (fig. 5B). These results are consistent with what one would expect under multicollinearity; that is, the data are insufficient to estimate separate contributions between the two covariates.

Multi-Year Data Simulation with Subsampling

The multiyear data simulation suggests that identifying a signal of the two countervailing covariates will be challenging. Of 96 simulated datasets, only 6 had both 90 percent posterior intervals which contained the true values and did not overlap zero. That is, these simulations had neither Type S nor Type M error. Two of these occurred in simulations of 2-year studies (36 total simulations, fig. 6) and four occurred in simulations of 3-year studies (60 total simulations, fig. 7). All 6 occurred in simulations where the subsample rate was set to 30 fish per sample. All but one occurred in simulations where the year-type combinations included more than one year-type and only one year-type class had more than one replication successfully estimate both parameters (“Low-Low-Moderate”). The models were somewhat better at distinguishing the effect of copepod infection rate alone (18 of 96). When considering only Type S error, 8 out of 96 models got the direction but not magnitude of the effect correct for both covariates, whereas 28 out of 96 models got the direction but not the magnitude correct for the copepod effect alone.