Size‐selective fishing and the potential for fisheries‐induced evolution in lake whitefish

Abstract The long‐term evolutionary effects of fishing on maturation schedules can depend on gear type, the shape of the gear type's size‐selectivity function, and the size and age structure of a population. Our goal was to better understand how environmentally induced differences in somatic growth influence the evolutionary effects of size‐selective fisheries, using lake whitefish (Coregonus clupeaformis) in Lake Huron as a case study. Using a state‐dependent optimization model of energy allocation parameterized for lake whitefish, we show that fishing with gill nets (bell‐shaped selectivity) and trap nets (sigmoid‐shaped selectivity) can be potent agents of selection on size thresholds for maturity. Compared to trap nets and large mesh (114 mm) gill nets, small mesh (89 mm) gill nets are better able to buffer populations from fishing‐induced evolution by safeguarding large, fecund fish, but only when overall fishing mortality is low and growth rates sufficiently fast such that fish can outgrow vulnerable size classes. Regardless of gear type, and all else being equal, high fishing mortality in combination with low growth rates is expected to intensify the long‐term evolutionary effects of fishing.


| INTRODUC TI ON
A central tenet of life history theory is that optimal age and size at maturity reflect a balance between the benefits of growing larger and achieving higher reproductive potential and the costs of dying before realizing that reproductive potential (Roff, 1984;Stearns, 1992). Thus, fishing mortality, which can target a large range of ages and sizes, is expected to be an important selection factor on life histories in exploited fish populations, leading to evolutionary changes in growth rate, size and age thresholds for maturation, and fecundity (Dunlop, Enberg, Jørgensen, & Heino, 2009;Heino, Pauli, & Dieckmann, 2015;Jørgensen et al., 2007;Kuparinen & Merilä, 2007). Indeed, evidence has emerged that fishing mortality of stocks, both marine and freshwater, can be sufficiently high for evolutionary changes in life history traits to occur at trackable, ecological timescales (Edeline et al., 2007;Heino et al., 2015;Jørgensen et al., 2007;Nusslé, Bornand, & Wedekind, 2009;Olsen et al., 2004).
For example, Jørgensen et al. (2009) demonstrated that trawl nets, which have sigmoid-shaped selectivity curves, consistently favor early maturation in Atlantic cod (Gadus morhua). In contrast, gill nets, which have bell-shaped size selectivity curves, can select for fish that grow through vulnerable size classes before maturing. This outcome has the potential to reduce fisheries-induced selection for early maturation. For any specific fishery stock, however, this outcome likely will depend on somatic growth rate and the potential for fish to grow through vulnerable stages quickly. Our goal was to better understand how the growth regime influences the evolutionary effects of size-selective fisheries in lake whitefish (Coregonus clupeaformis).
Lake whitefish populations have long experienced fishing mortality and variable growth regimes. In the Laurentian Great Lakes, lake whitefish support subsistence and commercial fisheries, and historically faced high exploitation rates in commercial gill net and trap net fisheries (Brenden, Brown, Ebener, Reid, & Newcomb, 2013;Ebener, Brenden, & Jones, 2010;Ebener et al., 2008). Ebener et al. (2008) provide an historical account of lake whitefish fisheries in the Great Lakes (also see Eberts et al., 2017). High exploitation rate in combination with other factors led to the collapse of lake whitefish fisheries throughout the Great Lakes from 1955 to 1970. Subsequent management actions then led to the recovery of lake whitefish populations and the reestablishment of large commercial fisheries. However, following the invasion of zebra and quagga mussels (Dreissena polymorpha and D. bugensis) and loss of Diporeia spp. amphipods beginning in the 1990s, lake whitefish underwent significant changes in diet and associated declines in growth and recruitment in affected regions (Gobin, Lester, Cottrill, Fox, & Dunlop, 2015;Pothoven & Madenjian, 2008;Pothoven, Nalepa, Schneeberger, & Brandt, 2001;Pothoven et al., 2006;Rennie, Sprules, & Johnson, 2009;Rennie & Verdon, 2008). Particularly in some regions of Lake Huron, food web changes in combination with high sea lamprey (Petromyzon marinus) mortality have contributed to very low fishing effort and yields in recent years (Caroffino & Lenart, 2017).
Prevailing growth and mortality regimes drive patterns of life history variation among lake whitefish populations. Across their range in North America, lake whitefish populations show extensive variation in growth, mortality, and age and size at maturity (Beauchamp, Collins, & Henderson, 2004). For example, maximum total length ranges from 470 to 818 mm and maximum weight ranges from 1,140 to 11,460 g among 419 populations (Rennie & Verdon, 2008). Among populations, faster prereproductive growth correlates with earlier age at maturity, which is consistent with life history theory (Charnov, 1993). Within the Laurentian Great Lakes, large-scale spatial variation in size thresholds for maturation among basins and lakes (e.g., smaller size thresholds in Lake Superior) has been hypothesized to result from local adaptation to large-scale spatial variation in prevailing growth and mortality regimes (Wang, Höök, Ebener, Mohr, & Schneeberger, 2008).
The idea that environmental variation in growth potential interacts with size-selective fisheries to affect the evolution of maturation traits is not new, yet this interaction has not been the focus of modeling studies. Certainly, all models of fisheries-induced evolution of maturation traits include somatic growth rate as an input variable. Growth has been formulated in different ways in these models, but typically is not manipulated to generate predictions. In Atlantic cod, several modeling studies have considered the evolutionary and demographic consequences of size-selective fisheries and different gear types under a specified growth regime (Hutchings, 2009;Jørgensen et al., 2009;Zimmerman & Jørgensen, 2017). In ecogenetic models of fisheries-induced evolution of maturation traits in lake whitefish, somatic growth rate is assumed to be negatively density-dependent according to a specified functional form. This allows for feedbacks between ecological and evolutionary processes in response to alternative forms and intensities of size-selective fishing (Dunlop, Eikeset, & Stenseth, 2015;Wang & Höök, 2009). Eco-genetic modeling also has been used to explore how different intensities of density-dependent growth affect population dynamics and stock productivity (Gobin, Lester, Fox, & Dunlop, 2016). Intrinsic growth rate also has been treated as a decision variable in models of size-selective fisheries Eikeset et al., 2016).
To better understand how environmental variation in growth might influence the evolutionary response to size-selective fisheries, we focussed on two lake whitefish stocks in Lake Huron with long-term assessment data. Our first objective was to characterize and compare the stock-specific growth regimes and maturation schedules of lake whitefish. We reasoned that if stocks differ in their growth and maturation traits, they might also differ in their evolutionary response to size-selective fishing. Our second objective was to model the influence of fishing mortality rate, different gear types, and growth regimes on optimal maturation schedules of lake whitefish. This was done using a state-dependent optimization model of energy allocation that was informed by observational data. In our model, we vary the type of fishing gear, fishing mortality rate, and somatic growth rate.

| Lake whitefish data
Our first objective was to characterize cohort-specific maturation traits, prematuration growth rates, and their relationship in Lake Huron lake whitefish. We used long-term monitoring data from the Ontario Ministry of Natural Resources and Forestry's Offshore Index Assessment (fall) program (Speers & Cottrill, 2008 We pooled the two sites in the Main Basin (MB) because the distance between sites (~150 km) was small relative to the presumed scale of movement in these fishes. We assume the MB population is reproductively isolated from southern Georgian Bay, but the population structure and movement of lake whitefish among basins are not completely understood, and some degree of mixing cannot be ruled out (Eberts et al., 2017). The MB site is contained with quota management area (QMA) OH-4 and OH-5; the GB site is contained within QMA GB-4 (see map in Gobin et al., 2015). Lake whitefish commercial fisheries are managed primarily through quota management, but also through gear restrictions (e.g., Caroffino & Lenart, 2017).
Previous studies have characterized growth and maturation traits using these or subsets of these data, but were not in the form we needed. Wang et al. (2008) analyzed variation in maturation schedules among multiple sites throughout the Laurentian Great Lakes, but did not analyze trends in prematuration growth rates or show the relationship between growth and maturation (i.e., the plastic effects of growth rate on maturation timing). Gobin et al. (2015) showed trends in juvenile growth rate in MB, but not in GB, and also did not directly relate cohort-specific growth to maturation traits.

| Growth and maturation traits
Prereproductive growth rates and maturation traits were characterized separately for each site and cohort. This cohort-based approach generally follows that of Morgan and Colbourne (1999). All data manipulation and data analyses were performed in SAS/STAT software version 9.3 (SAS Institute Inc. 2011). To estimate prereproductive growth rate, we fit linear models to length-at-age data for ages 1-3.
Cohort-specific slopes and 95% confidence limits were retained if slopes were significant (≠0) at α = 0.01. For females only, logistic regression was used to estimate the age and length at 50% maturity (A 50 and L 50 ) using procedure PROBIT in SAS. Females were used to be consistent with our optimization model, and because they typically mature at an older age and larger size than males (Beauchamp et al., 2004;Wang et al., 2008). Cohort-specific estimates of A 50 and L 50 were retained if the intercepts and slopes were significant (≠0) at α = 0.01. Fiducial confidence limits were calculated for A 50 and L 50 using procedure PROBIT in SAS. Annual variation in prereproductive growth rate, A 50 , and L 50 were visualized using loess curves in package ggplot in R (R Core Team 2016).
We used general linear models to quantify the relationship between prematuration growth rate and maturation traits (age and length at 50% maturity). For each maturation trait, the model included cohort-specific growth rate, the factor site, and the site × growth rate interaction. The interaction was dropped before assessing main effects (α = 0.05). A general linear model was justified on the basis of a causal relationship between maturation traits (y) and prematuration growth rate (x) (Smith, 2009). Given error in x, however, we also estimated the slope correction factor as var(x)/ (var(x)-var(u)), where var(u) is the estimate of variance in x attributed to technique error and sampling variance (Smith, 2009). In our case, var(x) was the variance in prereproductive growth rate (i.e., the x values). We used the median of the variance in the estimate of slope from the length versus age relationships to represent var(u).

| Catch per unit effort
The presence of large fish has the potential to buffer populations from the undesirable evolutionary effects of fishing (Hixon, Johnson, & Sogard, 2014;Law, 2007). Thus, we extracted information on the abundance of large fish in each year. For this analysis, we used catch counts rather than fish size data, because not all fish were measured (a maximum of 50 fish are measured per deployment). Our index of abundance was catch per unit effort (CPUE) in the two largest mesh sizes (114 and 127 mm). One unit of effort was one 24 hr (approximately) multimesh gill net set, which always included one 50-m panel each of 114-and 127-mm gill net. These mesh sizes have maximum selectivities for 513 and 573 mm lake whitefish, respectively, and rarely catch fish <400 mm (Zhao & Morbey, 2017).
Following Gobin et al. (2015), we calculated the annual geometric mean of ln(CPUE + 1) and the associated 95% confidence limits. This also involved adjusting catch rates of multifilament nets (used in years ≤ 1994) by 1.8 (Collins, 1979). Annual variation in CPUE was visualized using loess curves in package ggplot in R.

| Probabilistic maturation reaction norms
For each site (GB and MB), we characterized the probabilistic maturation reaction norm (pMRN) to assess site differences. A univariate MRN describes how age or size at maturity varies with an environmental factor (e.g., growth rate) within a single genotype (Falconer & Mackay, 1996). A bivariate MRN, which is especially popular for fish, considers the plasticity of both size and age at maturity in relation to growth, and is typically shown as a line on a bivariate plot of size versus age. Variation in growth rates can cause variation in age and size at maturation along the MRN but does not cause the estimated MRN itself to change (Stearns, 1992). An extension of the MRN concept involves accounting for the inherent between-individual variability in the timing of maturation by describing it as probabilistic. A pMRN describes the probability that an immature individual of a given age and size will mature in the following growing season (Heino, Dieckmann, & Godø, 2002). The pMRN shows the age-specific maturation thresholds (i.e., the maturation rules) and can also be shown on a bivariate plot of size versus age.
We estimated pMRNs (L p50,a values) for females following the procedures described in full by Barot, Heino, O'Brien, and Dieckmann (2004). Sites were analyzed separated using all the data available.
In step (1), we calculated the age-and length-specific probabilities of maturity conditional on being alive (age-specific ogives, o[a, L F ]). Following Ficker, Mazzucco, Gassner, Wanzenböck, and Dieckmann (2014), we first compared different logistic regression models of the maturity ogives o(a, L F ) to choose among different combinations of the predictor variables age as a continuous variable (a), age as a factor (a F ), and fork length (L F ). We considered five models: (1) a + ε; (2) To obtain reliable estimates from logistic regression, we only considered ages where the number of mature fish was ≥10, the number of immature fish was ≥10, and the total number of fish was ≥100. We chose the model with the best performance (lowest Bayesian information criterion with all variance inflation factors <10).
In step (2), we used estimates of age-specific annual growth rate (ΔL F [a]) and the ogives from step (1) This ensured that the relationship between logit(m) and L F was linear (based on visual examination, at L F values associated with m > 0.5, logit[m] began to decelerate). To estimate confidence intervals around L p50,a , bootstrapped resampling of the observed data was performed with 1,000 replicates. Our approach to estimate pMRNs differed slightly from Wang et al. (2008) in that we estimated age-specific growth rates from fitted von Bertalanffy models, and, as described above, applied constraints on samples sizes and values of m(a, L F ) to give more reliable estimates.
To estimate pMRNs, all years were pooled because there were insufficient data for a cohort-based approach. Thus, annual variation or trends in growth were not taken into account in Equation 1. This introduces uncertainty in the comparison of pMRNs between sites. As a result, we also tested for site differences in the maturity ogives, o(a, L F ), while allowing for random intercepts among years.
This was done using logistic regression in the procedure GLIMMIX in SAS, which permits the specification of random effects in logistic regression. The model included age, length, site as a fixed factor, and year as a random factor.

| Dynamic optimization model: overview
We used a phenotypic approach to model the growth and maturation of lake whitefish as a state-dependent energy allocation problem using a dynamic optimization model. In effect, a phenotypic approach to optimization let us determine what energy allocation strategies might be beneficial under different selection regimes, without any consideration of genetic inheritance mechanisms (for a comparison of different modeling approaches, see Kokko, 2007). Dynamic optimization models have been used widely to study behavioral or physiological decisions that depend on current state and time (e.g., time of day, day of year, or year of life) given their consequences on cumulative survival and/or reproduction at a future time (McNamara, Houston, & Collins, 2001). The model we developed closely parallels  model of energy allocation in Atlantic cod. In this model, a female allocates incoming energy to growth, storage, and reproduction in such a way to maximize the number of eggs produced over a lifetime in a stochastic environment. Females also experience different forms of mortality. Variation in lifetime patterns of growth and maturation emerge as the outcome of state-dependent energy allocation decisions. This model has been used to explore the phenomenon of skipped spawning (Jørgensen, Ernande, Fiksen, & Dieckmann, 2006) and the fitness consequences of size-selective fisheries for a variety of life history traits . This type of model was wellsuited to explore the interaction between size-selective fisheries and environmental variation in growth potential. We note that this approach does not explicitly consider ecological feedbacks that could influence somatic growth during the process of adaptation, such as negative density dependence via intraspecific competition for food. The dynamic optimization model contains a bioenergetics submodel. This describes an individual's metabolism as well as its capacity to consume, store, and grow as functions of its current size and energy stores. These functions are fully explained in . For lake whitefish, the context was a simple lake environment: Water temperature and prey availability were specified for each month of the year to reflect observed seasonal patterns.
The dynamic optimization model determines how excess energy should be used. If a fish's energy balance is positive after finding some amount of food, then the extra energy may be allocated to structural growth in length or to energy stores for future use in respiration or spawning. Every spawning month (e.g., assumed to be November for lake whitefish), the fish can convert all or none of its stored energy to eggs. Every month the fish is subject to some probability of surviving until the next month, which is the added effects of natural mortality (assumed to be invariant with respect to size and age) and size-selective fishing mortality. Because of mortality, a trade-off emerges between growth and storage. Allocating energy to stores in preparation for spawning is safer, in that some offspring will be produced early in life. The riskier strategy of taking time to grow larger before maturing is more rewarding in the end, but only if the fish can survive to reproduce.
An individual's fitness is the total number of offspring it produces over its lifetime. The dynamic programming algorithm finds the optimal allocation of net energy intake to storage (u) versus growth that maximizes the total number of eggs produced. The algorithm starts at the hypothetical last month (t = T) of the fish's life when the optimal choice for allocation (u) can only be "store as much as possible for the last chance to spawn." Going backwards in time, 1 month at a time, expected fitness (V) is calculated for all the combinations of states the modeled fish could be in: for every starting fork length (L F = 10-100 cm, in 2 cm increments), every starting energy stores (E = 0%-100% of storage capacity, in 10% increments), every food intake possible that month (ϕ = 20%-100%, in 10% increments), every allocation decision (u = 0%-100% of extra intake energy converted to stores, with the remainder going to structural growth, in 10% increments), and in November every spawning decision (yes/ no): is a function of fecundity b, survival probability S to the next age, the conditional probability of finding food P and food intake ϕ in the next time step, and V in the next time step. The L F (t + 1) and E(t + 1) functions are determined by the bioenergetics equations.
Reproductive output via the spawning of eggs (b) can only occur in November. V is assigned as the maximum value over all possible values of u. The allocation decision (u) resulting in the highest expected fitness is stored in an optimal decision matrix. The optimal spawning decision is also stored.
Individuals are then allowed to grow and interact with the optimal decision matrices, assuming a variable environment. Environmental stochasticity is modeled as variability in food intake around an expected proportion χ of maximum consumption C max . Thus, the forward simulation (FS) produces a random series of autocorrelated values for food intake and individual "fish" respond to these circumstances. Forward simulation reveals the optimal growth and maturation schedule.

| Dynamic optimization model: parameterization
The dynamic optimization model was applied to Lake Huron lake whitefish and was programmed in C++ using QT Creator 3.1.1 (https://www.qt.io/ide/). This involved some modification to elements of the Atlantic cod model. First, in the bioenergetics submodel, we used the respiration and maximum consumption functions for bloater (Coregonus hoyi) from Rudstam, Bindowski, and Miller (1994) and applied to lake whitefish by Madenjian et al. (2006). Respiration (R) is a power function of body weight and also depends on temperature. In this function, we initially used a coefficient of 0.00138 as recommended by Madenjian, Pothoven, and Kao (2013). Monthspecific values of temperature (T, °C) were obtained from a study of lake whitefish in Lake Michigan . These values were 1.5, 0.5, 0.5, 2.0, 4.5, 6.0, 7.5, 8.0, 11.0, 9.5, 6.5, and 4.0°C, respectively, for January to December. Respiration was converted into J g −1 month −1 by multiplying R by 13,556 J g O 2 −1 and 30.5 day/ month.
We used the temperature-dependent consumption function (C max ) for Coregonus hoyi from Binkowski and Rudstam (1994). C max is a power function of body weight and also depends on temperature according to the function F(T) from Kitchell, Stewart, and Weininger (1977). Maximum consumption was converted into energy available for allocation by accounting for the multiplicative effects of egestion (25% of C max ), excretion (10%), and the energetic cost of processing assimilated energy (specific dynamic action [SDA] = 17%). To convert to J g −1 month −1 , available energy was multiplied by 30.5 day/month and month-specific values of prey quality (J/g prey) from Madenjian et al. (2006). These values were 2,716 J/g for December to March, 2,449 J/g for April to June, 2,838 J/g for July and August, and 2,984 J/g for September to November.
Second, we assumed size selectivity fishing mortality from trap nets and gill nets, the two major types of fishing gear used to catch lake whitefish in Lake Huron. Based on data from the Ontario Ministry of Natural Resources and Forestry's commercial catch sampling program (Gile & Milne, 2006), we looked at the use of these two gear types in our sites. We selected the set of samples (i.e., deploy- Size selectivity functions for gill nets and trap nets were recently estimated for Lake Huron lake whitefish (Zhao & Morbey, 2017).
Size-selective mortality from gill nets follows a double logistic function and depends on mesh size (m in cm): where S m (L F ) is selectivity as a function of L F , α 1 = 6.23, β 1 = 4.10, α 2 = 2.87, and β 2 = 5.09. C is a constant (=1.2398) such that maximum S m = 1. We initially set m to 11.4 cm, which currently dominates in the commercial gill net fisheries in Lake Huron (maximum selectivity = 51.3 cm). Size-selective mortality from trap nets follows a logistic function: where S t (L F ) is trap net selectivity as a function of L F , α = 0.32 and β = 47.1 cm. This is an asymptotic function, where selectivity reaches 0.99 at 61.5 cm. Trap nets can vary in their construction, but sizes and dimensions of trap net components were not systematically recorded in the commercial catch sampling program.
A few additional modifications were made. We assumed no additional energy costs associated with movement from summer feeding areas to spawning areas, and instead assumed an additional mortality cost associated with spawning (M S ). We also applied a consumption multiplier (p) to permit a reduction in C max for spawners in November. We also assumed that excess energy would be stored in muscle only, rather than in muscle and liver as in Atlantic cod.
Additional parameters were required for functions relating stored energy to weight (ε, K max , K min , L STD , ρ E ), the conversion of net energy intake to growth and stores (δ growth , δ store , and ρ S ), and the conversion of stored energy to eggs (κ 5 ;Appendix 1). For some of these, we were able to estimate or approximate values using data for Lake Huron lake whitefish. When data were not available, we kept the original values used by . We assumed a value of 0.2 year −1 for instantaneous natural mortality rate (M), which led to a large range of ages in our modeled population in the absence of fisheries and spawning mortality. Under constant age-specific mortality, M = 0.2 corresponds to a 50% reduction in cohort size by age 3.5 years and a 95% reduction in cohort size by age 15 years.
After generating the decision matrices for optimal allocation u and spawning, forward simulation was used to extract information about size at maturity. In each FS, we let individuals (n = 50) grow and interact with the decision matrices. All environmental and bioenergetics parameters were the same as those used in the backward iteration phase. The output of a FS included age (converted into years), L F (converted to mm), and spawning decision (to spawn or not). Simulated individuals followed different growth trajectories owing to environmental stochasticity, but matured at similar lengths regardless of their experienced growth rate. The transition from being immature to mature was discrete with spawning occurring nearly every year thereafter ( Figure 1). In effect, the optimal maturation reaction norm was flat and could be characterized by the mean length at maturation (L m ). For each simulation, we calculated mean L m for the 50 individuals (the standard deviation of L m was small, ranging from about 1-3 mm).

| Model calibration
Given inevitable uncertainty in the parameters values used in the DP, we adjusted some parameter values so that maturation occurred at realistic sizes for lake whitefish. Prior to calibration, the DP systematically underestimated size at maturity relative to observed L 50 values. To increase the benefit of delayed maturation, we first increased the costs of spawning (M s and p). These were parameters for which we had no a priori information, and costs of migration and spawning are important components of the Atlantic cod model . We also adjusted the coefficient of the respiration function, which seemed to have a large effect on the model output.

| Model sensitivity
We assessed the sensitivity of optimal L m to fish growth rate, fishing mortality rate, and gear type. To induce differences in growth, we varied the energetic cost of processing food (SDA) from 0.085 (50% of 0.17) to 0.34 (200% of 0.17) by increments of 0.0425. Based on simulated growth trajectories, increasing SDA simulated a reduction in age 1-3 growth rate from 76 mm/year (when SDA = 0.085) to 55 mm/year (when SDA = 0.34). Other parameters affecting energy intake (e.g., χ) are expected to affect F I G U R E 1 Output of the state-dependent energy allocation model for 500 simulated lake whitefish using baseline parameters. (a) Growth trajectories for the first 10 simulated fish, showing the seasonal pattern of length increment and the months of November when the optimal decision was to spawn (solid symbols). (b) The relationship between age at maturation and estimated age 1-3 growth rate (mm/year) for all simulated fish. (c) The relationship between length at maturation (mm) and estimated age 1-3 growth rate (mm/year) for all simulated fish. The linear regression lines are shown in panels (b) and (c) maturation schedules in a similar way as SDA. Instantaneous fishing mortality (F) for maximally selected fish was arbitrarily set to 0, 0.1, 0.2, or 0.3 year −1 to show a range of responses. Gear type was set to be gill net or trap net. For gill nets, size selectivity was a function of mesh size, and we set mesh size to 89 mm (3.5 inches; maximum selectivity at 403 mm) or 114 mm (4.5 inches; maximum selectivity at 513 mm; Zhao & Morbey, 2017). Trap nets of different variants have not been described or modeled for lake whitefish fisheries.

| Growth, maturation traits, and CPUE
Lake whitefish showed some differences in growth rate between the two Lake Huron sites. In both sites, prematuration growth rate was variable among cohorts and showed a declining trend across the time series (Figure 2a). In MB, the 1983 and 1984 cohorts had the highest growth rates; in GB, the 1981 and 1982 cohorts had the F I G U R E 2 Patterns in the growth, maturation traits, and catch per unit effort of lake whitefish in Georgian Bay (GB; solid symbols) and the Main Basin of Lake Huron (MB; open symbols). On the left, cohort-specific values are shown for (a) age 1-3 growth rate, (b) age at 50% maturity, and (c) length at 50% maturity. In (d), catch per unit effort in 114-mm and 127-mm gill nets is shown for each year. Also shown in (e) are the relationship between age at 50% maturity and growth rate and (f) the relationship between length at 50% maturity and growth rate. Loess curves are shown in panels (a), (b), (c), and (d), whereas the linear regression lines are shown in panels (e) and (f). The 95% confidence limits are shown in panels (a) and (d); the 95% fiducial limits are shown in panels (b) and (c) Cohort Growth (mm/year) Year ln (CPUE+1) Growth (mm/year) Cohort L 50 (mm) Age at 50% maturity (A 50 ) was variable among cohorts and showed an increasing trend in both sites (Figure 2b). Length at 50% maturity (L 50 ) was also variable among cohorts, but did not display any obvious long-term trend in either site (Figure 2c). Variability in

| Probabilistic maturation reaction norms
The estimated pMRNs suggest site differences in maturation thresholds. In the estimation of pMRNs, the best-performing statistical models of the maturity ogives o(a, L F ) used age and length as continuous variables, meaning that fish were more likely to be mature if they were older and larger for their age. L p50,a values could be estimated for ages 2-6 in GB and ages 3-6 in MB, and showed that pMRNs were higher in GB than in MB (Figure 3). The slopes of the pMRNs were slightly negative in both sites (Figure 3).
When testing for site differences in the maturity ogives, there were overall effects of age (β = 0.616, t 20,592 = 13.0, p < .0001), length (β = 0.022, t 20,592 = 26.5, p < .0001), and site (β GB vs. MB = −0.784, t 71 = −3.4, p = .001; model intercept = −12.196). Thus, for any age and length combination, GB females had a lower probability of being mature than MB females. This is consistent with larger size thresholds for maturation in GB than in MB. For example, at age 4, the model estimates a length at 50% maturity of 470.4 mm in GB and 435.4 mm in MB.

| Dynamic optimization model
Calibration of the DP was required so that optimal L m better approximated observed L 50 values. When we included additional costs associated with spawning (M S = 0.3 and p = .25) and a lower value for the respiration coefficient (0.001 instead of 0.00138), optimal L m was 476 ± 1.6 mm when F = 0 and SDA = 0.17. With additional costs associated with spawning, lake whitefish delay spawning and instead allocate resources to growth to achieve a higher body size. A lower respiration coefficient helped to increase growth rate. Other explanations for the precalibration mismatch between model and data, such as a minimum required size for sexual maturation, likely would produce a similar effect to that of higher M S and p.
Fishing mortality and fish growth rate both affected optimal L m .
When controlling for gear type and SDA, increasing fishing mortality from 0 to 0.3 year −1 invariably reduced optimal L m ( Figure 4). When controlling for gear type and F, lowering growth rate by increasing SDA invariably reduced optimal L m ( Figure 4). Higher SDA favors lake whitefish that mature at a lower size threshold because slower growth lowers the fecundity advantage associated with delaying maturity.
The magnitude of the effect of increasing fishing mortality depended on gear type and fish growth rate. For trap nets and 114 mm gill nets, increasing fishing mortality reduced optimal L m across all values of SDA, but to a greater degree at lower values of SDA (fast growth; Figure 4a,b). Optimal L m was always lower than maximum selectivity for these two gear types (513 mm for 114 mm gill nets; >615 mm for trap nets). When fishing with 89 mm gill nets, fishing mortality more strongly reduced optimal L m at high values of SDA (slow growth) and had less effect on optimal L m at low values of SDA (fast growth; Figure 4c). Buffering occurred at low values of SDA because fish could quickly outgrow the sizes when they were most vulnerable to gill nets and mature at a larger size. Optimal L m was higher than maximum selectivity (403 mm) for a wide range of low F and low SDA values.
F I G U R E 3 Probabilistic maturation reaction norms for female lake whitefish from Georgian Bay (GB; solid symbols) and the Main Basin of Lake Huron (MB; open symbols). At each age a, L p50,a is shown along with the 95% confidence limits from bootstrap resampling. MB lake whitefish had smaller threshold sizes for maturing than those from GB, and in both sites, the pMRN was negatively sloped

| D ISCUSS I ON
We used a dynamic optimization model to explore the long-term evolutionary effects of different growth and mortality regimes on maturation size thresholds in lake whitefish. When we simulated a low growth selection regime by increasing the energetic cost of assimilating food (specific dynamic action or SDA), the model predicted a smaller optimal size threshold for maturation. When we simulated size-selective fishing using realistic functions for gill nets and trap nets, an increase in fishing mortality also favored a smaller optimal size threshold for maturation. These general predictions were not unexpected and do not differ from those generated by deterministic life history models (e.g., Roff, 1984). Moreover, growth and mortality regimes account for the local adaptation of maturation traits among populations of lake whitefish (Beauchamp et al., 2004), and we expect both to contribute to contemporary adaptation of maturation schedules.
In the model, the long-term evolutionary effects of fishing with different gear types depended on the growth regime. When fishing with trap nets, the optimal response to higher fishing mortality was to mature at a smaller size, and this effect was magnified under a fast growth regime. With trap nets, smaller fish are less vulnerable to fishing owing to the sigmoid-shaped selectivity of this gear type (Zhao & Morbey, 2017). Fish that mature at smaller sizes have lower reproductive output, but they also avoid growing into vulnerable size classes where their survival would be low. In contrast to trap nets, gill nets present two size refuges for fish, and the presence of an upper size refuge has the potential to buffer a population from the long-term evolutionary effects of fishing Law, 2007). However, this was not an outcome for Lake Huron lake whitefish when fishing with 114-mm gill nets, the mesh size currently dominating in the commercial fishery. Buffering was a possible outcome with the smaller mesh gill net (89 mm), but only when growth potential was high. Under these conditions, fish could quickly outgrow the period of vulnerability to harvest and then mature at a large size with high fecundity. When growth potential is low, however, fishing with the 89-mm gill net favors maturity at smaller sizes than would occur with a larger mesh size.
The growth rates experienced in the two sites have implications for the evolutionary outcome of size-selective fisheries. Within both Lake Huron sites, growth rates were variable among cohorts and generally declined across the time series. Trends of declining growth have been reported at multiple sites in the Laurentian Great Lakes and have been attributed to changing food webs (Ebener, 2013;Fera, Rennie, & Dunlop, 2015;Gobin et al., 2015;Lumb, Johnson, Cook, & Hoyle, 2007;Rennie et al., 2009). Through plasticity, a reduction in growth was associated with an increase in age at 50% maturity and a decrease in length at 50% maturity. The estimated pMRNs had shallow negative slopes in both sites. Overall, most of the estimated size thresholds for maturation (and the cohort-specific F I G U R E 4 The predicted, long-term evolutionary effects of mortality and growth regimes based on a state-dependent, dynamic optimization model of energy allocation for Lake Huron lake whitefish for trap nets (a), 114-mm gill nets (b), and (c) 89-mm gill nets. Fishing mortality rate (F) is shown on the x-axis, and the energetic cost of assimilating food (specific dynamic action or SDA) is shown on the y-axis. The contours show predicted lengths at maturity (mm), which were estimated from locally weighted smoothing

Gill net 89 mm
We expect the two sites in Lake Huron to differ in their evolutionary sensitivity to size-selective fisheries. Overall, GB fish grew faster and matured at a younger age and larger size than MB fish. Whereas the slope of the pMRN was similar between sites, its elevation was somewhat higher in GB fish. A higher elevation pMRN in GB could be the result of local adaptation to higher growth rates in GB, although other selection factors could be important. According to the dynamic optimization model, the evolutionary effects of increasing F should differ depend on growth rate and gear type. In a population already adapted to low growth (e.g., MB), fisheries-induced selection as a consequence of fishing with 114 mm gill nets is expected to be weaker (the contour lines are farther apart with increasing F) than in a population already adapted to high growth (e.g., GB). This pattern reverses, however, if fishing with 89-mm gill nets. This suggests that a strategy of using smaller mesh gill nets to protect large, fecund fish might be more effective in GB than in MB.
Dynamic optimization models determine optimal phenotypes at evolutionary equilibrium given a defined selection environment and do not make specific predictions about short-term evolutionary rates. The predicted outcomes for maturation size thresholds would only be realized if selection differentials were consistently high over multiple generations, and if there was sufficient heritability of statebased energy allocation and spawning decisions. In Lake Huron, lake whitefish populations experience highly variable selection environments both in terms of growth rate (e.g., Gobin et al., 2015) and fishing mortality (e.g., Ebener et al., 2008 (positively related to fish size) and the probability of survival given an attack (negatively related to fish size; Bence et al., 2003). In addition, lake whitefish can be targeted by recreational fishing. Thus, evolutionary trajectories of maturation thresholds might proceed at highly variable rates in nature.
Certain factors have the potential to buffer or counteract evolutionary change in maturation schedules. For example, fisheriesinduced evolution can be buffered by negative density dependence in somatic growth rate that is the outcome of intraspecific competition for limiting resources Gobin et al., 2016).
On the one hand, Healey (1980) provides experimental evidence of fisheries-induced, negative density-dependent growth at the whole lake level in small lakes, and Lorenzen and Enberg (2002) discuss more generally the importance of considering negative densitydependent growth as a population regulation mechanism in fish populations. On the other hand, Gobin et al. (2016) show an overall positive relationship between lake whitefish biomass and juvenile growth in QMA OH4-5. They further suggest that site-specific changes in lake whitefish density may be driven more so by variation in habitat quality, which could determine growth rate and carrying capacity, than to variation in fishing mortality. Positive densitydependent processes may also play a role. For example, in a spatially and temporally heterogeneous environment, habitat selection by mobile predators could lead to positive correlations between habitat quality, somatic growth rate, and density. Further research about density-dependent behavior and movement would help to reduce uncertainty about density effects at the different spatial scales relevant to lake whitefish.
Selection also may target intrinsic growth rate or one of the bioenergetic determinants of prematuration growth rate (e.g., maximum consumption, C max ). For example, long-term fishing pressure lead to the evolution of slower growth rates in Coregonus palaea in a small lake, where fishing is regulated by a minimum size limit (Nusslé et al., 2009). In experimental studies of captive Menidia menidia, harvesting of large fish leads to the evolution of slower growth rates, whereas harvesting of small fish leads to the evolution of faster growth rates (Conover & Munch, 2002). Evolutionary models show that fishing mortality can favor higher or lower growth rate, depending on the size selectivity of fishing mortality and other factors Eikeset et al., 2016).
Concern about fisheries-induced evolution has pervaded discourse on fisheries management and poses important questions about whether (and how) fisheries-induced evolution should be controlled or reversed (e.g., Dunlop et al., 2015;Heino et al., 2013).
A frequent recommendation is to protect large fecund fish, which could be done, for example, by fishing using a gear type with bellshaped size selectivity (Hixon et al., 2014;Law, 2007). Small mesh gill nets in particular can favor delayed maturation, extend the age and size distribution of the population, and result in good fisheries yield in the long term Zimmerman & Jørgensen, 2017). However, fishing with small mesh gill nets may intensify evolution toward smaller size thresholds for maturation under poor growth or high mortality regimes. Whether or not this happens will depend on the particular fishery and the biology of the target species. In unpredictable and uncertain selection regimes, reliance on small mesh gill nets to favor delayed maturation might be a risky strategy . A more cautious option for protecting large, fecund fish would be to keep overall mortality rates low, either through further controls of fishing or sea lamprey predation. For example, Hutchings (2009) advocates for fishing below a threshold level of fishing mortality (F evol ) that would favor delayed over early maturity. These values of F evol can be calculated for different stocks with different growth rates, mortality regimes, and gear types.
Our model provides insight about the potential for gill nets to buffer lake whitefish populations from the long-term evolutionary effects of fishing. Compared to trap nets, gill nets are expected to provide an upper size refuge for fish. Fish that grow through vulnerable stages to mature at a larger size will be more fecund over many years and could potentially help to buffer populations from low recruitment. However, this is not a general rule. When considering the specific case of Lake Huron whitefish, we predicted very similar evolutionary effects of fishing with trap nets and 114-mm gill nets.
Buffering only occurred with the smaller mesh gill net and only in a high growth, low mortality environment.

ACK N OWLED G EM ENTS
This project was supported by the Great Lakes Fishery Commission

CO N FLI C T O F I NTE R E S T
None declared.