Middle Columbia River Fish Indexing Analysis 2013

Occupancy and Species Richness

Occupancy which is the probability that a particular species was present at a site was estimated from the temporal replication of detection data (Kery, 2010; Kery and Schaub, 2011, pp. 238-242 and 414-418), i.e., each site was surveyed multiple times within a season. A species was considered to have been detected if one or more individuals of the species were caught or counted. Its important to note that the model estimates the probability that the species was present at a given (or typical) site in a given (or typical) year as opposed to the probability that the species was present in the entire study area. The estimated occupancies for multiple species were summed to give the expected species richnesses.

Key assumptions of the occupancy model include:

• Occupancy (probability of presence) varies with discharge regime and season.
• Occupancy varies randomly with site, year, and the interaction between site and year.
• Sites are closed, i.e., the species is present or absent at a site for all the sessions in a particular season of a year.
• Observed presence is described by a bernoulli distribution.

Count and Species Evenness

The count data were analysed using an overdispersed Poisson model (Kery, 2010; Kery and Schaub, 2011, pp. 168-170,180 and 55-56). Unlike Kery (2010) and Kery and Schaub (2011), which used a log-normal distribution to account for the extra-Poisson variation, the current model used a gamma distribution with identical shape and scale parameters because it has a mean of 1 and therefore no overall effect on the expected count. The model does not distinguish between the abundance and observer efficiency, i.e., it estimates the count which is the product of the two. As such it is necessary to assume that changes in observer efficiency are negligible in order to interpret the estimates as relative abundance. The shannon index of evenness (\(E\)) was calculated using the following formula where \(S\) is the number of species and \(p_i\) is the proportion of the total count belonging to the ith species.

\[ E = \frac{-\sum p_i \log(p_i)}{ln(S)}\]

Key assumptions of the count model include:

• Count density (count/km) varies with discharge regime, season and river kilometre.
• Count density (count/km) varies randomly with site, year, and the interaction between site and year.
• The relationship between count density and river kilometre (distribution) varies with discharge regime and season.
• The relationship between count density and river kilometre (distribution) varies randomly with year.
• Expected counts are the product of the count density (count/km) and the length of bank sampled.
• Sites are closed, i.e., the predicted count at a site is constant for all the sessions in a particular season of a year.
• Observed counts are described by a Poisson-gamma distribution.

Catch

The catch data were analysed using the same overdispersed Poisson model as the count data to provide estimates of relative abundance.

Site Fidelity

The extent to which sites are closed, i.e., fish remain at the same site between sessions, was evaluated from a binomial “t-test” (Kery, 2010, pp. 211-213). The “t-test” estimated the probability that intra-annual recaptures were caught at a different site as previously encountered.

Key assumptions of the site fidelity model include:

• Site fidelity varies with season.
• Observed site fidelity is described by a bernoulli distribution.

Abundance

The catch data were also analysed using a capture-recapture-based binomial mixture model (Kery, 2010; Kery and Schaub, 2011, pp. 253-257 and 134-136, 384-388) to provide estimates of capture efficiency and absolute abundance. To maximize the number of recaptures the model grouped all the sites into a supersite for the purposes of estimating the number of marked fish but analysed the total captures at the site level.

Key assumptions of the abundance model include:

• Lineal density (fish/km) varies with discharge regime, season and river km.
• Lineal density (fish/km) varies randomly with site, year and the interaction between site and year.
• The relationship between density and river kilometre (distribution) varies with discharge regime and season.
• The relationship between density and river kilometre (distribution) varies randomly with year.
• Efficiency (probability of capture) varies randomly by session within season and year.
• Marked and unmarked fish have the same probability of capture.
• There is no tag loss, mortality or misidentification of fish.
• Sites are closed.
• The number of fish captured are described by binomial distributions.

Capture Efficiency

In order to estimate the capture efficiency independent of abundance a recapture-based binomial model (Kery, 2010; Kery and Schaub, 2011, pp. 253-257 and 134-136,384-388) was fitted to just the marked fish. To maximize the number of recaptures the model grouped all the sites into a supersite.

Key assumptions of the efficiency model include:

• Efficiency (probability of capture) varies randomly by session within season and year.
• There is no tag loss, mortality or misidentification of fish.
• The supersite is closed.
• The number of marked fish caught is described by a binomial distribution.

Growth

Annual growth was estimated from the inter-annual recaptures using the Fabens method (Fabens, 1965) for estimating the von Bertalanffy growth curve (Bertalanffy, 1938).

Key assumptions of the growth model include:

• The growth coefficient varies with discharge regime.
• The growth coefficient varies randomly with year.
• Observed growth (change in length) is normally distributed.

Condition

Condition was estimated via an analysis of weight-length relations (He et al. 2008).

Key assumptions of the condition model include:

• Weight varies with length, discharge regime and season.
• Weight varies randomly with site, year and the interaction between site and year.
• Weight is log-normally distributed.

Length

Mean length was estimated from the measured lengths.

Key assumptions of the length model include:

• Length varies with discharge regime and season.
• Length varies randomly with site, year and the interaction between site and year.
• Length is log-normally distributed.

Biomass

The biomass was calculated from the posterior distributions for the Length, Condition and Abundance analyses.

Multivariate Analyses

In order to examine the relationships between environmental variables and the fish indexing metrics multivariate analyses were performed. More specifically the trimonthly mean discharge and elevation and the trimonthly mean absolute hourly discharge change were analysed to get their five primary eigenvectors. Next the correlations between the trimonthly environmental time series and the eigenvectors was quantified using a Bayesian model. The same model was also used to quantify the correlations between the eigenvectors and fish indexing time series related to the management hypotheses. Significant relationships were indicated using time series/eigenvector connectivity plots were nodes are connected if the relationship is significant and positive relationships are in black and negative ones in red.

Model Code

The JAGS model code, which uses a series of naming conventions, is presented below.

model {

  bOccupancy ~ dnorm(0, 5^-2)

  bOccupancySeason[1] <- 0
  for(i in 2:nSeason) {    
    bOccupancySeason[i] ~ dnorm(0, 5^-2)
  }

  bOccupancyRegime[1] <- 0
  for(i in 2:nRegime) {    
    bOccupancyRegime[i] ~ dnorm(0, 5^-2)
  }

  sOccupancyYear ~ dunif(0, 5)
  for (yr in 1:nYear) {
    bOccupancyYear[yr] ~ dnorm(0, sOccupancyYear^-2)
  }

  sOccupancySite ~ dunif(0, 5)
  sOccupancySiteYear ~ dunif(0, 5)
  for (st in 1:nSite) {
    bOccupancySite[st] ~ dnorm(0, sOccupancySite^-2)
    for (yr in 1:nYear) {
      bOccupancySiteYear[st, yr] ~ dnorm(0, sOccupancySiteYear^-2)
    }
  }

  for (i in 1:length(Year)) {

    logit(eOccupancy[i]) <- bOccupancy
      + bOccupancyRegime[Regime[i]] + bOccupancySeason[Season[i]] 
      + bOccupancySite[Site[i]] + bOccupancyYear[Year[i]] 
      + bOccupancySiteYear[Site[i],Year[i]]

    eObserved[i] <- eOccupancy[i]

    Observed[i] ~ dbern(eObserved[i])
  }
}

model {

  bDensity ~ dnorm(0, 5^-2)
  bDistribution ~ dnorm(0, 5^-2)

  bDensityRegime[1] <- 0
  bDistributionRegime[1] <- 0 
  for(i in 2:nRegime) {    
    bDensityRegime[i] ~ dnorm(0, 5^-2)
    bDistributionRegime[i] ~ dnorm(0, 5^-2)
  }

  bDensitySeason[1] <- 0
  bDistributionSeason[1] <- 0
  for(i in 2:nSeason) {    
    bDensitySeason[i] ~ dnorm(0, 5^-2)
    bDistributionSeason[i] ~ dnorm(0, 5^-2)
  }

  sDensityYear ~ dunif(0, 2)
  sDistributionYear ~ dunif(0, 2)
  for (i in 1:nYear) {
    bDensityYear[i] ~ dnorm(0, sDensityYear^-2)
    bDistributionYear[i] ~ dnorm(0, sDistributionYear^-2)
  }

  sDensitySite ~ dunif(0, 5)
  sDensitySiteYear ~ dunif(0, 2)
  for (i in 1:nSite) {
    bDensitySite[i] ~ dnorm(0, sDensitySite^-2)
    for (j in 1:nYear) {
      bDensitySiteYear[i, j] ~ dnorm(0, sDensitySiteYear^-2)
    } 
  }

  bDispersion ~ dgamma(0.1, 0.1)

  for (i in 1:length(Year)) {

      eDistribution[i] <- bDistribution
      + bDistributionRegime[Regime[i]] + bDistributionSeason[Season[i]]
      + bDistributionYear[Year[i]]

    log(eDensity[i]) <- bDensity 
      + eDistribution[i] * RiverKm[i] 
      + bDensityRegime[Regime[i]] + bDensitySeason[Season[i]] 
      + bDensitySite[Site[i]] + bDensityYear[Year[i]] 
      + bDensitySiteYear[Site[i],Year[i]]

    eCount[i] <- eDensity[i] * SiteLength[i] * ProportionSampled[i] 

    eDispersion[i] ~ dgamma(bDispersion, bDispersion)

    Count[i] ~ dpois(eCount[i] * eDispersion[i])  
  }
  tAbundance <- bDensityRegime[2]
  tDistribution <- bDistributionRegime[2]
}

model {

  bDensity ~ dnorm(0, 5^-2)
  bDistribution ~ dnorm(0, 5^-2)

  bDensityRegime[1] <- 0
  bDistributionRegime[1] <- 0 
  for(i in 2:nRegime) {    
    bDensityRegime[i] ~ dnorm(0, 5^-2)
    bDistributionRegime[i] ~ dnorm(0, 5^-2)
  }

  bDensitySeason[1] <- 0
  bDistributionSeason[1] <- 0
  for(i in 2:nSeason) {    
    bDensitySeason[i] ~ dnorm(0, 5^-2)
    bDistributionSeason[i] ~ dnorm(0, 5^-2)
  }

  sDensityYear ~ dunif(0, 2)
  sDistributionYear ~ dunif(0, 2)
  for (i in 1:nYear) {
    bDensityYear[i] ~ dnorm(0, sDensityYear^-2)
    bDistributionYear[i] ~ dnorm(0, sDistributionYear^-2)
  }

  sDensitySite ~ dunif(0, 5)
  sDensitySiteYear ~ dunif(0, 2)
  for (i in 1:nSite) {
    bDensitySite[i] ~ dnorm(0, sDensitySite^-2)
    for (j in 1:nYear) {
      bDensitySiteYear[i, j] ~ dnorm(0, sDensitySiteYear^-2)
    } 
  }

  bDispersion ~ dgamma(0.1, 0.1)

  for (i in 1:length(Year)) {
    eDistribution[i] <- bDistribution
      + bDistributionRegime[Regime[i]] 
      + bDistributionSeason[Season[i]]
      + bDistributionYear[Year[i]]

    log(eDensity[i]) <- bDensity 
      + eDistribution[i] * RiverKm[i] 
      + bDensityRegime[Regime[i]] + bDensitySeason[Season[i]] 
      + bDensitySite[Site[i]] + bDensityYear[Year[i]] 
      + bDensitySiteYear[Site[i],Year[i]]

    eCatch[i] <- eDensity[i] * SiteLength[i] * ProportionSampled[i] 

    eDispersion[i] ~ dgamma(bDispersion, bDispersion)

    Catch[i] ~ dpois(eCatch[i] * eDispersion[i])  
  }
  tAbundance <- bDensityRegime[2]
  tDistribution <- bDistributionRegime[2]
}

model {
  bMoved ~ dnorm(0, 5^-2)

  bMovedSeason[1] <- 0
  for(i in 2:nSeason) {    
    bMovedSeason[i] ~ dnorm(0, 5^-2)
  }

  for (i in 1:length(Season)) {

    logit(eMoved[i]) <- bMoved + bMovedSeason[Season[i]]

    Moved[i] ~ dbern(eMoved[i])
  }
}

model {

  bEfficiency ~ dnorm(0, 5^-2)
  bDensity ~ dnorm(0, 5^-2)
  bDistribution ~ dnorm(0, 5^-2)

  bDensityRegime[1] <- 0
  bDistributionRegime[1] <- 0
  for(i in 2:nRegime) {    
    bDensityRegime[i] ~ dnorm(0, 5^-2)
    bDistributionRegime[i] ~ dnorm(0, 5^-2)
  }

  bEfficiencySeason[1] <- 0
  bDensitySeason[1] <- 0
  bDistributionSeason[1] <- 0
  for(i in 2:nSeason) {
    bEfficiencySeason[i] ~ dnorm(0, 5^-2)
    bDensitySeason[i] ~ dnorm(0, 5^-2)
    bDistributionSeason[i] ~ dnorm(0, 5^-2)
  }

  sDensityYear ~ dunif(0, 2)
  sDistributionYear ~ dunif(0, 2)
  for (i in 1:nYear) {
    bDensityYear[i] ~ dnorm(0, sDensityYear^-2)
    bDistributionYear[i] ~ dnorm(0, sDistributionYear^-2)
  }

  sDensitySite ~ dunif(0, 5)
  sDensitySiteYear ~ dunif(0, 2)
  for (i in 1:nSite) {
    bDensitySite[i] ~ dnorm(0, sDensitySite^-2)
    for (j in 1:nYear) {
      bDensitySiteYear[i, j] ~ dnorm(0, sDensitySiteYear^-2)
    } 
  }

  sEfficiencySessionSeasonYear ~ dunif(0, 5)
  for (i in 1:nSession) {
    for (j in 1:nSeason) {
      for (k in 1:nYear) {
        bEfficiencySessionSeasonYear[i, j, k] ~ dnorm(0, sEfficiencySessionSeasonYear^-2)
      }
    }
  }

  for(i in 1:length(EffIndex)) {

    logit(eEff[i]) <- bEfficiency 
        + bEfficiencySeason[Season[EffIndex[i]]]
        + bEfficiencySessionSeasonYear[Session[EffIndex[i]],
                                       Season[EffIndex[i]], 
                                       Year[EffIndex[i]]]

    Marked[EffIndex[i]] ~ dbin(eEff[i], Tagged[EffIndex[i]])
  }

  for (i in 1:length(Year)) {

    logit(eEfficiency[i]) <- bEfficiency 
        + bEfficiencySeason[Season[i]]
        + bEfficiencySessionSeasonYear[Session[i], Season[i], Year[i]]

    eDistribution[i] <- bDistribution
      + bDistributionRegime[Regime[i]] 
      + bDistributionSeason[Season[i]]
      + bDistributionYear[Year[i]]

    log(eDensity[i]) <- bDensity 
      + eDistribution[i] * RiverKm[i] 
      + bDensityRegime[Regime[i]] 
      + bDensitySeason[Season[i]] 
      + bDensitySite[Site[i]] 
      + bDensityYear[Year[i]] 
      + bDensitySiteYear[Site[i], Year[i]]

    eSamplingEfficiency[i] <- min(eEfficiency[i] * ProportionSampled[i], 0.9)

    eAbundance[i] <- max(round(eDensity[i] * SiteLength[i]), MinAbundance[i])

    Catch[i] ~ dbin(eSamplingEfficiency[i], eAbundance[i])
  }
  tAbundance <- bDensityRegime[2]
  tDistribution <- bDistributionRegime[2]
}

model {

  bEfficiency ~ dnorm(0, 5^-2)

  bEfficiencySeason[1] <- 0
  for (i in 2:nSeason) {
    bEfficiencySeason[i] ~ dnorm(0, 5^-2)
  }

  sEfficiencySessionSeasonYear ~ dunif(0, 5)
  for (i in 1:nSession) {
    for (j in 1:nSeason) {
      for (k in 1:nYear) {
        bEfficiencySessionSeasonYear[i, j, k] ~ dnorm(0, sEfficiencySessionSeasonYear^-2)
      }
    }
  }

  for(i in 1:length(Year)) {
    logit(eEfficiency[i]) <- bEfficiency 
        + bEfficiencySeason[Season[i]]
        + bEfficiencySessionSeasonYear[Session[i], Season[i], Year[i]]

    Marked[i] ~ dbin(eEfficiency[i], Tagged[i])
  }
}

model {

  bK ~ dnorm (0, 5^-2)

  bKRegime[1] <- 0
  for(i in 2:nThreshold) {
    bKRegime[i] ~ dunif(-100, 100)
  }

  sKYear ~ dunif (0, 5)
  for (i in 1:nYear) {
    bKYear[i] ~ dnorm(0, sKYear^-2)
    log(eK[i]) <- bK + bKRegime[step(i - Threshold) + 1] + bKYear[i]
  }

  bLinf ~ dunif(100, 1000)
  sGrowth ~ dunif(0, 100)

  for (i in 1:length(Year)) {
    eGrowth[i] <- (bLinf - LengthAtRelease[i])
                  * (1 - exp(-sum(eK[Year[i]:(Year[i] + Years[i] - 1)])))

    Growth[i] ~ dnorm(eGrowth[i], sGrowth^-2)
  }
  tGrowth <-  bKRegime[2]
} 

model {

  bWeight ~ dnorm(5, 5^-2)
  bWeightLength ~ dnorm(3, 5^-2)

  bWeightRegime[1] <- 0
  for(i in 2:nRegime) {    
    bWeightRegime[i] ~ dnorm(0, 5^-2)
  }

  bWeightSeason[1] <- 0
  for(i in 2:nSeason) {    
    bWeightSeason[i] ~ dnorm(0, 5^-2)
  }

  sWeightYear ~ dunif(0, 5)
  for(yr in 1:nYear) {
    bWeightYear[yr] ~ dnorm(0, sWeightYear^-2)
  }  

  sWeightSite ~ dunif(0, 5)
  sWeightSiteYear ~ dunif(0, 5)
  for(st in 1:nSite) {
    bWeightSite[st] ~ dnorm(0, sWeightSite^-2)
    for(yr in 1:nYear) {
      bWeightSiteYear[st, yr] ~ dnorm(0, sWeightSiteYear^-2)
    }
  } 

  sWeight ~ dunif(0, 5)

  for(i in 1:length(Year)) {

    eWeightIntercept[i] <- bWeight
        + bWeightRegime[Regime[i]] 
        + bWeightSeason[Season[i]]
        + bWeightYear[Year[i]] + bWeightSite[Site[i]] 
        + bWeightSiteYear[Site[i],Year[i]]

    eWeightSlope[i] <- bWeightLength

    log(eWeight[i]) <- eWeightIntercept[i] + eWeightSlope[i] * Length[i]

    Weight[i] ~ dlnorm(log(eWeight[i]), sWeight^-2)
  }
  tCondition <- bWeightRegime[2]
}

model {

  bLength ~ dnorm(5, 5^-2)

  bLengthRegime[1] <- 0
  for(i in 2:nRegime) {    
    bLengthRegime[i] ~ dnorm(0, 5^-2)
  }

  bLengthSeason[1] <- 0
  for(i in 2:nSeason) {    
    bLengthSeason[i] ~ dnorm(0, 5^-2)
  }

  sLengthYear ~ dunif(0, 5)
  for(yr in 1:nYear) {
    bLengthYear[yr] ~ dnorm(0, sLengthYear^-2)
  }  

  sLengthSite ~ dunif(0, 5)
  sLengthSiteYear ~ dunif(0, 5)
  for(st in 1:nSite) {
    bLengthSite[st] ~ dnorm(0, sLengthSite^-2)
    for(yr in 1:nYear) {
      bLengthSiteYear[st, yr] ~ dnorm(0, sLengthSiteYear^-2)
    }
  } 

  sLength ~ dunif(0, 5)

  for(i in 1:length(Year)) {

    log(eLength[i]) <- bLength
      + bLengthRegime[Regime[i]] 
      + bLengthSeason[Season[i]] 
      + bLengthYear[Year[i]] + bLengthSite[Site[i]] 
      + bLengthSiteYear[Site[i],Year[i]]

    Length[i] ~ dlnorm(log(eLength[i]), sLength^-2)
  }
}

model {

  sValue ~ dunif(0, 2)

  for(k in 1:nEigen) {
    sWeight[k] ~ dunif(0, 5)
  }

  for (i in 1:nSeries) {
    for (k in 1:nEigen) {
      Weight[i, k] ~ dnorm(0, sWeight[k]^-2)
    }
    for (t in 1:nYear) {
      Fit[i, t] <- inprod(Weight[i,], Eigen[,t])
      Value[i, t] ~ dnorm(Fit[i, t], sValue^-2)
    }
  }
} 

Model Parameters

The posterior distributions for the fixed (Kery and Schaub 2011 p. 75) parameters in each model are summarised below.

Figures

Species Richness

Species Evenness

Occupancy - Burbot

Occupancy - Kokanee

Occupancy - Lake Whitefish

Occupancy - Northern Pikeminnow

Occupancy - Rainbow Trout

Occupancy - Redside Shiner

Occupancy - Sculpin

Count - Bull Trout

Count - Burbot

Count - Mountain Whitefish

Count - Northern Pikeminnow

Count - Rainbow Trout

Count - Suckers

Catch - Adult BT

Catch - Adult CSU

Catch - Adult MW

Catch - Adult RB

Catch - Juvenile BT

Catch - Juvenile MW

Site Fidelity - Adult BT

Site Fidelity - Adult CSU

Site Fidelity - Adult MW

Site Fidelity - Adult RB

Site Fidelity - Juvenile BT

Site Fidelity - Juvenile MW

Abundance - Adult BT

Abundance - Adult CSU

Abundance - Adult MW

Abundance - Juvenile MW

Capture Efficiency - Adult BT

Capture Efficiency - Adult CSU

Capture Efficiency - Adult MW

Capture Efficiency - Adult RB

Capture Efficiency - Juvenile BT

Capture Efficiency - Juvenile MW

Growth - Bull Trout

Growth - Mountain Whitefish

Condition - Adult BT

Condition - Adult MW

Condition - Adult RB

Condition - Juvenile BT

Condition - Juvenile MW

Biomass - Adult BT

Biomass - Adult MW

Biomass - Juvenile MW

Multivariate Analysis - Eigen Vector

Multivariate Analysis - Environmental

Multivariate Analysis - Indexing