To simplify the following presentation of population dynamics relationships, we supress subscripts for reach and species. It is to be understood that all age-structure accounting relationships (for animals age 1 yr and older) are done at system-year scale (once per year for each species), and all recruitment relationships (egg/hatchling to age 1 yr) at reach-month scale.

The age-structure population dynamics accounting is based on the Deriso(1980)-Schnute (1987) delay-difference equation structure, which is used widely in fisheries assessment (see its derivation in Hilborn and Walters 1992). For animals age 1 yr to an age at maturity k, numbers at each age are propagated using a simple survival relationship N(a + 1,tS_jN(a,t)

where S_j is an annual juvenile survival rate. For age k and older animals, we propagate total numbers N (summed over all ages, ignoring age-specific differences in survival rate) and adult biomass B using the delay-difference relations:

where t is time in years; S_a is adult annual survival rate (constant except for harvesting/control effects); R(t + 1) is age k recruitment in number of animals; A, ρ are Ford-Brody body growth model parameters; and w_k is mean body weight (mass) at maturity.

It is not necessary to propagate these equations at reach spatial scale because we do not assume reach-specific differences in age 1+ growth and survival rates (except for time-varying harvest effects on S_a for rainbow trout in the reach just below Glen Canyon dam, where most of the rainbow population is found). The age 1 numbers entering the juvenile age accounting each year (leading at age k to total recruits R) are calculated as a sum, over reaches, of recruitment rates by reach.

Adults using each reach (a proportion of N,B based initially on historical data and then updated in relation to simulated recruitment success by reach over years) are assumed to produce eggs/hatchlings E(t) proportional to biomass, E(t) = eB(t) where e = eggs/adult biomass. Derivations in Walters and Korman (1998) imply that we can predict the number of age-1 recruits resulting from E(t), while accounting for monthly changes in predation risk related to food supply and growth, by a Beverton-Holt (1957) function of the form:

where S_egg is egg survival rate, assumed to be constant; m is month of age (after hatching) m = 1…12; v_m is relative vulnerability of m-month-old juveniles to predation; S_m(t) is maximum (low-density) survival rate over the month; K₁,K₂ are scaling constants that depend on units of measurement (see below); P(m) is an index of total predation risk per time foraging; F(m) is a weighted (by diet composition) index of food density; H_r(m) is a function defining relative (0-1) variation in the risk ratio effect (P/F effect on survival) with changes in physical habitat factors; and H_f(m) is a function defining relative (0-1) variation in usable foraging habitat (foraging arena size) with changes in physical habitat factors.

There are two key components in (A13.3): trophic linkage to predation risk and food availability via the risk ratio P(m)/F(m), and linkage to physical habitat factors that affect both predation risk and food competition via the habitat functions H_r and H_f.

The P and F functions are weighted sums over predators and food types of predators, and food-specific relative risks/opportunities. In initializing each simulation for baseline year 1993, we scale these functions so that P = 1, F = 1 if predator and food abundances identical to 1993 values reoccur in other simulation years. Thus, for example, if we have assumed that a particular predator accounts for 20% of the juvenile mortality in the 1993 baseline situation, we take its contribution to P in other years to be 0.2B(t)/B(1993), where B(t) is predator biomass in year t. Contributions to F are calculated the same way. We attempted a spreadsheet trophic modeling exercise to see if we could replace such relative predation and feeding rate calculations with dimensional mortality rate parameters based on absolute estimates of predator abundances, feeding rates, and diet compositions. This exercise revealed really discouraging gaps in basic abundance and feeding data, particularly for key exotic species. The relative P and F calculations allow us at least to define clear alternative hypotheses about possible predation and food impacts, without pretending that we can quantify each of the component calculations required to justify any particular hypothesis.

In the first (1993 baseline) simulation year of each scenario, we calculate the equation A13.3 H, P, and F monthly factors for every indicator species, and then use these to calculate K₁ and K₂ so that N(1,12) (age 1 recruits) will just balance assumed natural mortalities of older animals, provided that habitat, predation, and food conditions remain constant. K₁ is estimated over reaches as a single, nonspatial factor, using observed relative abundances by reach to weight the spatially varying H, P, and F factors. But the K₂ "calibration" is done by default on a reach-by-reach basis, so the K₂ parameter can represent habitat capacity effects in addition to those modeled explicitly with H_f. This default calculation can be turned off by model users (initialization calculates single abundance-weighted average K₂ over reaches), to see if explicitly recognized habitat relationships included in the calculation of H_f can explain observed spatial distribution patterns (this two-stage optional calculation was necessary, particularly in early model development, in order that we could "drive" the submodels for particular species with reasonable spatial relative abundance estimates for predator/competitor species, and work through the habitat linkage parameterization in some systematic way across species).

In fact, K₁ represents a "compensation parameter", in the sense that the product S_m(m) over months represents maximum survival rate of eggs to age 1 yr. Model users can either specify this maximum survival rate directly, or can assume the maximum survival rate to be K*(R₀/E₀), where R₀/E₀ represents the survival rate needed to produce a population equilibrium at the baseline (1993) egg production rate E₀. In the second approach, K* represents the ratio of maximum survival rate at very low population density to the survival rate for a balanced population. There are no empirical data on K* for Grand Canyon vertebrates (no stock-recruitment data); based on comparative studies of fish recruitment relationships (Myers and Barrowman 1996), we have advised users to assume K* values on the order of 3-10.

Each "habitat linkage function" H_r(m) and H_f(m) is calculated as a product of effects over habitat factor values for month m. (Habitat factor values available in the simulation for each reach, and month include temperature, turbidity, max-min stage, wetted area, and warm tributary + backwater area). Using the interface shown in Appendix 8, model users can sketch functions of the form h(V_k), where 0< h <1, and V_k is the value of the kth habitat factor in month m (e.g., V₁ is the max-min stage, V₅ is turbidity); V scales are set relative to 1993 baseline values. Each H is then just the product h(V₁)h(V₂)… over those factors V for which the model user has sketched some relationship (h = 1 for "inactive" V's). By defining H as a product of h(V) effects (rather than, say, a sum of effects), we assume that each habitat factor has a strong effect independent of other factors (a high f value for one habitat variable cannot "compensate" for a low value for another variable, e.g., having a large habitat area cannot compensate for poor thermal conditions). In the same model interface, users can sketch the v_m function (0-1 relative values) defining how sensitivity to both H and P/R varies with age over the first year of life; we initially allowed definition of a separate v_m function for each habitat factor, but we found this greatly complicated the functional specifications without substantially improving our ability to represent alternative hypotheses about the importance of various habitat factors.