Skip to contents

salmonMSE utilizes an age-structured model. The population is tracked by age and year but various dynamics correspond to the salmon life stages as described below.

Definitions

Definition of variable names and the corresponding slots in either the input (SOM) or output (SMSE) objects in salmonMSE.

Name Definition Type Class Slot
NOS\textrm{NOS} Natural origin spawners Natural production SMSE NOS
HOS\textrm{HOS} Hatchery origin spawners Hatchery SMSE HOS
HOSeff\textrm{HOS}_\textrm{eff} Effective number of HOS, spawning output discounted by γ\gamma Hatchery SMSE HOSeff
FryNOS\textrm{Fry}^\textrm{NOS} Fry production by natural origin spawners, assumed to be equal to egg production Natural production SMSE Fry_NOS
FryHOS\textrm{Fry}^\textrm{HOS} Fry production by hatchery origin spawners, assumed to be equal to egg production Hatchery SMSE Fry_HOS
SmoltNOS\textrm{Smolt}^\textrm{NOS} Smolt production by natural origin spawners, density-dependent Natural production SMSE Smolt_NOS
SmoltHOS\textrm{Smolt}^\textrm{HOS} Smolt production by hatchery origin spawners, density-dependent Hatchery SMSE Smolt_HOS
CsmoltC_\textrm{smolt} Carrying capacity of smolts (Beverton-Holt stock-recruit parameter) Natural production + Habitat SOM capacity_smolt
SmaxS_\textrm{max} Spawning output that maximizes smolt production (Ricker stock-recruit parameter) Natural production + Habitat SOM Smax
κ\kappa Maximum recruitment production rate (stock-recruit parameter), units of recruit per spawner Natural production + Habitat SOM kappa
ϕ\phi Unfished egg production rate, units of egg per smolt Natural production SOM phi
rr Maturity at age, i.e., recruitment rate Natural production SOM p_mature
uPTu^\textrm{PT} Pre-terminal fishery harvest rate Harvest SOM u_preterminal
uTu^\textrm{T} Terminal fishery harvest rate Harvest SOM u_terminal
mm Mark rate of hatchery fish Harvest SOM m
δ\delta Mortality from catch and release (proportion) Harvest SOM release_mort
vv Relative vulnerability by age to the fishery Harvest SOM vulPT, vulT
Fec\textrm{Fec} Fecundity of spawners (eggs per female) Natural production SOM fec
Fecbrood\textrm{Fec}^\textrm{brood} Fecundity of broodtake (eggs per female) Natural production SOM fec_brood
pfemalep^\textrm{female} Proportion of female spawners in broodtake and spawners Natural production + Hatchery SOM p_female
SAR\textrm{SAR} Smolt-to-adult recruit survival Natural production + Hatchery - -
MM Instantaneous natural mortality of juvenile corresponding to either the freshwater or marine environment by age class Natural production + Hatchery SOM
MjuvNOS,MjuvHOSMjuv_NOS, Mjuv_HOS
NOB\textrm{NOB} Natural origin broodtake Hatchery SMSE NOB
HOB\textrm{HOB} Hatchery origin broodtake Hatchery SMSE HOB
syearlings_\textrm{yearling} Survival of hatchery eggs to yearling life stage Hatchery SOM s_egg_smolt
ssubyearlings_\textrm{subyearling} Survival of hatchery eggs to subyearling life stage Hatchery SOM s_egg_subyearling
pyearlingp_\textrm{yearling} Proportion of hatchery releases as yearling (vs. subyearling) Hatchery Internal state variable -
sprespawns_\textrm{prespawn} Survival of adult broodtake in hatchery Hatchery SOM s_prespawn
nyearlingn_\textrm{yearling} Target number of hatchery releases as yearlings Hatchery SOM n_yearling
nsubyearlingn_\textrm{subyearling} Target number of hatchery releases as subyearlings Hatchery SOM n_subyearling
pNOBp_\textrm{NOB} Proportion of the total broodtake of natural origin (vs. hatchery origin) Hatchery Internal state variable -
pHOBp_\textrm{HOB} Proportion of the total broodtake of hatchery origin (vs. natural origin) Hatchery Internal state variable -
pmaxNOBp^\textrm{NOB}_\textrm{max} Maximum proportion of the natural origin broodtake from the escapement, i.e., NOB/NOS ratio Hatchery SOM pmax_NOB
NOR\textrm{NOR} Natural origin return Natural production SMSE Return_NOS
HOR\textrm{HOR} Hachery origin return Hatchery SMSE Return_HOS
phatcheryp^\textrm{hatchery} Proportion of hatchery origin escapement to hatchery, available for broodtake Hatchery SOM phatchery
premovalHOSp^\textrm{HOS}_\textrm{removal} Proportion of hatchery origin escapement removed from spawning grounds, not available for broodtake Hatchery SOM premove_HOS
γ\gamma Reduced reproductive success of HOS (relative to NOS) Hatchery SOM gamma
z\bar{z} Mean phenotypic value of cohort in natural and hatchery environments Fitness Internal state variable and SOM zbar_start
θ\theta Optimal phenotypic value for natural and hatchery environments Fitness SOM theta
σ2\sigma^2 Variance of phenotypic traits in population Fitness SOM fitness_variance
Ω\Omega Selection strength: fitness variance relative to phenotypic variance Fitness SOM selection_strength
h2h^2 Heritability of phenotypic traits Fitness SOM heritability
W\bar{W} Population fitness in the natural and hatchery environments Fitness SMSE fitness
i\ell_i Relative fitness loss at the i-th life stage (egg, fry, smolt) Fitness SOM rel_loss
PNI\textrm{PNI} Proportionate natural influence Fitness SMSE PNI
mm Mark rate of hatchery fish (affects fishery retention of hatchery fish relative to natural fish) Harvest SOM m
pNOSp_\textrm{NOS} Proportion of natural origin spawners (vs. effective HOS) Population dynamics Internal state variable -
pHOSeffp_\textrm{HOSeff} Proportion of effective hatchery origin spawners (vs. NOS) Population dynamics SMSE NOS, HOS_effective
pHOScensusp_\textrm{HOScensus} Proportion of hatchery origin spawners (vs. NOS) Population dynamics SMSE NOS, HOS

Natural production

First, we consider natural production in the absence of fitness effects from the hatchery improvement.

Egg production

From the spawners (NOS and HOS) of age aa in year yy, the corresponding egg production of the subsequent generation is calculated as:

EggyNOS=aNOSy,a×pfemale×FecaEggyHOS=aHOSeffy,a×pfemale×Feca\begin{align} \textrm{Egg}^\textrm{NOS}_y &= \sum_a\textrm{NOS}_{y,a} \times p^\textrm{female} \times \textrm{Fec}_a\\ \textrm{Egg}^\textrm{HOS}_y &= \sum_a\textrm{HOS}_{\textrm{eff}y,a} \times p^\textrm{female} \times \textrm{Fec}_a \end{align}

where HOSeff=γ×HOS\textrm{HOS}_{\textrm{eff}} = \gamma \times \textrm{HOS} and the superscript denotes the parentage of the progeny.

Fry

Fry production is assumed to be equal to egg production, i.e., Fryy+1NOS=EggyNOS\textrm{Fry}^\textrm{NOS}_{y+1} = \textrm{Egg}^\textrm{NOS}_y and FryyHOS=EggyHOS\textrm{Fry}^\textrm{HOS}_y = \textrm{Egg}^\textrm{HOS}_y.

Smolts

Survival to the smolt stage is density-dependent. With the Beverton-Holt stock-recruit relationship, the age-1 smolt production is

Smolty+1NOS=α×Fryy+1NOS1+β(Fryy+1NOS+Fryy+1HOS+ny+1sub)Smolty+1HOS=α×Fryy+1HOS1+β(Fryy+1NOS+Fryy+1HOS+ny+1sub)\begin{align} \textrm{Smolt}^\textrm{NOS}_{y+1} &= \frac{\alpha \times \textrm{Fry}^\textrm{NOS}_{y+1}}{1 + \beta(\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1} + n^\textrm{sub}_{y+1})}\\ \textrm{Smolt}^\textrm{HOS}_{y+1} &= \frac{\alpha \times \textrm{Fry}^\textrm{HOS}_{y+1}}{1 + \beta(\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1} + n^\textrm{sub}_{y+1})} \end{align}

where α=κ/ϕ\alpha = \kappa/\phi, β=α/Csmolt\beta = \alpha/{C_\textrm{smolt}}, the unfished egg per smolt ϕ=a(i=1a1exp(MiNOS)(1ri))×ra×pfemale×Feca\phi = \sum_a\left(\prod_{i=1}^{a-1}\exp(-M^\textrm{NOS}_i)(1-r_i)\right)\times r_a \times p^\textrm{female} \times \textrm{Fec}_a, with rar_a as the maturity at age.

The density-independent component of the survival equation is controlled by α\alpha and the density-dependent component of survival is controlled by β\beta and scaled by the total number of fry, as well as subyearling hatchery releases (see Hatchery section), in the system.

If there is knife-edge maturity, i.e., all fish mature at the same age, then the equation simplifies to ϕ=SAR×pfemale×Fec\phi = \textrm{SAR} \times p^\textrm{female} \times \textrm{Fec}, with SAR\textrm{SAR} as the marine survival (between 0-1).

With the Ricker stock-recruit relationship, smolt production is

Smolty+1NOS=α×Fryy+1NOS×exp(β[Fryy+1NOS+Fryy+1HOS+ny+1sub])Smolty+1HOS=α×Fryy+1HOS×exp(β[Fryy+1NOS+Fryy+1HOS+ny+1sub])\begin{align} \textrm{Smolt}^\textrm{NOS}_{y+1} &= \alpha \times \textrm{Fry}^\textrm{NOS}_{y+1}\times\exp(-\beta[\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1} + n^\textrm{sub}_{y+1}])\\ \textrm{Smolt}^\textrm{HOS}_{y+1} &= \alpha \times \textrm{Fry}^\textrm{HOS}_{y+1}\times\exp(-\beta[\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1} + n^\textrm{sub}_{y+1}]) \end{align}

where α=κ/ϕ\alpha = \kappa/\phi and β=1/Smax\beta = 1/{S_\textrm{max}}.

Habitat improvement

Habitat improvement can improve either the productivity parameter, capacity parameter, or both.

Base terms α\alpha and β\beta are used in the historical period of the model, e.g., conditioning from data. In the projection, the stock recruit parameters α\alpha^\prime and β\beta^\prime are updated with the improvement parameters, IκI_\kappa or ICI_C, as specified in the operating model.

The corresponding parameters in the projection is:

α=κϕ×Iκβ={αCsmolt×ICBeverton-Holt1Smax×ICRicker\begin{align} \alpha^\prime &= \dfrac{\kappa}{\phi}\times I_\kappa\\ \beta^\prime &= \begin{cases} \dfrac{\alpha^\prime}{C_\textrm{smolt} \times I_C} & \textrm{Beverton-Holt}\\ \dfrac{1}{S_\textrm{max} \times I_C} & \textrm{Ricker}\\ \end{cases} \end{align}

Comparison of the change in two density-dependent stock-recruit functions if the improvement parameter is applied to either the compensation ratio (kappa) or capacity parameter. The base value of kappa is 3, Beverton-Holt capacity is 1000, and Ricker Smax is 500. The dotted line is the one-to-one unfished replacement line (corresponding to 1/ϕ1/\phi).

Hatchery production

Hatchery production is controlled by several sets of variables specified by the analyst.

The first consideration is to specify the target number of annual releases of sub-yearlings ntargetsubyearlingn^\textrm{subyearling}_\textrm{target} and yearlings ntargetyearlingn^\textrm{yearling}_\textrm{target}. Going backwards, the corresponding number of eggs needed to reach the target number depends on the egg survival to those life stages. Finally, the corresponding number of broodtake is dependent on the brood fecundity and hatchery survival of broodtake.

Additional considerations are the composition (natural vs. hatchery origin) of the broodtake. To minimize genetic drift of the population due to hatchery production, it is desirable to maintain a high proportion of natural origin broodtake. This is controlled by ptargetNOBp^\textrm{NOB}_\textrm{target}, the desired proportion of natural broodtake relative to all broodtake, but can be exceeded if there is insufficient escapement of natural origin fish.

Second, it is also desirable to maintain high spawning of natural origin fish. This is controlled by pmaxNOBp^\textrm{NOB}_\textrm{max}, the maximum allowable proportion of the natural origin escapement to be used as broodtake. This value is never exceeded.

The following equations then generate the annual broodtake and hatchery production from the state variables given these constraints.

Broodtake

The annual target egg production for the hatchery is calculated from the target releases as

Eggtarget,broodtake=ntargetyearlingsyearling+ntargetsubyearlingssubyearling \textrm{Egg}_\textrm{target,broodtake} = \dfrac{n^\textrm{yearling}_\textrm{target}}{s^\textrm{yearling}} + \dfrac{n^\textrm{subyearling}_\textrm{target}}{s^\textrm{subyearling}}

where ss is the corresponding survival term from the egg life stage.

From the escapement in year yy, some proportion pbroodtakep^\textrm{broodtake} is used as broodtake:

NOBy,a=pybroodtake,NOB×NORy,aescapementHOBy,a=pybroodtake,HOB×phatchery×HORy,aescapement\begin{align} \textrm{NOB}_{y,a} &= p^\textrm{broodtake,NOB}_y \times \textrm{NOR}^\textrm{escapement}_{y,a}\\ \textrm{HOB}_{y,a} &= p^\textrm{broodtake,HOB}_y \times p^\textrm{hatchery} \times \textrm{HOR}^\textrm{escapement}_{y,a} \end{align}

The proportion of the available hatchery fish for broodtake is controlled by phatcheryp^\textrm{hatchery}, which can include fish swimming back to the hatchery or removed from spawning grounds.

The realized hatchery egg production is

EggyNOB=aNOBy,a×sprespawn×pfemale×FecabroodEggyHOB=aHOBy,a×sprespawn×pfemale×Fecabrood\begin{align} \textrm{Egg}_\textrm{y}^\textrm{NOB} &= \sum_a \textrm{NOB}_{y,a} \times s^\textrm{prespawn} \times p^\textrm{female} \times \textrm{Fec}^\textrm{brood}_a\\ \textrm{Egg}_\textrm{y}^\textrm{HOB} &= \sum_a \textrm{HOB}_{y,a} \times s^\textrm{prespawn} \times p^\textrm{female} \times \textrm{Fec}^\textrm{brood}_a \end{align}

where egg production is subject to a survival term sprespawns^\textrm{prespawn} for the broodtake.

The proportion pybroodtakep^\textrm{broodtake}_y is solved annually to satisfy the following conditions:

aNOBy,aaNOBy,a+aHOBy,a=ptargetNOB\dfrac{\sum_a\textrm{NOB}_{y,a}}{\sum_a\textrm{NOB}_{y,a} + \sum_a\textrm{HOB}_{y,a}} = p^\textrm{NOB}_\textrm{target}

0<pybroodtake,HOB<10 < p^\textrm{broodtake,HOB}_y < 1

0<pybroodtake,NOB<pmaxNOB0 < p^\textrm{broodtake,NOB}_y < p^\textrm{NOB}_\textrm{max}

EggyNOB+EggyHOB=Eggbroodtake\textrm{Egg}_\textrm{y}^\textrm{NOB} + \textrm{Egg}_\textrm{y}^\textrm{HOB} = \textrm{Egg}_\textrm{broodtake}

The target ratio ptargetNOBp^\textrm{NOB}_\textrm{target} ensures that there is a sufficiently high proportion of natural origin fish in the broodtake. The maximum removal rate of natural origin fish pmaxNOBp^\textrm{NOB}_\textrm{max} ensures that there is high abundance of natural origin spawners.

The total egg production in a given year can fail to reach the target if there is insufficient natural origin escapement. In this case, the NOB take is set to the maximum removal rate (pybroodtake,NOB=pmaxNOBp^\textrm{broodtake,NOB}_y = p^\textrm{NOB}_\textrm{max}), and the remaining deficit in egg production is met using HOB.

Smolt releases

After the total hatchery egg production is calculated, the production of yearlings and subyearlings is calculated to ensure the annual ratio is equal to the target ratio. To do so, the parameter pyegg,yearlingp^\textrm{egg,yearling}_y is solved subject to the following conditions:

Eggbrood,y=EggyNOB+EggyHOB\textrm{Egg}_\textrm{brood,y} = \textrm{Egg}_\textrm{y}^\textrm{NOB} + \textrm{Egg}_\textrm{y}^\textrm{HOB}

ny+1yearling=pyegg,yearling×Eggbrood,y×syearlingn^\textrm{yearling}_{y+1} = p^\textrm{egg,yearling}_y \times \textrm{Egg}_\textrm{brood,y} \times s^\textrm{yearling}

ny+1subyearling=(1pyegg,yearling)×Eggbrood,y×ssubyearlingn^\textrm{subyearling}_{y+1} = (1 - p^\textrm{egg,yearling}_y) \times \textrm{Egg}_\textrm{brood,y} \times s^\textrm{subyearling}

nyyearlingnysubyearling+nyyearling=ntargetyearlingntargetsubyearling+ntargetyearling\frac{n^\textrm{yearling}_y}{n^\textrm{subyearling}_y + n^\textrm{yearling}_y} = \frac{n^\textrm{yearling}_\textrm{target}}{n^\textrm{subyearling}_\textrm{target} + n^\textrm{yearling}_\textrm{target}}

From the total broodtake, the smolt releases is calculated as

Smolty+1Rel=ny+1yearling+α×ny+1subyearling1+β(Fryy+1NOS+Fryy+1HOS+ny+1subyearling) \textrm{Smolt}^\textrm{Rel}_{y+1} = n^\textrm{yearling}_{y+1} + \frac{\alpha \times n^\textrm{subyearling}_{y+1}}{1 + \beta(\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1} + n^\textrm{subyearling}_{y+1})}

or

Smolty+1Rel=ny+1yearling+α×ny+1subyearling×exp(β(Fryy+1NOS+Fryy+1HOS+ny+1subyearling)) \textrm{Smolt}^\textrm{Rel}_{y+1} = n^\textrm{yearling}_{y+1} + \alpha \times n^\textrm{subyearling}_{y+1} \times \exp(-\beta(\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1} + n^\textrm{subyearling}_{y+1}))

Subyearlings are subject to density-dependent survival in competition with natural production of fry.

Pre-terminal fishery

Let Ny,ajuvN^\textrm{juv}_{y,a} be the juvenile abundance in the population and Ny,a=1juv,NOS=SmoltyNOS+SmoltyHOSN^\textrm{juv,NOS}_{y,a=1} = \textrm{Smolt}^\textrm{NOS}_y + \textrm{Smolt}^\textrm{HOS}_y and Ny,a=1juv,HOS=SmoltRelN^\textrm{juv,HOS}_{y,a=1} = \textrm{Smolt}^\textrm{Rel}. The superscript for the smolt variable corresponds to the parentage while the superscript for NN denotes the origin of the current cohort.

Harvest uPTu^\textrm{PT} in the pre-terminal (PT\textrm{PT}) fishery, assuming no mark-selective fishing, is modeled as a seasonal process. The kept catch KK is

Ky,aNOS,PT=(1exp(vaPTFyPT))Ny,ajuv,NOSKy,aHOS,PT=(1exp(vaPTFyPT))Ny,ajuv,HOS\begin{align} K^\textrm{NOS,PT}_{y,a} &= \left(1 - \exp(-v^\textrm{PT}_a F^\textrm{PT}_y)\right)N^\textrm{juv,NOS}_{y,a}\\ K^\textrm{HOS,PT}_{y,a} &= \left(1 - \exp(-v^\textrm{PT}_a F^\textrm{PT}_y)\right)N^\textrm{juv,HOS}_{y,a}\\ \end{align}

The instantaneous fishing mortality is solved to meet the following condition

uPT=aKy,aNOS,PT+aKy,aHOS,PTaNy,ajuv,NOS+aNy,ajuv,HOS u^\textrm{PT} = \dfrac{\sum_a K^\textrm{NOS,PT}_{y,a} + \sum_a K^\textrm{HOS,PT}_{y,a}}{\sum_a N^\textrm{juv,NOS}_{y,a} + \sum_a N^\textrm{juv,HOS}_{y,a}}

Recruitment and maturity

The recruitment is calculated from the survival of juvenile fish after pre-terminal harvest and maturation:

NORy,a=Ny,ajuv,NOSexp(vaFyPT)ry,aHORy,a=Ny,ajuv,HOSexp(vaFyPT)ry,a\begin{align} \textrm{NOR}_{y,a} &= N^\textrm{juv,NOS}_{y,a}\exp(-v_aF^\textrm{PT}_y)r_{y,a}\\ \textrm{HOR}_{y,a} &= N^\textrm{juv,HOS}_{y,a}\exp(-v_aF^\textrm{PT}_y)r_{y,a} \end{align}

The juvenile abundance in the following year consists of fish that did not mature and subsequently survived natural mortality MM:

Ny+1,a+1juv,NOS=Ny,ajuv,NOSexp([vaFyPT+My,aNOS])(1ry,a)Ny+1,a+1juv,HOS=Ny,ajuv,HOSexp([vaFyPT+My,aHOS])(1ry,a)\begin{align} N^\textrm{juv,NOS}_{y+1,a+1} &= N^\textrm{juv,NOS}_{y,a}\exp\left(-[v_aF^\textrm{PT}_y + M^\textrm{NOS}_{y,a}]\right)(1 - r_{y,a})\\ N^\textrm{juv,HOS}_{y+1,a+1} &= N^\textrm{juv,HOS}_{y,a}\exp\left(-[v_aF^\textrm{PT}_y + M^\textrm{HOS}_{y,a}]\right)(1 - r_{y,a}) \end{align}

Natural mortality is specified by age class. Accordingly, this mortality corresponds to either the freshwater or marine survival depending on age class.

Terminal fishery

Assuming no mark-selective fishing, the retained catch of the terminal (T\textrm{T}) fishery is calculated from the harvest rate similarly as with the pre-terminal fishery:

Ky,aNOS,T=(1exp(vaTFyT))NORy,aKy,aHOS,T=(1exp(vaTFyT))HORy,a\begin{align} K^\textrm{NOS,T}_{y,a} &= \left(1 - \exp(-v^\textrm{T}_a F^\textrm{T}_y)\right)\textrm{NOR}_{y,a}\\ K^\textrm{HOS,T}_{y,a} &= \left(1 - \exp(-v^\textrm{T}_a F^\textrm{T}_y)\right)\textrm{HOR}_{y,a}\\ \end{align}

subject to

uT=aKy,aNOS,T+aKy,aHOS,TaNORy,a+aHORy,a u^\textrm{T} = \dfrac{\sum_a K^\textrm{NOS,T}_{y,a} + \sum_a K^\textrm{HOS,T}_{y,a}}{\sum_a \textrm{NOR}_{y,a} + \sum_a \textrm{HOR}_{y,a}}

Escapement and spawners

The escapement consists of the survivors of the terminal fishery:

NORy,aescapement=NORy,aexp(vaFyT)HORy,aescapement=HORy,aexp(vaFyT)\begin{align} \textrm{NOR}^\textrm{escapement}_{y,a} &= \textrm{NOR}_{y,a}\exp(-v_aF^\textrm{T}_y)\\ \textrm{HOR}^\textrm{escapement}_{y,a} &= \textrm{HOR}_{y,a}\exp(-v_aF^\textrm{T}_y) \end{align}

The abundance of natural origin spawners consists of the escapement reduced by the broodtake.

NOSy,a=(1pybroodtake,NOB)NORy,aescapement \textrm{NOS}_{y,a} = (1 - p^\textrm{broodtake,NOB}_y) \textrm{NOR}^\textrm{escapement}_{y,a}

The hatchery origin spawners is the escapement, reduced by the proportion phatcheryp^\textrm{hatchery} returning to the hatchery, either by swim-in facilities or direct removal. A second removal term premovalHOSp^\textrm{HOS}_\textrm{removal} removes fish from the spawning grounds, these animals are not available for broodtake.

HOSy,a=HORy,aescapement(1phatchery)(1premovalHOS) \textrm{HOS}_{y,a} = \textrm{HOR}^\textrm{escapement}_{y,a} (1 - p^\textrm{hatchery}) (1 - p^\textrm{HOS}_\textrm{removal})

Fitness effects on survival

Reproductive success of first generation hatchery fish has been observed to be lower than their natural counterparts, and is accounted for in the γ\gamma parameter (see review in Withler et al. 2018).

Through genetic and epigenetic factors, survival of hatchery juveniles in the hatchery environment selects for fish with a phenotype best adapted for that environment, and likewise for juveniles spawned in the natural environment. Since these traits are heritable, the fitness of the natural population can shift away from the optimum for the natural environment towards that of the hatchery environment on an evolutionary time scale, i.e., over a number of generations, when hatchery fish are allowed to spawn.

As described in Ford 2002 and derived in Lande 1976, the fitness loss function WW for an individual with phenotypic trait value zz in a given environment is

W(z)=exp((zθ)22ω2) W(z) = \exp\left(\dfrac{-(z-\theta)^2}{2\omega^2}\right)

where θ\theta is the optimum for that environment and ω2\omega^2 is the fitness variance.

If the phenotypic trait value zz in the population is a random normal variable with mean z\bar{z} and variance σ2\sigma^2, then the mean fitness of the population in generation gg is W(z)=W(z)f(z)dz\bar{W}(z) = \int W(z) f(z) dz, where f(z)f(z) is the Gaussian probability density function. The solution is proportional to

W(z)exp((zθ)22(ω+σ)2) \bar{W}(z) \propto \exp\left(\dfrac{-(\bar{z}-\theta)^2}{2(\omega+\sigma)^2}\right)

The mean phenotype z\bar{z} is calculated iteratively, where the change Δz\Delta\bar{z} from generation g1g-1 to gg is

Δz=zgzg1=(zg1zg1)h2zg=zg1+(zg1zg1)h2\begin{align} \Delta\bar{z} &= \bar{z}_g - \bar{z}_{g-1} = (\bar{z}^\prime_{g-1} - \bar{z}_{g-1})h^2\\ \bar{z}_g &= \bar{z}_{g-1} + (\bar{z}^\prime_{g-1} - \bar{z}_{g-1})h^2\\ \end{align}

where h2h^2 is the heritability of zz and zg1\bar{z}^\prime_{g-1} is the trait value after applying the fitness function, defined as:

zg1=1Wg1Wg1(z)×zf(z)dz=zg1ω2+θσ2ω2+σ2\begin{align} \bar{z}^\prime_{g-1} &= \dfrac{1}{\bar{W}_{g-1}}\int W_{g-1}(z)\times zf(z)dz\\ &= \dfrac{\bar{z}_{g-1}\omega^2 + \theta \sigma^2}{\omega^2 + \sigma^2} \end{align}

Let zg1(θ)\bar{z}^\prime_{g-1}(\theta) be a function that returns the mean trait value after selection in an environment with optimum value θ\theta. With a hatchery program, the mean trait value of the progeny in the natural environment is a weighted average of the mean trait value in natural and hatchery origin spawners, with selection in the natural environment, i.e., with optimum trait value θnatural\theta^\textrm{natural}:

zgnatural=(1pg1HOSeff)×(zg1natural+[zg1natural(θnatural)zg1natural]h2)+pg1HOSeff×(zg1hatchery+[zg1hatchery(θnatural)zg1hatchery]h2)\begin{align} \bar{z}^\textrm{natural}_g = & (1 - p^\textrm{HOSeff}_{g-1}) \times \left(\bar{z}^\textrm{natural}_{g-1} + [\bar{z}^{\prime\textrm{natural}}_{g-1}(\theta^\textrm{natural}) - \bar{z}^\textrm{natural}_{g-1}] h^2\right) +\\ & p^\textrm{HOSeff}_{g-1} \times \left(\bar{z}^\textrm{hatchery}_{g-1} + [\bar{z}^{\prime\textrm{hatchery}}_{g-1}(\theta^\textrm{natural}) - \bar{z}^\textrm{hatchery}_{g-1}] h^2\right) \end{align}

where pHOSeff=γ×HOS/(NOS+γ×HOS)p^\textrm{HOSeff} = \gamma\times\textrm{HOS}/(\textrm{NOS} + \gamma\times\textrm{HOS}).

Similarly, the mean trait value in the hatchery environment zghatchery\bar{z}^\textrm{hatchery}_g is a weighted average of the mean trait value of the natural and hatchery broodtake, with selection in the hatchery environment, i.e., with optimum trait value θhatchery\theta^\textrm{hatchery}:

zghatchery=pg1NOB×(zg1natural+[zg1natural(θhatchery)zg1natural]h2)+(1pg1NOB)×(zg1hatchery+[zg1hatchery(θhatchery)zg1hatchery]h2)\begin{align} \bar{z}^\textrm{hatchery}_g = & p^\textrm{NOB}_{g-1} \times \left(\bar{z}^\textrm{natural}_{g-1} + [\bar{z}^{\prime\textrm{natural}}_{g-1}(\theta^\textrm{hatchery}) - \bar{z}^\textrm{natural}_{g-1}] h^2\right) +\\ & (1 - p^\textrm{NOB}_{g-1}) \times \left(\bar{z}^\textrm{hatchery}_{g-1} + [\bar{z}^{\prime\textrm{hatchery}}_{g-1}(\theta^\textrm{hatchery}) - \bar{z}^\textrm{hatchery}_{g-1}] h^2\right) \end{align}

where pNOB=NOB/(NOB+HOB)p^\textrm{NOB} = \textrm{NOB}/(\textrm{NOB} + \textrm{HOB}).

The fitness variance ω2\omega^2 and phenotype variance σ2\sigma^2 are assumed constant between the two environments. The fitness variance is parameterized relative to the variance of the phenotype, i.e., ω=Ωσ\omega = \Omega \sigma where Ω\Omega is the “selection strength”.

The mean fitness of generation gg in the natural environment is then:

Wgnatural=exp((zgnaturalθnatural)22(ω+σ)2) \bar{W}^\textrm{natural}_g = \exp\left(\dfrac{-(\bar{z}^\textrm{natural}_g-\theta^\textrm{natural})^2}{2(\omega+\sigma)^2}\right)

Mixed brood-year return

If a mixed-brood year return in year yy across several ages aa produces the smolt cohort in year y+1y+1, then the mean trait value in the progeny is calculated from a weighted average by brood year and age class fecundity:

zy+1natural=apy,aNOS×(zya+1natural+[zya+1natural(θnatural)zya+1natural]h2)+apy,aHOSeff×(zya+1hatchery+[zya+1hatchery(θnatural)zya+1hatchery]h2)\begin{align} \bar{z}^\textrm{natural}_{y+1} = & \sum_a p^\textrm{NOS}_{y,a} \times \left(\bar{z}^\textrm{natural}_{y-a+1} + [\bar{z}^{\prime\textrm{natural}}_{y-a+1}(\theta^\textrm{natural}) - \bar{z}^\textrm{natural}_{y-a+1}] h^2\right) +\\ & \sum_a p^\textrm{HOSeff}_{y,a} \times \left(\bar{z}^\textrm{hatchery}_{y-a+1} + [\bar{z}^{\prime\textrm{hatchery}}_{y-a+1}(\theta^\textrm{natural}) - \bar{z}^\textrm{hatchery}_{y-a+1}] h^2\right) \end{align}

zy+1hatchery=apy,aNOB×(zya+1natural+[zya+1natural(θhatchery)zya+1natural]h2)+apy,aHOB×(zya+1hatchery+[zya+1hatchery(θhatchery)zya+1hatchery]h2)\begin{align} \bar{z}^\textrm{hatchery}_{y+1} = & \sum_a p^\textrm{NOB}_{y,a} \times \left(\bar{z}^\textrm{natural}_{y-a+1} + [\bar{z}^{\prime\textrm{natural}}_{y-a+1}(\theta^\textrm{hatchery}) - \bar{z}^\textrm{natural}_{y-a+1}] h^2\right) +\\ & \sum_a p^\textrm{HOB}_{y,a} \times \left(\bar{z}^\textrm{hatchery}_{y-a+1} + [\bar{z}^{\prime\textrm{hatchery}}_{y-a+1}(\theta^\textrm{hatchery}) - \bar{z}^\textrm{hatchery}_{y-a+1}] h^2\right) \end{align}

where

py,aNOS=Feca×NOSy,aaFeca(NOSy,a+γ×HOSy,a)p^\textrm{NOS}_{y,a} = \dfrac{\textrm{Fec}_a \times \textrm{NOS}_{y,a}}{\sum_a\textrm{Fec}_a(\textrm{NOS}_{y,a} + \gamma \times \textrm{HOS}_{y,a})}

py,aHOSeff=Feca×γ×HOSy,aaFeca(NOSy,a+γ×HOSy,a)p^\textrm{HOSeff}_{y,a} = \dfrac{\textrm{Fec}_a \times \gamma \times \textrm{HOS}_{y,a}}{\sum_a\textrm{Fec}_a(\textrm{NOS}_{y,a} + \gamma \times \textrm{HOS}_{y,a})}

py,aNOB=Fecabrood×NOBy,aaFecabrood(NOBy,a+HOBy,a)p^\textrm{NOB}_{y,a} = \dfrac{\textrm{Fec}^\textrm{brood}_a \times \textrm{NOB}_{y,a}}{\sum_a\textrm{Fec}^\textrm{brood}_a(\textrm{NOB}_{y,a} + \textrm{HOB}_{y,a})}

py,aHOB=Fecabrood×HOBy,aaFecabrood(NOBy,a+HOBy,a)p^\textrm{HOB}_{y,a} = \dfrac{\textrm{Fec}^\textrm{brood}_a \times \textrm{HOB}_{y,a}}{\sum_a\textrm{Fec}^\textrm{brood}_a(\textrm{NOB}_{y,a} + \textrm{HOB}_{y,a})}

Effective proportions, i.e., weighting by age-class fecundity, accounts for older age classes that are more fecund and more likely to contribute to the production of next cohort.

Fitness loss

Fitness can reduce survival in the egg, fry, and immature life stages:

Fryy+1NOS=EggyNOS×(Wy+1nat.)eggFryy+1HOS=EggyHOS×(Wy+1nat.)egg\begin{align} \textrm{Fry}^\textrm{NOS}_{y+1} &= \textrm{Egg}^\textrm{NOS}_y \times (W^\textrm{nat.}_{y+1})^{\ell_\textrm{egg}}\\ \textrm{Fry}^\textrm{HOS}_{y+1} &= \textrm{Egg}^\textrm{HOS}_y \times (W^\textrm{nat.}_{y+1})^{\ell_\textrm{egg}} \end{align}

Smolty+1NOS=α×Fryy+1NOS1+β(Fryy+1NOS+Fryy+1HOS)Smolty+1HOS=α×Fryy+1HOS1+β(Fryy+1NOS+Fryy+1HOS)\begin{align} \textrm{Smolt}^\textrm{NOS}_{y+1} &= \frac{\alpha^{\prime\prime} \times \textrm{Fry}^\textrm{NOS}_{y+1}}{1 + \beta^{\prime\prime} (\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1})}\\ \textrm{Smolt}^\textrm{HOS}_{y+1} &= \frac{\alpha^{\prime\prime} \times \textrm{Fry}^\textrm{HOS}_{y+1}}{1 + \beta^{\prime\prime} (\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1})} \end{align}

My,aNOS=log(exp(My,abase,NOS)×(Wyanat.)smolt) M^\textrm{NOS}_{y,a} = -\log(\exp(-M^\textrm{base,NOS}_{y,a}) \times (W^\textrm{nat.}_{y-a})^{\ell_\textrm{smolt}})

with i\ell_i is the proportion of the fitness loss apportioned among the three life stages, ii=1\sum_i \ell_i = 1, and density-dependent parameters α=κϕ×Iκ×(Wy+1nat.)fry\alpha^{\prime\prime} = \frac{\kappa}{\phi} \times I_\kappa \times (W^\textrm{nat.}_{y+1})^{\ell_\textrm{fry}}, β=α/[Csmolt×IC×(Wy+1nat.)fry]\beta^{\prime\prime} = \alpha^{\prime\prime}/[C_\textrm{smolt} \times I_C \times (W^\textrm{nat.}_{y+1})^{\ell_\textrm{fry}}].

In the marine environment, age-specific natural survival is reduced proportional to the fitness loss term and modeled as a cohort effect.

With the Ricker density-dependent survival,

Smolty+1NOS=α×Fryy+1NOS×exp(β[Fryy+1NOS+Fryy+1HOS])Smolty+1HOS=α×Fryy+1HOS×exp(β[Fryy+1NOS+Fryy+1HOS])\begin{align} \textrm{Smolt}^\textrm{NOS}_{y+1} &= \alpha^{\prime\prime} \times \textrm{Fry}^\textrm{NOS}_{y+1} \times \exp(-\beta^{\prime\prime} [\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1}])\\ \textrm{Smolt}^\textrm{HOS}_{y+1} &= \alpha^{\prime\prime} \times \textrm{Fry}^\textrm{HOS}_{y+1} \times \exp(-\beta^{\prime\prime} [\textrm{Fry}^\textrm{NOS}_{y+1} + \textrm{Fry}^\textrm{HOS}_{y+1}]) \end{align}

with β=1/[Smax×IC×(Wy+1nat.)fry]\beta^{\prime\prime} = 1/[S_\textrm{max} \times I_C \times (W^\textrm{nat.}_{y+1})^{\ell_\textrm{fry}}].

PNI

PNI (proportionate natural influence) is an approximation of the rate of gene flow from the hatchery to the natural environment, calculated for the progeny in year y+1y+1 from the parental composition of year yy:

PNIy+1=apy,aNOBapy,aNOB+apy,aHOSeff \textrm{PNI}_{y+1} = \dfrac{\sum_a p^{\textrm{NOB}}_{y,a}}{\sum_a p^{\textrm{NOB}}_{y,a} + \sum_a p^{\textrm{HOSeff}}_{y,a}}

Generally, a combination of minimizing hatchery releases, increasing natural broodtake, and reducing the number hatchery origin spawners maintains high PNI, i.e., low rate of gene flow from the hatchery to natural environment.

Wild salmon

With single brood-year returns, the proportion of wild salmon, natural origin spawners whose parents were also natural spawners, can be calculated as

pgWILD=(1pgHOScensus)×(1pg1HOScensus)2(1pg1HOScensus)2+2γ×pg1HOScensus(1pg1HOScensus)+γ2(pg1HOScensus)2 p^\textrm{WILD}_g = (1 - p^\textrm{HOScensus}_g) \times \dfrac{(1 - p^\textrm{HOScensus}_{g-1})^2} {(1 - p^\textrm{HOScensus}_{g-1})^2 + 2 \gamma \times p^\textrm{HOScensus}_{g-1}(1 - p^\textrm{HOScensus}_{g-1}) + \gamma^2 (p^\textrm{HOScensus}_{g-1})^2}

where pHOScensus=HOS/(HOS+NOS)p^\textrm{HOScensus} = \textrm{HOS}/(\textrm{HOS} + \textrm{NOS}).

The first term is the proportion of natural spawners in the current generation gg.

The ratio comprising the second term discounts the proportion of the current generation to include natural spawners whose parents were both natural spawners. Assuming non-assortative mating, the three terms in the denominator gives the composition of generation gg whose parents who are both natural origin, mixed origin (one parent in natural origin and the other is hatchery origin), and both hatchery origin.

To generalize for mixed-brood year return, we calculate the probability weighted across brood-years and age class fecundity:

pyWILD=aNOSy,aa(NOSy,a+HOSy,a)×(apya,aNOScensus)2(apya,aNOScensus)2+2γ×(apya,aNOScensus)(apya,aHOScensus)+γ2(apya,aHOScensus)2 p^\textrm{WILD}_y = \sum_a \dfrac{\textrm{NOS}_{y,a}}{\sum_{a'}(\textrm{NOS}_{y,a'} + \textrm{HOS}_{y,a'})} \times \dfrac{(\sum_{a'} p^\textrm{NOScensus}_{y-a,a'})^2} {(\sum_{a'} p^\textrm{NOScensus}_{y-a,a'})^2 + 2 \gamma \times (\sum_{a'}p^\textrm{NOScensus}_{y-a,a'})(\sum_{a'}p^\textrm{HOScensus}_{y-a,a'}) + \gamma^2 (\sum_{a'}p^\textrm{HOScensus}_{y-a,a'})^2}

where

py,aNOScensus=Feca×NOSy,aaFeca(NOSy,a+HOSy,a)p^\textrm{NOScensus}_{y,a} = \dfrac{\textrm{Fec}_a \times \textrm{NOS}_{y,a}}{\sum_a{\textrm{Fec}_a (\textrm{NOS}_{y,a}} + \textrm{HOS}_{y,a})}

py,aHOScensus=Feca×HOSy,aaFeca(NOSy,a+HOSy,a)p^\textrm{HOScensus}_{y,a} = \dfrac{\textrm{Fec}_a \times \textrm{HOS}_{y,a}}{\sum_a{\textrm{Fec}_a (\textrm{NOS}_{y,a}} + \textrm{HOS}_{y,a})}

The probability of finding a wild salmon in year yy is the sum of probabilities of finding a wild salmon over all ages. For each age aa, the first ratio is the probability of finding a natural spawner in year yy. The second ratio is the probability of mating success from two parental natural spawners in year yay-a using a Punnett square, assuming non-assortative mating across age and origin. The summation across dummy age variable aa' calculates the total proportion of spawners in a given year.

Effective proportions, i.e., weighting by age-class fecundity, in the parental composition accounts for older age classes that are more fecund and more likely to contribute to the production of offspring.

Mark-selective fishing

If the mark rate mm of hatchery fish is greater than zero, then mark-selective fishing is implemented for both the pre-terminal and terminal fisheries. The mark rate is a proxy for retention and the harvest rate uharvestu^\textrm{harvest} corresponds to the ratio of the kept catch and abundance. The exploitation rate uexploitu^\textrm{exploit} is calculated from kept catch and dead releases. Exploitation rates differ between hatchery and natural origin fish because there is no retention of the latter.

Let the instantaneous fishing mortality for kept catch and released catch be

Fkept=mEFrel.=(1m)δE\begin{align} F^\textrm{kept} &= mE\\ F^\textrm{rel.} &= (1 - m)\delta E \end{align}

where δ\delta is the proportion of released fish that die.

EE is an index of fishing effort, also referred to as the encounter rate by the fishery, that links together FkeptF^\textrm{kept} and Frel.F^\textrm{rel.}. Intuitively, fishing effort can increase in a mark-selective fishery compared to a non-selective fishery. For example, if the mark rate is 20 percent, then the fishing effort could be 500 percent higher than in a non-selective fishery in order to attain the kept quota or bag limit. Additional catch and release mortality then occurs for un-marked fish, according to δ\delta.

In the pre-terminal (PT\textrm{PT}) fishery, EE is solved to satisfy the following equation for hatchery fish:

uharvest,HOS,PT=aKy,aHOS,PTaNy,ajuv,HOS u^\textrm{harvest,HOS,PT} = \dfrac{\sum_aK^\textrm{HOS,PT}_{y,a}}{\sum_a N^\textrm{juv,HOS}_{y,a}}

where the kept catch KK is

Ky,aHOS,PT=Fykept,PTFykept,PT+Fyrel,PT(1exp(vaPT[Fkept,PT+Frel,PT]))Ny,ajuv,HOSK^\textrm{HOS,PT}_{y,a} = \dfrac{F^\textrm{kept,PT}_y}{F^\textrm{kept,PT}_y + F^\textrm{rel,PT}_y}\left(1 - \exp(-v^\textrm{PT}_a[F^\textrm{kept,PT} + F^\textrm{rel,PT}])\right)N^\textrm{juv,HOS}_{y,a}.

The exploitation rate for natural origin fish is calculated from dead discards. The exploitation rate for hatchery origin fish is calculated from kept catch and dead discards:

uyexploit,NOS,PT=a(1exp(vaFyrel.,PT))Ny,ajuv,NOSaNy,ajuv,NOSuyexploit,HOS,PT=a(1exp(va[Fykept,PT+Fyrel.,PT]))Ny,ajuv,HOSaNy,ajuv,HOS\begin{align} u^\textrm{exploit,NOS,PT}_y &= \dfrac{\sum_a (1 - \exp(-v_a F^\textrm{rel.,PT}_y))N^\textrm{juv,NOS}_{y,a}}{\sum_a N^\textrm{juv,NOS}_{y,a}}\\ u^\textrm{exploit,HOS,PT}_y &= \dfrac{\sum_a (1 - \exp(-v_a[F^\textrm{kept,PT}_y + F^\textrm{rel.,PT}_y]))N^\textrm{juv,HOS}_{y,a}}{\sum_a N^\textrm{juv,HOS}_{y,a}} \end{align}

Similarly, in the terminal fishery, the fishing effort satisfies the equation

uharvest,HOS,T=aKy,aHOS,TaHORy,a u^\textrm{harvest,HOS,T} = \dfrac{\sum_aK^\textrm{HOS,T}_{y,a}}{\sum_a \textrm{HOR}_{y,a}}

with the corresponding exploitation rates:

uyexploit,NOS,T=a(1exp(vaFyrel,T))NORy,aaNORy,auyexploit,HOS,T=a(1exp(va[Fykept,T+Fyrel.,T]))HORy,aaHORy,a\begin{align} u^\textrm{exploit,NOS,T}_y &= \dfrac{\sum_a (1 - \exp(-v_a F^\textrm{rel,T}_y))\textrm{NOR}_{y,a}}{\sum_a \textrm{NOR}_{y,a}}\\ u^\textrm{exploit,HOS,T}_y &= \dfrac{\sum_a (1 - \exp(-v_a[F^\textrm{kept,T}_y + F^\textrm{rel.,T}_y]))\textrm{HOR}_{y,a}}{\sum_a \textrm{HOR}_{y,a}} \end{align}