Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Spatial and length-based models for management of migratory transboundary species : application to Pacific… Wor Lima, Catarina 2018

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


24-ubc_2018_may_wor_lima_catarina.pdf [ 4MB ]
JSON: 24-1.0365575.json
JSON-LD: 24-1.0365575-ld.json
RDF/XML (Pretty): 24-1.0365575-rdf.xml
RDF/JSON: 24-1.0365575-rdf.json
Turtle: 24-1.0365575-turtle.txt
N-Triples: 24-1.0365575-rdf-ntriples.txt
Original Record: 24-1.0365575-source.json
Full Text

Full Text

Spatial and length-based models for management of migratorytransboundary species: application to Pacific hake (Merlucciusproductus)byCatarina Wor LimaB.S. Fisheries Engineering, Universidade Federal Rural de PernambucoM.Sc. Marine Science, Virginia Institute of Marine Science - College of William & MaryA THESIS SUBMITTED IN PARTIAL FULFILLMENTOF THE REQUIREMENTS FOR THE DEGREE OFDoctor of PhilosophyinTHE FACULTY OF GRADUATE AND POSTDOCTORAL STUDIES(Zoology)The University of British Columbia(Vancouver)April 2018c© Catarina Wor Lima, 2018AbstractThis dissertation develops new modeling tools to provide new scientific perspectives on migratorytransboundary fish populations. I particularly focus on two main issues: (1) the interaction betweenage/size based migratory movement, spatial availability, and fisheries exploitation rates, and (2)time-varying fisheries selectivity caused by size based migration and cohort targeting. I use Pacifichake (Merluccius productus) as a case study. Pacific hake occurs off the Pacific coast of the U.S.A.and Canada and is characterized by ontogenetic migratory movement (older fish migrate furthernorth), strong recruitment events, and time-varying selectivity due to targeting of strong cohorts. Inthis dissertation, I present two new modeling approaches, and explore the effects of spatial structureon management outcomes using a closed-loop evaluation. First, I use a Lagrangian approach todevelop a migration model that describes the Pacific hake dynamics including seasonal migrations,fisheries dynamics, and cohort targeting. Second, I introduce a new stock assessment method thatbypasses the requirement of estimating selectivity by using catch at length and growth parametersto produce estimates of exploitation rate at age. This method produces mixed results because oflow precision in selectivity estimates. Third, I evaluate the impacts of harvest control rules on theoutcomes experienced by Canada and the U.S.A while sharing the Pacific hake resource. I use themigration model described above in a closed-loop simulation to evaluate the long-term impact of61 harvest control rules. The results indicate that there are differences in performance of harvestcontrol rules between the two nations when maximizing potential long-term yield and log yields.This is a result of the reduced availability of the resource in Canadian waters as the overall harvestrate increases. Caps on allowable catch may help to avoid reduced availability issues. I believe thatthe results and conclusions presented in this dissertation can inform the future management andmodeling of Pacific hake. In addition, the methods presented here could be used for management ofother resources subject to time-varying selectivity and other transboundary stocks managed underagreements that do not consider spatial management explicitly.iiLay SummaryWhen fisheries resources are shared by two or more nations, tracking the spatial range of the fishand avoiding management actions that severely change this range becomes important. I developedmodels to aid in the management of shared fisheries using Pacific hake as an example. Pacific hakeis fished by the U.S.A. and Canada, and the extent of their annual south-north migrations dependson the age/size of the individual fish, with only larger and older fish present in Canadian waters.I present two new modelling tools: a model that imitates annual Pacific hake migrations; and amodel that provides consistent estimates of fish abundance regardless of changes in the locationof the fish or fishing vessels. I use these models to explore the impacts that fish movement canhave on fisheries management. The contributions presented here can inform future management ofmigratory and shared resources.iiiPrefaceThis thesis is part of a project entitled: “Developing management procedures robust to variability instock productivity arising through trophic interactions and persistent environmental changes” whichwas part of a larger collaboration initiative, the Canadian Fisheries Research Network (CFRN),funded by the Natural Sciences and Engineering Research Council of Canada (NSERC). The ob-jective of all projects under the CFRN was to investigate research questions that are relevant toCanadian fisheries while simultaneously fostering collaboration between members of industry, gov-ernment and academia. The questions examined in this thesis were designed in partnership withthe industry representatives of the Canadian Pacific hake fisheries industry: Brian Mose (Deep SeaTrawlers Association) and Bruce Turris (Canadian Groundfish Research and Conservation Society)and the representative from Fisheries and Oceans Canada and who was also a committee member:Nathan Taylor.A version of Chapter 2 has been published: Wor, C., McAllister, M.K., Martell, S.J.D., Taylor,N.G., and Walters, C.J. 2017. A Lagrangian approach to model movement of migratory species.Canadian Journal of Fisheries and Aquatic Sciences doi:10.1139/cjfas-2017-0093.I developed the model, performed analyses and wrote the paper. My co-authors helped withmodel development, analyses and edited drafts.A version of Chapter 3 has been submitted for publication:Wor, C., van Poorten, B.T., and Licandeo, R., C.J. Walters Stock Reduction Analysis usingcatch at length data: Length-SRA.I developed the model, performed analyses and wrote the paper. My co-authors helped withmodel development, analyses and edited drafts.ivContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiLay Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivContents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiiAcknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiDedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiii1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Management of transboundary stocks . . . . . . . . . . . . . . . . . . . . . . . . . 31.2 Size segregation and fishing selectivity . . . . . . . . . . . . . . . . . . . . . . . . 51.3 Pacific hake - fisheries ecology and management . . . . . . . . . . . . . . . . . . . 61.4 Objectives and dissertation structure . . . . . . . . . . . . . . . . . . . . . . . . . 92 A Lagrangian approach to model movement of migratory species . . . . . . . . . . . 112.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132.3 Movement model framework . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132.3.1 Application to Pacific hake and simulation-estimation procedure . . . . . . 212.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 252.4.1 Model dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 252.4.2 Simulation-estimation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 262.5 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30v3 Stock Reduction Analysis using catch at length data: Length-SRA . . . . . . . . . . 373.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 373.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 393.2.1 Stock reduction analysis with catch at length data - length-SRA . . . . . . 393.2.2 Simulation-evaluation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 443.2.3 Misspecification of growth parameters . . . . . . . . . . . . . . . . . . . . 463.2.4 Real data examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 463.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 473.3.1 Simulation-evaluation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 473.3.2 Misspecification of growth parameters . . . . . . . . . . . . . . . . . . . . 513.3.3 Real data examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 533.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 544 Evaluation of harvest control rules for transboundary stocks . . . . . . . . . . . . . 604.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 604.2 Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 634.2.1 The simulation model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 644.2.2 The 40:10 harvest control rule . . . . . . . . . . . . . . . . . . . . . . . . 694.2.3 Linear harvest control-rules . . . . . . . . . . . . . . . . . . . . . . . . . . 694.2.4 Performance metrics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 704.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 734.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 825 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 875.1 Research summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 885.2 Future research directions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94viList of TablesTable 2.1 Indexes and variable definitions . . . . . . . . . . . . . . . . . . . . . . . . . . 15Table 2.2 Variable definitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16Table 2.3 Population dynamics model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17Table 2.4 Lagrangian movement model . . . . . . . . . . . . . . . . . . . . . . . . . . . 19Table 2.5 Fisheries dynamics model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21Table 2.6 Pacific hake model dimensions and parameter values . . . . . . . . . . . . . . . 23Table 2.7 Simulation-estimation scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . 25Table 3.1 Indexes, variable definition, and values used in simulation-evaluation . . . . . . 41Table 3.2 Variable definition for operating model, MSY quantities, and values used insimulation-evaluation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42Table 3.3 population dynamics for Length-SRA and operating model . . . . . . . . . . . 43Table 3.4 Management quantities and operating model . . . . . . . . . . . . . . . . . . . 44Table 3.5 Likelihood functions and penalties . . . . . . . . . . . . . . . . . . . . . . . . 45Table 3.6 Simulation-estimation scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . 46Table 3.7 Scenarios for testing misspecification of L∞ . . . . . . . . . . . . . . . . . . . . 46Table 4.1 List of scenarios for closed loop simulations . . . . . . . . . . . . . . . . . . . 67Table 4.2 Pacific hake operating model dimensions and parameter values . . . . . . . . . 68Table 4.3 Equations used to calculate performance metrics for evaluation of harvest con-trol rules. All quantities were averaged across simulation runs. . . . . . . . . . . 72viiList of FiguresFigure 2.1 Diagram illustrating the sine function used to describe the cyclic movementdynamics used in the Lagrangian model. The position (y-axis) of an individualchanges in cyclic waves as time (x-axis) progresses. In this diagram, we assumethat a cycle lasts 12 time steps and that movement occurs between values of 30and 60 along the spatial scale. . . . . . . . . . . . . . . . . . . . . . . . . . . 18Figure 2.2 Map of the study area. The dashed lines indicate the division between fishinggrounds. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22Figure 2.3 Illustration of the distribution of fish of age 5, in the month of July for thesingle and multiple groups versions in the absence of fishing. The multiplegroups version output is shown for all groups individually and combined. Thevertical dashed line represents the border between territories. . . . . . . . . . . 27Figure 2.4 Monthly representation of differences between spatial distribution of biomassfor single and multiple groups versions. Vertical bars indicate relative fishingeffort. Continuous and dashed lines represent vulnerable biomass for fish at age1 and 5 throughout one year migration cycle, one panel for month). Differentshades of gray indicate multiple and single group versions. The vertical dashedline represents the border between territories. Both fish biomass and effort werere-scaled from 0 to 1 for plotting. . . . . . . . . . . . . . . . . . . . . . . . . . 28Figure 2.5 % Relative error for estimation of key parameters of the movement model. Sce-narios are indicate in the upper right corner of each plot and the scenarios com-position is indicated in the plot titles. SG = single group version , MG= multiplegroups version and FG = fishing grounds. . . . . . . . . . . . . . . . . . . . . 29Figure 2.6 Median biomass distribution for simulation and estimation models for eachmonth of the year for scenario 12 - multiple groups, 5 areas, τ=2.0. . . . . . . . 31viiiFigure 3.1 Relative proportional error for main parameters for all scenarios considered inthe simulation-evaluation. Boxplots center lines indicate the median estimate.Lower and upper hinges indicate first and third quartiles. Upper and lowerwhiskers are given by the maximum and minimum values within the intervalsgiven by the hinge value +/- 1.5 · inter-quartile range (distance between the firstand third quartiles). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48Figure 3.2 Relative proportional error for main parameters for all scenarios considered inthe simulation-evaluation. Boxplots center lines indicate the median estimate.Lower and upper hinges indicate first and third quartiles. Upper and lowerwhiskers are given by the maximum and minimum values within the intervalsgiven by the hinge value +/- 1.5 · inter-quartile range (distance between the firstand third quartiles). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49Figure 3.3 Simulated and realized selectivity estimates for a set of years within simulation-evaluation time series. The estimated solid lines indicate median, 2.5% and97.5% quantiles for the derived selectivities. . . . . . . . . . . . . . . . . . . . 50Figure 3.4 Simulated and realized exploitation rate at length Ul,t when L∞ is misspeci-fied. Results shown for the last four years of simulation-evaluation time series.Boxplots center lines indicate the median estimate. Lower and upper hingesindicate first and third quartiles. Upper and lower whiskers are given by themaximum and minimum values within the intervals given by the hinge value+/- 1.5 · inter-quartile range (distance between the first and third quartiles). . . 52Figure 3.5 Fit to index of abundance, historical catches, and Yieldtarget and Utarget es-timates for Pacific hake and jack mackerel. Observed indexes of abundanceare shown in open circles, closed dots in Yieldtarget and Utarget panels indicatemodel estimates. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54Figure 3.6 Realized selectivity at age patterns across years for Pacific hake and jack mack-erel. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55Figure 4.1 Logistic functions used in the three movement scenarios considered in this study. 66Figure 4.2 Illustration of the relationship between relative spawning biomass and TAC fortwo example linear harvest control rules and the 40:10 harvest control rule usedin this study. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71Figure 4.3 Historical and median and 95% intervals for projected catches for three exam-ple harvest control rules under the "no strong recruitment" scenario. Harvestcontrol rules include the 40:10 rule with cap, and two linear harvest controlrules with biomass threshold (SBt/SBo) of 0.1 and exploitation rates of 0.1 and0.5. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74ixFigure 4.4 Comparison of relative performance (normalized to the maximum value by na-tion) of average log yield (log utility) and average yield between nations (trade-offs). Quantities were normalized by the maximum observation per nation andscenario for comparison. Nation 1 = U.S.A and Nation 2 = Canada. Colorsindicate harvest rates for linear harvest control rules and point shape indicatebiomass threshold relative to unfished spawning biomass (SBt/SBo). Text indi-cates performance of the 40:10 harvest control rule. . . . . . . . . . . . . . . . 76Figure 4.5 Median Average Annual Variability in yield (AAV). Surfaces indicate perfor-mance for linear harvest control rules. White circle indicates performance levelcomparable to the 40:10 harvest control rule. Minimum value along the surfaceis indicated with black square. Recruitment scenarios are indicated on the rows.Nation 1 = U.S.A and Nation 2 = Canada.Total is aggregate result. . . . . . . . 79Figure 4.6 % closures for total fisheries over the 60 years of simulation. Surfaces indi-cate performance for linear harvest control rules. White circle indicates per-formance level comparable to the 40:10 harvest control rule. Minimum valuealong the surface is indicated with a black square. Recruitment scenarios are in-dicated on the columns: A - no strong recruitments, B - one strong recruitmentper decade,C - two strong recruitments per decade. . . . . . . . . . . . . . . . 80Figure 4.7 % of time that biomass is above 40% of average unfished levels. Surfaces in-dicate performance for linear harvest control rules. White circle indicates per-formance level comparable to the 40:10 harvest control rule. Minimum valuealong the surface is indicated with a black square. Recruitment scenarios are in-dicated on the columns: A - no strong recruitments, B - one strong recruitmentper decade,C - two strong recruitments per decade. . . . . . . . . . . . . . . . 80Figure 4.8 Comparison of relative performance of average log yield (log utility) and aver-age yield between nations (trade-offs) for the three alternative movement sce-narios. Quantities were normalized by the maximum observation per nationand scenario for comparison. Nation 1 = U.S.A and Nation 2 = Canada. Colorsindicate harvest rates for linear harvest control rules and point shape indicatebiomass threshold relative to unfished spawning biomass (SBt/SBo). Text indi-cates performance of the 40:10 harvest control rule. . . . . . . . . . . . . . . . 81xAcknowledgmentsI would like to thank all my committee members: Steve Martell, Murdoch McAllister, Carl Walters,Nathan Taylor, and Villy Christensen for all the support throughout the development of this thesis.In particular I would like to thank Carl Walters for his patience in teaching me new concepts, hismentoring and guidance. I also owe my gratitude to Brett van Poorten who provided excellentadvice and mentoring in the development of Chapter 3. Additional thanks go to Steve Martell forthe internships at the International Pacific Halibut Commission and all the time spent teaching mehow to (really) use a computer and the secrets behind catch at age models.I would like to thank the Natural Sciences and Engineering Research Council of Canada (NSERC)for funding the Canadian Fisheries Research Network (CFRN) that included the research projectunder which this thesis was developed. My deep gratitude goes to all members of the CanadianFisheries Research Network. Being part of the CFRN was an immense privilege and I an amaz-ing learning experience on collaborative work and interdisciplinarity. The CFRN group showedme the value of collaboration in fisheries management and science. Special thanks goes to thegroup of students that worked with me on the Fisheries Evaluation Framework: Eric Angel, AllanDebertin, Danielle Edwards, Aaron Greenberg, Andrea Haas, Mike Hawkshaw, Sarah Hawkshaw,Robin Messenger, Dan Mombourquette and Courtenay Parlee. I hope to collaborate more with youin the years to come.My sincere thanks also goes to the group of students in my lab: Shannon Obradovich, RobertoLicandeo, Aaron Greenberg, Mike Hawkshaw, Sarah Hawkshaw, Danielle Edwards, Ben Nelson,Rachel Chudnow, Rachel Neuenhoff and Rachael Louton. You provided me with a huge amountxiof support that was essential in helping me complete this dissertation. Thank you for reading myearly/dreadful drafts, for the brainstorming sessions, for the company, for the kind words wheneverI needed them. You are all smart and amazing people.I would like to thank my family. My parents, Risolete Wor and Uka Lima, for providing all theemotional and financial support so that I could pursue my education and also for always believingin me and in my dreams. I also would like to thank Carole, Martin, Julian, Victoria, Daryl andKadyn: thank you for welcoming me into your family and for making me feel at home in Canada.Thank you to Dalek, for getting me away from my computer and out for walks. And last but notleast, I would like to thank my partner, Dominic, for all the patience, care, financial support andunconditional love.xiiTo Dominic, mainha and papai.xiiiChapter 1IntroductionFish species have evolved to optimize reproductive success, feeding rates and to avoid unfavourableconditions, such as high predation and cannibalism (McKeown, 1984). For many species this op-timization is achieved with migratory movement away from nursery grounds (for new recruits),and subsequently between feeding and spawning grounds (Harden Jones, 1968; McKeown, 1984).Migratory movement is indeed a common feature in the ecology of fishes, and has been a frequentsubject of study over the last centuries (see Morais and Daverat (2016) for a historic review of fishmigration literature). For iteroparous species, the migration routes are frequently traced betweenspawning and feeding grounds, and can occur repeatedly, in a cyclic manner, throughout an indi-vidual’s life (Secor, 2015; Walters and Martell, 2004). Examples of cyclic migration behavior arefound in tunas and tuna-like species (Nakamura, 1969; Merten et al., 2016), sharks (Hammerschlaget al., 2011), sardines (Lo et al., 2011), flatfishes (Hunter et al., 2003) and gadoids (Robichaud andRose, 2004; Alheit and Pitcher, 1995).Migration is frequently associated with segregation of subgroups within a population (Harden Jones,1968). For example, segregation can occur by sex (e.g., Okamura et al., 2014; Mucientes et al.,2009), size or age (e.g., Ressler et al., 2007; Chen et al., 2005). Size segregation is often a result offish migrating over larger distances; fish tend to migrate further and faster as they get larger leadingto segregation of individuals by size/age (Secor, 2015). Fish movement is of particular concern1when the migratory path of the exploited resource crosses political boundaries, so that if the fishstock is exploited, management may require transboundary management agreements (Munro, 2005;Sumaila, 2013). Management of transboundary stocks frequently assumes that the distribution ofthe resource is stable over time (Liu and Heino, 2013). Whenever changes in distribution occur,management agreements can become unstable and require renegotiations (Bjørndal and Ekerhovd,2014). One notable example of how changes in fish distribution can impact international manage-ment agreements is the case of Atlantic mackerel. The stock was initially managed through aninternational agreement between Norway, the EU and the Faroe Islands. However, around 2008 thestock started migrating into Icelandic waters which led to the development of a fishery for Atlanticmackerel in that country. The increase in Icelandic catches lead to conflict and consequent failureof the management agreement previously held (Hannesson, 2013; Spijkers and Boonstra, 2017).Size segregation can further the problem of transboundary management for migratory species.This is because there is a relationship between size composition and fishing mortality in a pop-ulation. If we assume that recruitment and natural mortality remain constant, populations underhigher fishing pressure will have truncated age structures, i.e., much lower proportions of olderand larger fish (Beverton and Holt, 1957). For this reason, management quantities, such as targetand threshold exploitation rates, will have an effect on the size distribution of the population. Thisimpacts the migration patterns of the exploited population that segregate by size/age; if larger/olderfish migrate further, the migration extent will be shortened by a truncated size/age structure. Thechange in migration extent may affect the nations sharing the resource to different degrees, andcould significantly diminish the availability of the resource to one or more of the fishing nations.In addition to affecting availability, migratory movement and size/age segregation can also leadto cohort targeting behavior by the fishing fleets. Cohort targeting occurs when fish segregate byage and the fishing fleet targets areas of high abundance. This leads to a disproportionate amountof effort being applied to the most abundant cohorts, either because the cohort is the outcomeof a strong recruitment event or because the cohort has been recently recruited and has not yetbeen fished. When cohorts are targeted, the fishery selectivity changes as the cohorts grow older,2resulting in temporal changes in fisheries selectivity. This is an important issue for many of themodern stock assessment methods that rely on estimates of selectivity to produce age specificfishing mortality rates and derived management quantities.Both migratory movement and size segregation are common for many exploited species aroundthe globe. For many species, size segregation comes as byproduct of the migratory movement withlarger individuals migrating further away from spawning grounds. A few examples of importantcommercial stocks that exhibit this behavior are: tunas (Nakamura, 1969), sharks (Camhi et al.,2009), Lake Erie walleye (Berger et al., 2012), Pacific halibut (Webster et al., 2013), sardines (Loet al., 2011) and Pacific hake (Bailey et al., 1982).This dissertation aims at identifying robust management procedures for transboundary fish pop-ulations that are subject to size/age segregation associated with migratory movement. Here, I de-fine the term “management procedure” as the series of choices made regarding the managementof exploited fisheries resources, including choices about data collection methods, stock assessmentmodels, target harvest rates and harvest control rules. I evaluate two main aspects that relate to themanagement of transboundary migratory populations: (1) changes in availability of a fisheries re-source to a nation due to management outcomes and (2) the issue of time-varying selectivity causedby changes in availability and cohort targeting by fishing fleets. I use the offshore stock of Pacifichake (Merluccius productus) as a case study, and focus the research questions on issues that arerelevant for Pacific hake management.1.1 Management of transboundary stocksMany fisheries around the world exploit stocks that are distributed across waters under the jurisdic-tion of two or more coastal nations, such stocks are called transboundary stocks (Sumaila, 1999;Miller and Munro, 2004). Several complications arise in the management of transboundary stocksbecause the optimal harvest strategy from the perspective of an individual nation usually differsfrom that of a group of nations exploiting a shared fisheries resource (Criddle and Strong, 2014).3Several studies have modeled the management of transboundary stocks under a game theoryapproach (Bailey et al., 2010). In this approach, nations involved in transboundary fisheries arerepresented by players who seek to maximize their benefits (utility) (Munro, 1979). The playerscan adopt cooperative or non-cooperative strategies with examples of both strategies found in fish-eries management (Munro, 2005). However, multiple studies show that although the individual gaincan be higher under a non-cooperative strategy, the overall utility is greater when all players opt fora cooperative strategy (e.g., Bailey et al., 2013; Ishimura et al., 2013). The adoption of cooperativestrategies is also frequently associated with better long-term management of the resource and in-creased fishery sustainability (Ishimura et al., 2013). Cooperative strategies are usually representedby treaties and agreements signed by two or more countries or by Regional Fisheries ManagementOrganizations. These agreements typically assume that stock distribution remains somewhat con-stant over time and that catch shares for each nation are proportional to the biomass of the stockoccupying that nation’s Exclusive Economic Zone (EEZ). This assumption does not always holdas changes in stock distribution often occur and lead to inadequacy of the management agreements(Liu and Heino, 2013).Variability in the distribution of exploited resources is not uncommon and can be caused bya wide variety of factors, including range contraction, i.e., as a population decreases, individualsconcentrate in optimal areas (e.g., Brodie et al., 1998), or changes in habitat characteristics suchas temperature or food availability (e.g., Rodríguez-Sánchez et al., 2002). The disparity betweenactual stock distributions and that assumed by management agreements can lead to ineffective man-agement regulations and sub-optimal outcomes for one or all nations sharing the resource (Millerand Munro, 2004; Bjørndal and Ekerhovd, 2014).Movement and migration patterns can both affect management outcomes and be affected bychanges in management procedures. Therefore, potential shifts in distribution need to be consid-ered when evaluating the performance of management procedures for transboundary stocks. How-ever, very few studies to date have included a spatial dimension when evaluating the effectivenessof fisheries management procedures. This effect is partly because spatial models are rarely used4for fisheries management (Berger et al., 2017a). However, this scenario has been changing in re-cent years due to increased spatial data availability, increases in computer power, and developmentof new data collection and modeling techniques (Goethel et al., 2011, 2016). These factors havecontributed to a better understanding of movement dynamics and spatial structures of fished pop-ulations. The improved understanding has also led to a surge in management related questionsassociated with movement and spatial distribution of exploited populations (Goethel and Berger,2017).1.2 Size segregation and fishing selectivitySelectivity is defined by a combination of two processes: vulnerability and availability. The vul-nerability process represents the contact selectivity, i.e., the proportion of individuals at a given ageor size retained by a fishing gear. The availability process represents the spatial dimension of thefisheries, i.e., the degree of overlap between the spatial distribution of the population being fishedand the spatial distribution of the fishing activity (Lee et al., 2017).Both vulnerability and availability of an exploited population can vary over time and thereforemodify the resulting selectivity. Vulnerability changes are generally associated with changes in thefishing fleet and gear type. Such changes can be identified if there is good knowledge of the historyof the fishery. Documentation on technological innovations, changes in target species or popula-tion groups, and development of new fleets are often found in stock assessment documents, legaldocuments, and industry reports. Availability changes as fish and fishers move. Fish movementcan be driven by migration and dispersal, and although ecological studies of fish movement areabundant, the inclusion of spatial distribution of fish in stock assessment models is comparativelyscarcer (Berger et al., 2017a). For this reason, in many cases, the effects of fish movement andspatial structure is modeled in stock assessment by considering spatial areas as different fishingfleets that are subject to different selectivities. Although this methodology can capture some of thevariance caused by migration and spatial structure, it can sometimes also lead to biased assessmentresults (Hurtado-Ferro et al., 2014).5The spatial distribution of fishing fleets is also non-static and affects availability. Many fishingfleets target a given size range, either for economic or regulatory reasons. For this reason, a fleettends to optimize its distribution in order to target areas of high abundance of the resource withintheir preferred size. For fisheries where the preferred size range is broad, fleets tend to target strongcohorts, especially if the fished resource tends to segregate by size or age. The cohort targetingbehavior will lead to changes in the fisheries selectivity over time, which in turn can be a difficultcharacteristic to model in stock assessments. Cohort targeting is a characteristic of the Pacific hakefisheries, which is consequently, a strong indication that the selectivity patterns for that fisherychanges as strong cohorts grow.Time varying selectivity is an important issue for current fisheries stock assessment (Gud-mundsson and Gunnlaugsson, 2012). In statistical catch at age models, now a widespread as-sessment method for many data rich fisheries (Methot and Wetzel, 2013), selectivity is used as amultiplier to fishing mortality, and therefore acts as a scaler to fishing mortality at age and, conse-quently, for fishing mortality of each cohort. If not appropriately considered in stock assessments,changes in selectivity over time can severely impact the estimates of fishing mortality and derivedmanagement quantities.In 2014, a special issue of the Fisheries Research journal was dedicated to the estimation ofselectivity, and the implications for stock assessment and fisheries management (Maunder et al.,2014). This special issue contains articles that evaluate the impacts of misrepresenting selectivityin assessments (e.g., Butterworth et al., 2014; Martell and Stewart, 2014), present methods fortreating selectivity estimation (e.g., Nielsen and Berg, 2014; Waterhouse et al., 2014) and reviewthe occurrence of various selectivity patterns in real data (e.g., Sampson, 2014).1.3 Pacific hake - fisheries ecology and managementPacific hake, also known as Pacific whiting, inhabits the waters off the west coast of North America,living within the California Current system (Lloris et al., 2005). Several stocks are encounteredthroughout the species distribution range (Chittaro et al., 2013; King et al., 2012), but a single and6offshore transboundary stock supports most of the Pacific hake commercial fisheries off the westcoast of the U.S.A. and Canada. This dissertation focuses on the offshore transboundary stock, andall references to Pacific hake are, therefore, referring to that stock.The recruitment of Pacific hake is highly variable, and very strong recruitment events happenwith some regularity (once or twice per decade) although not in a predictable pattern. The popu-lation age structure is heavily influenced by strong recruitment events, and the fishing fleets fromboth U.S.A. and Canada tend to target the strong cohorts. Maximum age is around twenty yearsand natural mortality is estimated to be around 0.23 year−1, but individuals older than fifteen yearsold are rarely encountered in the fisheries (Methot and Dorn, 1995). The offshore Pacific hakepopulation exhibits seasonal migratory behavior. Spawning occurs in offshore waters off southernCalifornia during the winter with fish migrating north between spring and fall to feed (Bailey et al.,1982; Ressler et al., 2007). The extent of the migrations is correlated with individual size. Largerfish, typically older than age four, migrate longer distances and are found to be more abundant thanyounger age classes in Canadian waters. Fish three years of age and younger tend to remain in U.S.waters off the coast of California and Oregon (Methot and Dorn, 1995; Ressler et al., 2007).The fishery occurs between May and November off the coast of northern California, Oregon,Washington (U.S.A.) and British Columbia (Canada), and is conducted almost exclusively by mid-water trawls. The fishing vessels operate in areas with bottom depth ranging between 100 to 500m.Until the 1990s, the Pacific hake fishery was strongly dominated by foreign fleets. The large scalefishery was started in 1966 with factory trawlers from the Soviet Union (Forrester et al., 1978)and expanded in the mid-1970s when factory trawlers from Poland, Federal Republic of Germany,the German Democratic Republic and Bulgaria joined the fishery. A joint venture fishery betweenU.S.A. trawlers and Soviet factory ships acting as motherships started in 1978. The U.S.A. nationalfishery expanded during the 1980s, and the fleet became entirely domestic by the early 1990s(Methot and Dorn, 1995; Ressler et al., 2007). In Canada, the domestic fleet also expanded afterthe mid-1980s but joint venture initiatives continued to occur until 2011. In recent years, most7Pacific hake production has been processed into headed and gutted products and surimi (Ressleret al., 2007; Nelson, 2014).The offshore Pacific hake stock is managed through a bilateral agreement between Canada andthe United States, known as the Pacific hake Treaty (United States State Department, 2004). Thistreaty and was initially ratified by the U.S. in 2006 but an error in the original text delayed its im-plementation until a new ratification in 2010. The treaty has been considered in force in Canadasince 2008. The agreement determines that a stock assessment should be conducted every year bya Joint Technical Committee (JTC). The JTC is composed of scientists appointed by each countryand independent members chosen by a private sector advisory panel. The stock assessment resultsare used in association with a defined harvest control rule (40:10 rule) to generate a recommen-dation for a coast-wide Total Allowable Catch (TAC), which is evaluated and adjusted by a JointManagement Committee (JMC). The JMC approved TAC is then shared between the U.S.A. andCanada following a fixed proportion: 73.88% and 26.12% of the TAC is allocated to the U.S.A.and Canada, respectively (United States State Department, 2004).Throughout the Pacific hake exploitation history there has been some concern about potentialchanges in stock distribution and stock biomass associated with gaps in the occurrence of strongrecruitment events. Whenever there is a strong cohort hiatus, as occurred between 1999 and 2010,the population age structure tends to become more truncated with most of its biomass being repre-sented by fish of younger age classes. Because the Pacific hake migratory movement is associatedwith fish size/age, the truncation of the population age composition has an effect on the populationrange and distribution (Hicks et al., 2016). This variability in stock distribution can particularlyimpact the northern distribution of the stock, resulting in lower availability of Pacific hake off thecoast of Washington and Canada. These changes have generated concern to the Canadian Pacifichake Advisory Panel regarding the efficiency of the current harvest control rule, and the consequentappropriateness of total allowable catch recommendations (Canadian Advisory Panel, 2013).A few studies have evaluated the performance of current and alternative management proce-dures considered for Pacific hake. Ishimura et al. (2005) evaluated the performance of the 40:108harvest control rule and a series of linear harvest control rules for the aggregate Pacific hake stock.They considered performance metrics relating to average yield (magnitude and standard deviation),probability of closures, and total biomass and found that the best results were obtained for low har-vest rates and low threshold biomass. They also found that the 40:10 harvest control rule performedwell in relation to the best linear harvest control rules. Punt et al. (2008) assessed the performanceof a set of harvest control rules a for groundfish stocks off the U.S.A. west coast, including Pacifichake. They evaluated the performance of such rules in relation to average catch and conservationmetrics. They found that threshold harvest control rules, i.e., rules that progressively adjust harvestrates as the biomass decreases (reach thresholds), tend to perform better in relation to conserva-tion objectives, but results are sensitive to assumed productivity and recruitment variability. Anadditional evaluation has been carried out by the Pacific hake JTC since 2014 (Taylor et al., 2014).However none of the evaluations listed above considered issues involving the spatial structure ofthe stock.1.4 Objectives and dissertation structureThis dissertation focuses on questions that are relevant to the management of transboundary stocks.Particularly, those stocks that are subject to migratory movement, ontogenetic segregation, andtime-varying selectivity. We use the Pacific hake resource as an example to illustrate the use ofthe modeling tools developed. The dissertation is organized into three modeling chapters and aconclusion. The three main chapters evaluate the following questions: (1, Chapter 2) How to modelcyclic ontogenetic migrations and fleet dynamics associated with cohort targeting? (2, Chapter3) Is it possible to overcome the requirement of estimating selectivity curves by using catch atlength and growth curve information? (3, Chapter 4) How do commonly considered harvest controlrules perform for the management of a transboundary fish stock when performance measures areevaluated separately for each fishing nation?In Chapter 2, I develop and implement a continuous migration model to describe the populationdynamics of the Pacific hake resource. The model simulates the cyclic migratory movement that is9characteristic of Pacific hake including cohort segregation. In this model, the population dynamicsof the resource are described with an age-structured model that allows for high variability in re-cruitment and strong recruitment events. This model is also coupled with a fleet dynamics modelthat simulates the tendency of fishing fleets to target areas of high abundance, which leads to thetargeting of strong cohorts and consequent time varying selectivity. Using a simulation-evaluationanalysis, I demonstrate how the model’s movement parameters can be estimated given commonlyavailable spatial catch at age composition data.In Chapter 3, I present a novel length based stock reduction analysis (Length-SRA) approach.This assessment model bypasses the requirement of estimating selectivity parameters, generatingno constraints on the occurrence of time-varying selectivity, similarly to what is obtained with aVirtual Population Analysis (VPA). I demonstrate the model performance under three exploita-tion trajectories, and with both time varying and time invariant selectivity patterns. I also use theapproach to estimate trends in selectivity and management quantities for two real data examples:Peruvian jack mackerel and Pacific hake.In Chapter 4, I use the movement model developed in Chapter 2 as an operating model in aclosed loop simulation to evaluate the performance of a series of linear harvest control rules and theharvest control rule currently used to manage the Pacific hake resource, the 40:10 harvest controlrule with a cap on maximum allowable catch. I evaluate five performance metrics that representmy interpretation of three potential fisheries objectives: high catches with limited variability, lowprobability of closures and maintenance of biomass above a threshold.The research presented in this dissertation offers opportunities for improving the understandingof the effects of spatial dynamics on fisheries assessment and management. The modeling toolsand framework presented here is informative for the future management of Pacific hake and othertransboundary migratory stocks.10Chapter 2A Lagrangian approach to model movementof migratory species2.1 IntroductionTwo main drivers of cyclic migration in marine fish species are seasonal availability of food andspawning behavior (Walters and Martell, 2004). In many iteroparous species, these drivers areresponsible for a continuous migration cycle between feeding and spawning areas (e.g., Hunteret al., 2003; Costa et al., 2012; Merten et al., 2016). Migration modeling is an important componentof fisheries science due to the common presence of migratory behavior in exploited species. In thisstudy, I present a novel method for modeling the migration of fisheries resources and the associatedfisheries dynamics that can arise from fish movement.The migration cycle observed in many iteroparous species can vary in extent (i.e., distance ortiming) for subgroups within a population, e.g., for different age or size groups (Ressler et al.,2007), sex (Okamura et al., 2014) or sub stocks (Carlson et al., 2014). This variability can leadto spatial segregation of subgroups within a fish population. When spatial segregation is present,subgroups are susceptible to distinct environmental and ecological drivers leading to differences innatural mortality, recruitment and to additional variability in spatial distribution (Ciannelli et al.,2008). In addition, spatial segregation within a population can cause spatial differences in vulner-11ability and fishing effort or mortality that are independent from fishing gear or fishing techniques(Martell and Stewart, 2014). If ignored, migratory movement and spatial segregation can lead tostrong bias in surveys and stock assessment which could result in unreliable management advice(McAllister, 1998; Waterhouse et al., 2014). It can also lead to failure of management strategies,particularly when considering space/time closures (Martell et al., 2000; Grüss et al., 2011). Theevaluation of impacts of fish movement on fisheries management often requires the use of migra-tion models and simulation studies. See Kerr and Goethel (2014) for a comprehensive review onthe application of simulation studies to evaluate the impacts of movement and migration on fishpopulation dynamics and fisheries management.Migration models are diverse in the fisheries literature (Goethel et al., 2011) and can be groupedbased on two numerical methods used to implement them: Eulerian and Lagrangian models. Theseterms are originally applied to fluid dynamics but have been used to describe migration modelingfor aquatic resources (Kerr and Goethel, 2014; Walters and Martell, 2004). The two approachesdiffer in the way that movement is measured. In the Eulerian approach space is divided in prede-termined areas and the movement rate of individuals across such areas is the variable of interest.Eulerian models have been coupled with stock assessment models and tagging studies (Carrutherset al., 2011; Methot and Wetzel, 2013; Sippel et al., 2015) and are useful when there is interestin tracking the net flux of individual across predetermined spatial areas. In the Lagrangian ap-proach, the movement of individuals is tracked through time and space and the movement tracksare the variable of interest (Walters et al., 1999). Lagrangian models rely on the assumption thatan underlying force drives movement. For example, this force could be environmental drivers forc-ing ichthyoplankton dispersal (Lett et al., 2008) or homing behavior driving salmonid runs (Caveand Gazey, 1994; Branch and Hilborn, 2010). Lagrangian models are not commonly applied toiteroparous species, despite the suitability of the approach to model migratory behavior. Migrationhypotheses, such as the cyclic movement between feeding and spawning areas, can be used as theunderlying pattern that drives the movement of individuals through space in a Lagrangian model.12In this paper, I describe a Lagrangian movement model designed for iteroparous species thatperform cyclic migrations throughout their lifetime. The objective of this model is to provide away of formalizing movement hypotheses into mathematical models. The resulting model can beused to summarize data and test the validity of alternative migration hypotheses and to representcomplex movement dynamics as an operating model in closed loop simulation studies. The modelis applied to the Pacific hake (Merluccius productus) as a study case. In particular, I focus on theoffshore Pacific hake stock that inhabits the California current system. This stock is believed toperform annual migration cycles between the spawning area off southern California and feedinggrounds along the West coast of North America, from Oregon to Southeast Alaska (Ressler et al.,2007). The migration cycle of Pacific hake is partially influenced by the age/length of individuals,with older/larger individuals reaching waters further north (Methot and Dorn, 1995; Ressler et al.,2007). The migration cycle is a key component to consider for management strategies for this stock.The agreement between Canada and the US defines an aggregate (i.e., non-spatial) harvest controlrule, but, given that prosecution of the fishery itself affects the mean age/size of fish, it also affectsthe distribution of the stock and hence, the distribution of the fishery’s benefits between the parties.In addition to the model description, I provide a description of data requirements to estimate themodel movement parameters. I also show how the model can be extended to incorporate covariatesrepresenting biological and environmental factors that alter the distribution and migration range ofthe populations being modeled.2.2 Methods2.3 Movement model frameworkI decompose the model into the following sections: population dynamics, movement, and fisheriesdynamics. This division is somewhat arbitrary as all three parts are interdependent, but the divisionis done to ease the description of the model. All three sections are structured by age, time, space(i.e., modeled area, fishing grounds and territories), and group. All model indexes, i.e., variablesused to designate the model dimensions, are presented in Table 2.1. The age dimension, denoted as13a, aggregates individuals by cohort, and spans from age one to a plus group A, which encompassesindividuals of age 20 and older. Individuals age is set to increase in the first month of each year,i.e., t = t0. The time dimension is denoted at t. I assumed monthly time steps, so all quantitieswithin the model were computed twelve times within a year/migrations cycle, but any step lengthcould be used. The indexing of time is cyclic, so that the range {t0, ..., tmax} is repeated every year.The first month of the migration cycle is also indexed with the year counters (y), so that y−1 = t0in the previous year. The space dimension is considered at three scales: area, fishing ground, andterritory. The variable area (r) refer to small and equally sized intervals of space, which are usedto discretize the variables of interest (biomass, fishing effort, etc.). The area range denotes thelimits of the modeled space. Fishing grounds and territories are larger than areas, but are containedwithin the area range (i.e., within the limits of the modeled space). The fishing grounds correspondto areas where a fleet is known to operate and territories correspond territorial waters of a nation orstate. Each territory may contain one or more fishing grounds, and the biomass within one territoryis only accessible to the fleets operating within that territory. Lastly, the group dimension is usedto represent parcels within an age cohort that move at different speeds.In the population dynamics section (Table 2.3) the processes relating to recruitment, survival,aging and growth are described. These processes are modeled for each group in the populationthrough time and space. The spatial section (Table 2.4) encompasses the description of movementof individuals at age and group through time. And finally, the fishing dynamics section (Table 2.5)describes the model effort distribution as a function of spatio-temporal effort scalers and the spatialdistribution of fish biomass (summed over all ages and groups). The spatial distribution of fishingeffort is used to generate of spatially explicit fishing mortality. In the movement section, groups aremodeled in two ways yielding two versions of the model: single group or multiple groups versions.All model variables are defined in Table 2.2.The population dynamics is composed of survival and recruitment processes (Table 2.3). Themodel assumes age-specific vulnerability and fecundity as well as Beverton-Holt recruitment oc-curring at age 1. Recruitment and aging are assumed to happen at the first time step within a14Table 2.1: Indexes and variable definitionsSymbol range Descriptiont {t = t0, ..., tmax} time steps within a migration cycley {y = yo, ...,Y} year indexa {a = 1, ...,A} age indexr {r = rini, ...,R} area indexk {k = 1, ...,K} fishing ground indexn {n = 1, ...,N} territory indexg {g = 1, ...,G} group indexRM range of modeled areadr interval between two adjacent areastmax ·Y total numbers of time steps modeledmigration cycle (Equation T3.2 - cases i-iii) and survival due to natural and fishing mortality arecalculated at monthly time steps (Equation T3.2 - cases ii-iv). Because of the continuous nature ofthe model, sometimes some portion of the population might be located outside the boundaries ofmodeled area and therefore outside the fishing grounds. When this happens, we assumed that indi-viduals are subject to natural mortality only (Equation T3.2 - second term on cases ii-iv). Spawningbiomass is combined over all areas and groups (Equation T3.3).The movement section (Table 2.4) assumes the individuals in a population are distributed alonga unidimensional gradient X , and that they perform annual cyclic migrations between spawningand feeding areas. This cyclic migration is modeled with a sine function in which the position ofindividuals change as a function of time (Figure 2.1). I have developed two alternative versions:single group and multiple groups. In both versions, the cyclic movement of individuals betweenspawning and feeding areas is modeled with a sine curve (Equation T4.1). The magnitude ofthe movement is determined by two parameters representing maximum and minimum positionsof the migration cycle (Xmax,a,g and Xmin), and the starting time step of the cycle is given by theparameter t0. The maximum and minimum positions of the cyclic movement, Xmax and Xmin, canbe modeled as a function of covariates, such as age, size or environmental drivers. Here, I modelthe maximum position as a logistic function of age (Equation T4.2) and fix the minimum positionas constant for all ages. Once the maximum and minimum positions are defined, the sine curveis used to calculate the mean position of individuals in each group of age a, in time step t (X¯a,t).15Table 2.2: Variable definitionsSymbol DescriptionPopulation dynamicsNa,t,g Numbers at age, time and groupM Annual natural mortalityva Vulnerability at ageS0 Maximum juvenile survival rateSBt Spawning biomass at time tβ Beverton & Holt density dependencepropg Proportion of total recruitments that recruits to group gwt Normally distributed recruitment errorPra,r,t,g Proportion of population of age a, group g and at time t in area rwa Weight at agefa Fecundity at agef Proportion of females in the populationΦ(x|µ,σ) Cumulative normal distribution with mean µ and standard deviationσ , evaluated at xSpatial dynamicsX¯a,t mean position at age and at time tX¯min minimum average positionX¯max,a maximum average position at aget0 time step at which individuals are at their minimum average positiona50 inflection point for maximum average position logistic functionσXmax standard deviation for maximum average position logistic functionσXa standard deviation for position at ageCV coefficient of variation for position at at agesingle group versionXa,t probability distribution for position at age and timemultiple groups versionX¯a,t,g mean position at age, time t and for group gQXa,t,g quantiles for mean group distributionδ distance between quantiles for mean group distributionσXa,g standard deviation for position at ageX¯a,t,g mean position at age and groupXa,t,g probability distribution for position at age and time for each groupFisheries dynamicsEr,t Effort in area r and at time tEy,k Maximum yearly effort scaler in fishing ground kEt,k Maximum time step effort scaler in fishing ground kAr,t Attractiveness of each area r and at time tλ Attractiveness powerV Br,t Vulnerable biomass by area and timeEpot,r Relative effort potential for each area rFr,t Fishing mortality rate in area r and at time tq Catchability coefficientCa,r,t Catch at age for each small area r and time t16Table 2.3: Population dynamics modelInitial year and time step - y = yoNa,t,g =S0·SB01+β ·SB0 · exp(wt) · propg, (i) t = t0 & a = 1Na−1,t,g · exp(−M/tmax), (ii) t = t0 & 1 < a < ANa−1,t,g·exp(−M/tmax)1−exp(−M/tmax) , (iii) t = t0 & a = A(T3.1)Age-schedule informationNa,t,g =S0·SBy−11+β ·SBy−1 · exp(wt) · propg, (i) t = t0 & a = 1r∈RM∑ (Na−1,t−1,g ·Pra−1,r,t−1,g · exp(−M/tmax−q ·Er,t−1 · va−1))+Na−1,t−1,g ·(1−r∈RM∑ Pra−1,r,t−1,g)· exp(−M/tmax), (ii) t = t0 & 1 < a < Ar∈RM∑(Na−1,t−1,g·Pra−1,r,t−1,g·exp(−M/tmax−q·Er,t ·va−1)1−exp(−M/tmax−q·Er,t−1·va))+Na−1,t−1,g·(1−r∈RM∑ Pra−1,r,t−1,g)·exp(−M/tmax)1−exp(−M/tmax) , (iii) t = t0 & a = Ar∈RM∑ (Na,t−1,g ·Pra,r,t−1,g · exp(−M/tmax +q ·Er,t−1 · va))+Na,t−1,g ·(1−r∈RM∑ Pra,r,t−1,g)· exp(−M/tmax), (iv) t0 < t ≤ tmax &1≤ a≤ A(T3.2)SBt =g∑a∑Na,t,g · fa ·wa · f (T3.3)propg ={1, single groupΦ(QXg +δ |µ = 0,σ = 1)−Φ(QXg−δ |µ = 0,σ = 1), multiple groups(T3.4)17The proportion of individuals in each area, i.e., small intervals of space, is given by the cumulativenormal distribution function given the mean position and standard deviation of the cohort (EquationT4.5) or group (Equation T4.11). X spatial gradientTime1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 2430354045505560t0 t0XmaxXmidXminXmidXmaxFigure 2.1: Diagram illustrating the sine function used to describe the cyclic movement dy-namics used in the Lagrangian model. The position (y-axis) of an individual changes incyclic waves as time (x-axis) progresses. In this diagram, we assume that a cycle lasts12 time steps and that movement occurs between values of 30 and 60 along the spatialscale.The differences between single and multiple groups versions are in how the individuals aredistributed around the mean position (X¯a,t). When a single group is present, the individuals areassumed to be normally distributed around the mean with variance σ2Xa (Equation T4.4). In themultiple groups version, each group’s position follows a group specific normal distribution andthe mean and variance of each group is a function of the overall mean and variance (EquationsT4.6-T4.9). The main difference between the two versions is that in the single group version, thedistribution of fish regenerates to a normal distribution at every time step, but with a new averageposition. Although individual fish might experience different fishing mortality depending on lo-cation, the regeneration assumption does not allow localized depletion to occur. Instead, fish areredistributed across their range in each time step following a normal density function. When multi-ple groups are considered, each group is distributed according to a group specific distribution that ismuch narrower than the overall distribution of individuals at an age class. This mechanism allows18Table 2.4: Lagrangian movement modelMovement informationX¯a,t = X¯min +(X¯max,a− X¯min) ·(0.5+0.5∗ sin(t · 2pitmax− t0 · 2pitmax −pi2))(T4.1)X¯max,a =11+ exp(−(a−a50)/σXmax)· e(vt∼N (0,σvt)) (T4.2)σXa = X¯max,a ·CV (T4.3)Single group versionXa,t ∼N (X¯a,t ,σXa) (T4.4)Pra,r,t =Φ(x = r+dr2|µ = X¯a,t ,σ = σXa)−Φ(x = r−dr2|µ = X¯a,t ,σ = σXa) (T4.5)Multiple groups versionX¯a,t,g = QXa,t,g ·σXa + X¯a,t (T4.6)QXg = δ · (g−G/2.0) (T4.7)δ = 6.0/G (T4.8)σXa,g =√σXa 2G2(T4.9)Xa,t,g ∼N (X¯a,t,g,σXa,g) (T4.10)Pra,r,t,g =Φ(x = r+dr2|µ = X¯a,t,g,σ = σXa,g)−Φ(x = r−dr2|µ = X¯a,t,g,σ = σXa,g) (T4.11)certain groups to experience different fishing pressures depending on where the group is located,which might lead to higher or lower fishing pressure over extended periods of time. Because thedistribution range of each group is narrower the regeneration only occurs within a narrow range,which allows local depletion to become apparent.In the fishing dynamics section (Table 2.5), spatial fishing effort allocation is done with a gravitymodel (Caddy, 1975) . These models assume that effort in each area is a function of the latent yearlyand monthly effort in a fishing territory and the attractiveness index of that area (Equation T5.3).The effort potential quantities Ey,k and Em,k are scaling matrices of dimensions Y x K and tmax x Krespectively. The values in these matrices range between 0 and 1. Values of 0 indicate that no effort19is allowed in that particular year or month, conversely, values of 1 indicate that full effort is allowedin that particular year or month. The attractiveness index can include factors such as fishing cost,target fish abundance and bycatch abundance (Caddy, 1975; Walters and Martell, 2004). I makeattractiveness (Equation T5.2) a function of vulnerable biomass V Br,t−1, the power parameter λand effort potential Epot,r (Equation T5.1). I use the vulnerable biomass in the previous time stepbecause I assume that effort distribution is guided by the biomass distribution observed in theprevious time step. The parameter λ is used to indicate if the attractiveness is directly proportionalto abundance (λ = 1), or if the fleet tends to disproportionately aggregate in high abundance areas(λ > 1). One example in which λ > 1, is when there is communication between fishing vesselswhen a high abundance area is located. In such cases, fishing effort would tend to aggregate in highabundance areas generating a patchy distribution of effort. The effort potential Epot,r parameterrepresents the avoidance factors for a given area, e.g., fishing costs and bycatch avoidance. Thisparameter can be used to represent a range of differences in fishing fleets, such as storage capacity,autonomy at sea, and distance between fishing grounds and home port. The Epot,r parameter canalso be used to represent areas that are avoided due to high bycatch occurrence. Avoidance areasaffect the ability to concentrate fishing effort at a given location, and therefore should be consideredin the modeling process. Effort is then multiplied by q, the effort scaler, resulting in the fishingmortality in that area. Lastly, catches are calculated for each time step using the Baranov catchequation (Equation T5.4).Process and observation random errors are incorporated in the model. The process random errorwas represented annual recruitment variability, annual variability in the maximum average position,and annual variability in the effort scaler q. All these variability components were modeled withlognormally distributed error. The age composition data, in numbers and aggregated by fishingground, are generated with multivariate logistic sampling error.20Table 2.5: Fisheries dynamics modelFisheries dynamicsV Br,t−1 =a∑g∑va ·Na,t−1,g ·wa ·Pra,r,t−1,g (T5.1)Ar,t = V Br,t−1r∈n∑ V Br,t−1λ ·Epot,r (T5.2)Er,t = Ey,k ·Et,k · Ar∑r Ar(T5.3)Ca,r,t =Fr,t · vaFr,t · va +M ·Na,r,t · (1− exp(−(Fr,t · va +M))) (T5.4)Fr,t = q · e(wx∼N (0,σwx)) ·Er,t (T5.5)2.3.1 Application to Pacific hake and simulation-estimation procedureThe offshore Pacific hake stock was used as an example to illustrate the model dynamics. In thePacific hake case, the fish are known to perform annual migration between the waters off SouthCalifornia and northern British Columbia. Therefore, I model the movement of hake in terms oflatitude degrees, from 30◦N to 60◦N (Figure 2.2). The driving population dynamics parameters forPacific hake were obtained from the 2015 stock assessment (Taylor et al., 2015). The parametersfor the effort dynamics were set to mimic the fisheries dynamics described in the stock assessmentdocument (Taylor et al., 2015). The parameters used in the simulation-estimation procedure arelisted in Table 2.6.I did a simulation experiment to evaluate the estimability of the movement parameters. I sim-ulated and estimated population movement dynamics using both single group and multiple groupsversions (20 groups). I simulated total catch and catch at age data with observation error, and usedthat data to estimate the models parameters.The model parameterization is divided into two categories: fixed (extracted from other models),and parameters that could be estimated, given seasonal catch at age data. The fixed parametersinclude all the population dynamics and fisheries capacity parameters. These parameters were therecruitment function parameters (Ro and h) and natural mortality (M). These parameters were21fishing ground 1fishing ground 2fishing ground 3fishing ground 4fishing ground 5CanadaU.S.A303540455055−135 −130 −125 −120 −115LongitudeLatitudeFigure 2.2: Map of the study area. The dashed lines indicate the division between fishinggrounds.22Table 2.6: Pacific hake model dimensions and parameter valuesSymbol value or range DescriptionDimensionst 1−12 Time steps within a migration cycley 1−30 Yearsa 1−20 Ager 30−60Areak 3 or 5 Fishing grounds42, 46, 48.5 and 51 Fishing ground boundaries in latitude degreesn 2 Territories48.5 Territory boundary in latitude degreesg 1−20 Groupsdr 1 Interval between two adjacent areastmax ·Y 360 Total numbers of time steps mod-eledPopulation dynamics parametersM 0.223 Annual natural mortalityS0 15.331 Maximum juvenile survival rateβ 5.422 Beverton & Holt density dependenceMovement parameterst0 1 Time step at which individuals are at their mini-mum average positionCV 0.1 Coefficient of variation for position at at agea50 3.0 Inflection point for maximum average position lo-gistic functionσXmax 1.5 Standard deviation for maximum average positionlogistic functionerror levelsτ 0.4 or 1.0 Standard deviation for the multivariate logistic er-ror around catch at ageσwx 0.08 Standard deviation for lognormal variation aroundthe maximum average positionσvt 0.1 Standard deviation for lognormal variation aroundthe effort scalerEffort parametersEy,k Constant for all years and equal to (1,1, 1, 0.2, 0.2) for fishing grounds from1 to 5Yearly effort scaler - fishing grounds 1-2 and 4-5were combined when only 3 fishing grounds wereconsideredEt,k (0.0 0.0 0.0 0.0 0.5 1.0 1.0 1.0 0.5 0.10.0 0.0) for fishing grounds 1-3 and(0.0 0.0 0.0 0.0 0.0 1.0 1.0 1.0 0.5 0.30.0 0.0) for fishing grounds 4-5Monthly effort scaler - fishing grounds 1-2 and 4-5were combined when only 3 fishing grounds wereconsideredq 1 Effort scaler23considered as a fixed input to the model and were assumed to be known without error in bothsimulation and estimation models. The estimable parameters are: t0, CV , a50, σXmax and q. Theseparameters were estimated with a multivariate logistic likelihood function fitted to simulated agecomposition data.A total of 12 simulation-evaluation scenarios were considered in this study (Table 2.7). I con-sidered two data aggregation scenarios with data reported from three or five large fishing grounds(aggregated over all areas within fishing ground). These two levels of data aggregation were con-sidered in order to explore the sensitivity of the model to levels of spatial aggregation of the agecomposition information. I initially considered these two levels of data aggregation for the singleand multiple groups version, however I got very low levels of model convergence (50% or less)when three fishing ground were considered using the multiple groups version. This fact lead me todrop the three fishing ground scenarios when the multiple groups version was used. The aggregatedcatch at age data was assumed to have multivariate logistic error with two levels of observation er-ror, i.e., the standard deviation (τ) being either (0.4 or 1.0). These two levels of variability in theage composition data were chosen in order to evaluate the sensitivity of the model to measurementand sampling errors. In addition, I considered two attractiveness scenarios (λ = 1 and λ = 2). Thedifferent levels of attractiveness are likely to change the spatial distribution of fishing effort, andwould likely impact the degree of depletion experienced by fish in different areas. This is par-ticularly relevant for the multiple group scenarios as the different λ values might affect the localdepletion patterns.The estimation model was identical to the simulation model, and parameters were estimatedwith a multivariate logistic likelihood function fitted to simulated age composition data and alognormal likelihood fitted to the simulated total catches. A total of 100 simulations were runfor each scenario-version combination. Simulation and estimation routines were performed us-ing ADMB (Fournier et al., 2012). The code and simulated data are available for download from 2.7: Simulation-estimation scenariosScenario number Model Version catch at age error Fishing grounds λ value1 single group low (τ = 0.4) 3 1.02 single group low (τ = 0.4) 3 2.03 single group high (τ = 1) 3 1.04 single group high (τ = 1) 3 2.05 single group low (τ = 0.4) 5 2.06 single group low (τ = 0.4) 5 1.07 single group high (τ = 1) 5 1.08 single group high (τ = 1) 5 2.09 multiple groups low (τ = 0.4) 5 1.010 multiple groups low (τ = 0.4) 5 2.011 multiple groups high (τ = 1) 5 1.012 multiple groups high (τ = 1) 5 2.02.4 Results2.4.1 Model dynamicsFigure 2.3 shows spatial distribution of biomass in the absence of fishing for the single and multiplegroups versions. If fishing is absent, the biomass spatial distributions are practically identical.However, the spatial distributions for the two models tend to change if fishing is present and theeffort is not homogeneously distributed throughout the unfished distribution of the resource (Figure2.4). For the multiple groups version, the spatial distribution of fish at each age tends to deviatefrom the initial normal distribution assumption as the fish grow older. This distortion is caused bythe fact that not all groups are subject to the same effort intensity, hence they encounter differentfishing mortality, and therefore depletion levels over their lifetime.In the Pacific hake example, effort is assumed to be concentrated towards the northern areas(bars on Figure 2.4). Therefore, when the multiple groups version of the model is considered(dark line on Figure 2.4), the groups that move to higher latitudes tend to be subject to strongerfishing pressure, and therefore become more depleted than groups that remain in lower latitudes.This effect cannot be detected for the younger ages (Figure 2.4, age 1). However, as cohort ages,the groups that move further (located to the right of the mean distribution for the entire cohort)start exhibiting perceptibly higher depletion levels compared to groups within the same cohort that25do not migrate as far (located to the right of the mean distribution for the entire cohort, Figure2.4, age 5). The higher depletion levels on the fish that move further also causes the mean in theoverall distribution of each cohort to shift to the south, which, over time would tend to diminish theavailability of fish to the fishing fleets in the northern areas.2.4.2 Simulation-estimationThe simulation-evaluation analysis showed that the five key parameters of the movement modelcan be estimated given spatial catch at age data, assuming that the model assumptions are satisfied.When data were simulated using the single group version (scenarios 1-8 - Table 2.7), it was possibleto predict the main movement parameters with practically no bias, i.e., parameter estimates werewithin 10% of true values (Figure 2.5). Variability in parameter estimates was lower when data werereported for five fishing grounds when compared to scenarios where data were reported only forthree fishing grounds. When the data reporting was assumed to occur for five fishing grounds, thedata was less aggregated and therefore more informative to the estimation of the movement patternalong the latitudinal migration route of the fish. As expected, higher values of τ , the standarddeviation for the catch at age error, resulted in higher variability in the parameter estimates. Nodifference was observed for scenarios with different values for the λ parameter. The parameter λallows effort to concentrate in areas of high abundance. It is likely that the effects different valuesof λ were not noticeable because of the regeneration assumption that is inherent to the single groupversion of the model. In other words, even though effort might have aggregated in certain areas, itdid not affect the overall distribution of fish for each cohort. 40 50 60LatitudeBiomass modelmultiplesinglesum_multipleFigure 2.3: Illustration of the distribution of fish of age 5, in the month of July for the singleand multiple groups versions in the absence of fishing. The multiple groups versionoutput is shown for all groups individually and combined. The vertical dashed linerepresents the border between territories.27Oct Nov DecJul Aug SepApr May JunJan Feb Mar30 40 50 60 30 40 50 60 30 40 50 600.000.250.500.751. (areas)Relative Biomass/Effortmodelmultiplesingleage15Figure 2.4: Monthly representation of differences between spatial distribution of biomass forsingle and multiple groups versions. Vertical bars indicate relative fishing effort. Con-tinuous and dashed lines represent vulnerable biomass for fish at age 1 and 5 throughoutone year migration cycle, one panel for month). Different shades of gray indicate mul-tiple and single group versions. The vertical dashed line represents the border betweenterritories. Both fish biomass and effort were re-scaled from 0 to 1 for plotting.28t0 CV a50 σXmax q-, 3 FG, τ = 0.4, λ = 1.0t0 CV a50 σXmax q-, 3 FG, τ = 0.4, λ = 2.0t0 CV a50 σXmax q-, 3 FG, τ = 1.0, λ = 1.0t0 CV a50 σXmax q-, 3 FG, τ = 1.0, λ = 2.0t0 CV a50 σXmax q-, 5 FG, τ = 0.4, λ = 1.0t0 CV a50 σXmax q-, 5 FG, τ = 0.4, λ = 2.0t0 CV a50 σXmax q-, 5 FG, τ = 1.0, λ = 1.0t0 CV a50 σXmax q-, 5 FG, τ = 1.0, λ = 2.0t0 CV a50 σXmax q-, 5 FG, τ = 0.4, λ = 1.0t0 CV a50 σXmax q-, 5 FG, τ = 0.4, λ = 2.0t0 CV a50 σXmax q-, 5 FG, τ = 1.0, λ = 1.0t0 CV a50 σXmax q-, 5 FG, τ = 1.0, λ = 2.0relative errorParameterFigure 2.5: % Relative error for estimation of key parameters of the movement model. Scenarios are indicate in the upper rightcorner of each plot and the scenarios composition is indicated in the plot titles. SG = single group version , MG= multiplegroups version and FG = fishing grounds.29When data were simulated with the multiple groups version (scenarios 9-12 - Table 3.6), bias inthe estimated parameters became more prominent (Figure 2.5). The effort scaler parameter, q, wasunderestimated in all scenarios by about 40%. In order to evaluate the impact of the parameter biason the predicted biomass distribution, I used plotted the median predicted biomass for scenario 12,one of the cases where the bias was most prominent (Figure 2.6). Despite the bias in parameterestimates for the multiple groups version the impact on the distribution of total biomass over timepredicted by the model was relatively small (Figure 2.6). Because of the lower effort predictionsdue to the underestimation of q, higher median biomass distributions are predicted by estimationmodel. However, little impact is seen in the overall proportion of biomass in each area.2.5 DiscussionThe Lagrangian movement model described in this paper is an alternative to traditional Eulerianapproaches commonly used to model the distribution of adult iteroparous fish (Goethel et al., 2011;Sippel et al., 2015). The Lagrangian approach shown here allows for the explicit consideration ofmigration hypotheses, such as the cyclic migration between spawning and feeding grounds, andexploration of potential impacts of covariates in shaping migration variability within a population.Traditional Eulerian models are generally represented by spatially discrete box or bulk transfermodels (Carruthers et al., 2015; Methot and Dorn, 1995) or continuous advection-diffusion mod-els (Sibert et al., 1999). Box models are simpler but require the predetermination of, usually large,spatial areas from which flow (i.e., movement) is measured. The determination of such areas can bechallenging, especially due to the assumption that fishing mortality is homogeneously distributedwithin each area and that flow in between boxes is mainly caused by migration or diffusive move-ment (Walters and Martell, 2004; Carruthers et al., 2011). Frequently, these large spatial areasare determined by political or management boundaries or historical division of data, and do notcorrespond to an ecologically relevant partition of the habitat of the species being studied. Suchartificial partitions can lead to violation of the box model assumptions. The Lagrangian model pre-sented here is continuous in space, and therefore does not require the definition of homogeneous301975197619771978197919801981198219831984198519861987198819891990199119921993199419951996199719981999200020012002200320042005200620072008200920102011201220135 10 15 20Age (years)YearPacific hake19801981198219831984198519861987198819891990199119921993199419951996199719981999200020012002200320042005200620072008200920102011201220132.5 5.0 7.5 10.0 12.5Age (years)Yearjack mackerelFigure 2.6: Median biomass distribution for simulation and estimation models for each monthof the year for scenario 12 - multiple groups, 5 areas, τ=2.0.31spatial areas. However, the model output can be aggregated in larger areas for comparison withhistorical data or for investigation of management questions relating to political boundaries.Advection-diffusion models are continuous in space, but are much more data intensive andgenerally require the availability of tagging and tracking data (Sibert et al., 1999; Costa et al., 2012).Similarly to advection-diffusion models, the Lagrangian model presented in this paper is continuousin space and time, allowing for predictions of biomass at any location in the time-space continuum.However, differently from the diffusion models, the approach shown here does not assume thatanimals move at random. Instead, the movement is assumed to be directed by an innate migrationhypothesis, frequently derived from observed seasonal size or age distribution of the species ofinterest (e.g., Ressler et al., 2007; Francis and Clark, 1998). The use of a migration assumptionreplaces the need to directly estimate advective terms to explain fish movement and, therefore, doesnot rely on tagging data to determine the direction of fish movement. This is an advantage becausefor many fish species tagging data is not readily available and tagging studies can be difficult tocarry out. This is particularly true for deep water species due the complications and high mortalityrates resulting from the barotrauma caused by bringing the individuals to the surface for tagging(Nichol and Chilton, 2006; Winter et al., 2007). In addition, tagging studies usually involve highoperational costs, which include tag deployment, recapture surveys and/or publicity campaigns toincrease voluntary tag reporting rates in the fisheries and, in some cases, high costs associatedwith the tags themselves (Pine et al., 2003). Despite the high costs and the large effort associatedwith tagging studies, the data quality is frequently diminished by violation of basic assumptionsof tagging experiments, such as changes in reporting rates over time, time-varying catchabilityand influence of the tag on fish behavior (Pine et al., 2003). These difficulties associated withtagging data frequently result in the lack of enough information in the data, despite extraordinaryexpense and effort to accurately resolve the movement dynamics of the species of interest withouta number of simplifying assumptions. Alternatively, the movement model presented here reliesmainly on seasonal catch and age composition data to estimate its movement parameters. Catchand age composition data are conventionally available for many temperate exploited fish species32(Melnychuk et al., 2017) and the collection of such data is already part of the management programsfor such species, i.e., used in stock assessments. Therefore, the Lagrangian model presented can beused as an alternative tool to formalize and test the movement hypotheses based on data, even whentagging or movement track data are not available. It is important to note, however, that the approachdescribed here is not meant to be a complete replacement for tagging studies. There are severalcaveats associated with the present approach, including the assumption of known population andfishing dynamic parameters. Such parameters are seldom known with certainty, and if the assumedvalues deviate significantly from the reality, the model outputs would likely be severely biased aswell.Appreciable differences between the two versions of the model, i.e., the single group and mul-tiple groups versions, were observed. These differences are associated with the way in whichmovement of each age class is treated. In the single group version, all individuals within a co-hort are considered equal and individuals are assumed to be normally distributed around the meanposition at time. In the multiple groups version each cohort is sliced into groups which have agroup-specific mean position at time. The single group model is simpler and, for this reason, muchmore computationally efficient. However, it relies on the assumption of regeneration of spatial dis-tribution at each time step. That is, fish of a given age are assumed to spread spatially followinga normal distribution at each time step regardless of possible local depletion due to concentrationof fishing effort. This might be a reasonable assumption if the distribution of effort is relativelyhomogeneous over the distribution of the age group being exploited. However, effort concentrationand local depletion is not uncommon, and has been observed for many fish stocks around the world(e.g., Maury and Gascuel, 2001; Bez et al., 2006). Alternatively, if fishing effort is known to berestricted to some areas, or is known to concentrate heavily in some areas, then using the multiplegroups version is more appropriate. I found that a minimum of 15 groups is required to appropri-ately reduce the effects of regeneration of spatial distribution. The multiple groups formulation willallow local depletion to occur and persist through the life of each cohort.33An important feature of the simulation approach presented is that it explicitly accounts for spa-tial fishing effort with the gravity model (Caddy, 1975). Because the effort dynamics are modeledexplicitly, it is possible to use the model to investigate questions relating to active avoidance ofspecific areas due to high bycatch occurrence or spatial aggregation in fishing effort. In the specificexample of Pacific hake, the U.S. fleet has a strong incentive to avoid areas where abundance ofbycatch species is high, despite potential high abundances of the target species. This effect can bemodeled by linking the Epot,r to the abundance and distribution of bycatch species. Spatial aggre-gation of fishing effort is also important because if fishers aggregate in areas of high abundance,it is likely that strong cohorts will be targeted disproportionately causing variation in selectivityacross time. Time varying selectivity can be confounded with fishing mortality estimates leadingto biases in the estimation of reference points and management targets (Punt et al., 2015; Martelland Stewart, 2014).I found that it is possible to estimate the driving movement parameters of the model using spatialcatch at age information collected throughout the migration cycle of the species being modeled.Data were generated and estimated using both single and multiple groups versions of the model.This procedure resulted in unbiased parameter estimates for the single group simulation scenariosand up to 40% median relative error in parameter estimates for multiple group scenarios. Thecauses for biased parameter estimates in the multiple groups version is unclear. However, thehigher data requirements of the multiple groups version is understandable given the higher degreeof complexity in the model. I have not tested the performance of this model assuming catch atage from more than five fishing grounds, but it is possible that less aggregated data would resultin better resolution of the model parameters. However it is also important to consider that lessaggregated data will probably have higher observation error levels and that might also impair theestimation of the model parameters.It is possible that the relative error estimates obtained in the simulation-estimation analysis areover optimistic because the same model structure was used for simulating data and for estimatingparameters. This similarity is likely to have improved the realized model fit. In addition, the pa-34rameter estimates obtained in this study are dependent on assumptions regarding the fixed modelparameters, i.e., the fisheries and population dynamics parameters. The true values of such param-eters are usually not well known, and estimates can change dramatically over time (Brooks et al.,2015). However, when we attempted to jointly estimate the movement and population dynamicsparameters, i.e., use the current model as a stock assessment model, there was a great amount ofconfounding between the estimates of productivity, recruitment deviations and the movement pa-rameters. Therefore, it is unlikely that the movement dynamics presented here could be integratedinto stock assessment models.A promising application of the Lagrangian model described here is its potential to be used asan operating model in closed loop simulations. Such simulations can be used to evaluate the ef-fects of management strategies for exploited fish populations (Giske et al., 2001; Sainsbury et al.,2000). The model can be used to represent the complex population dynamics of migratory species,as well as the variability in distribution of stocks due to intrinsic (e.g., growth, maturity) and ex-trinsic (e.g., environmental forcing, fishing effort) forces to the population. One advantage is thatthe mechanics of this model is different from that usually implemented in stock assessment models(e.g., Methot and Wetzel, 2013; Fournier et al., 1998), which would yield significant differencesbetween operating and estimation models in closed loop simulations. Similarities between operat-ing and estimation models lead to improved performance of estimation models, which in turn resultin overly favorable performance of management strategies (McClure et al., 2014). In reality it isunlikely that the estimation models capture all the processes that occur in nature.A couple model extensions were not included in this chapter, but could easily be added to themodel dynamics. These extensions are the addition of multidimensional movement and the use ofother mathematical functions to represent migration hypotheses. The model illustration presentedin this study only describes movement in a unidirectional basis, that is from spawning to feedinggrounds. This simplification of movement trajectory stems partly from the knowledge of the specieschosen as a case study. A unidimensional model is commonly used to describe the movement dy-namics of Pacific hake, and not much is known regarding the population trajectory and habitat use35in other dimensions (Ressler et al., 2007). However for many species, migration routes are morecomplex and involve simultaneous migration between spawning and feeding grounds, shallow anddeep water and between inshore and offshore grounds (e.g., Misund et al., 1998; Barbaro et al.,2009; Merten et al., 2016). The current approach could be easily extended to a multidimensionalapproach. The addition of new dimensions, however, would require the development of mathe-matical functions that describe movement in each dimension. The sine function presented here isa good candidate to describe cyclic movements, but linear, logistic or knife-edge functions couldalso be used if movement in a given dimension is thought to be permanent, as would be true formigration between nursery and rearing grounds.I believe that the model presented here is a useful approach to model movement of migra-tory fish species. I anticipate that the model can be used to examine the plausibility of differentmovement hypotheses and to explore the possible links between fish migration and ecosystem in-teractions. In addition, we suggest that the model is a good candidate to be used as an operatingmodel in closed loop simulations, especially when there is interest in evaluating the implication ofmigratory movement on management outcomes.36Chapter 3Stock Reduction Analysis using catch atlength data: Length-SRA3.1 IntroductionModern stock assessments typically attempt to fit population dynamics models to catch at ageand/or catch at length data, in hopes of extracting information from these data about age/size selec-tivity, cohort strength, and fishing mortality patterns (Methot and Wetzel, 2013; Hilborn and Wal-ters, 1992). Some assessment methods attempt to put aside the length frequency data, by convertingthese data to age compositions using age-from-length tables, perhaps using iterative methods to es-timate proportions of fish at age for each length interval (Kimura and Chikuni, 1987). In caseswhere age data are lacking, models like MULTIFAN-CL attempt to obtain estimates of selectivity,fishing mortality and population dynamics parameters only from size distribution data (Fournieret al., 1998). Combined with a few assumptions regarding the structure and variability in lengthat age, this procedure can even be used to attempt to recover information about changes in bodygrowth patterns if there is a strong age-class signal in the length frequency data (Fournier et al.,1998). It is typical for assessment results from length-based assessment models to show substantialdeviations between predicted and observed length distributions of catches, reflecting both sampling37variation in the length composition data and incorrect assumptions about stability of growth andselectivity patterns (Hilborn and Walters, 1992).Selectivity to fishing is the combination of two processes: vulnerability to the fishing gearand availability of the fished population in the area being fished (Beverton and Holt, 1957). Bothprocesses can vary over time and therefore modify the resulting selectivity. Although selectivityprocess can often be directly measured through gear experiments, availability is generally harderto measure as it depends on the size-based distribution of the exploited population and the spatialdistribution of the fishing fleet. Fish movement, size-structured changes in fish distribution, andchanges in fleet distribution, can all affect availability and consequently lead to selectivity changes.Changes in selectivity are not uncommon (Sampson and Scott, 2012) but are usually difficult totrack over time. This difficulty is associated with an inability to distinguish between changes infishing mortality and changes in selectivity in most age- and length-based stock assessment meth-ods. For this reason, many assessment methods rely on ad hoc parametric selectivity models thatmay or may not include changes over time (Maunder et al., 2014). If misspecified, such modelsmight lead to severe bias in fishing mortality estimates and other model parameters, which couldresult in misleading management advice (Martell and Stewart, 2014).Here, I introduce an alternative approach to assessment modeling that begins by assuming thatthe assessment model should exactly reproduce the observed catch at length composition. Thisapproach follows the dynamics of an age structured stock reduction analysis (SRA) (Walters et al.,2006; Kimura et al., 1984; Kimura and Tagart, 1982), which follows a “conditioned on catch”format, in which catch composition is assumed to be known without error. The observed catchesat age are then subtracted from modeled numbers at age to project numbers at age over time. Agood review of SRA-type models is provided in Thorson and Cope (2015). The assumption ofknown catch composition is analogous to the classical assumption in virtual population analysisthat reconstructed numbers at age should exactly match observed catch at age data (Hilborn andWalters, 1992). The suggested approach may have two key advantages over statistical catch at ageand/or catch at length models: (1) it does not require estimation of age or size selectivity schedules,38and (2) catch at length data are commonly available for every year, even when age compositionsampling has not been conducted.I have named this approach a Length-SRA assessment model. Here, I present the model formu-lation, demonstrate its performance with a simulation-evaluation analysis, and apply it to fisheriesdata from the Peruvian jack mackerel (Trachurus murphyi) and Pacific hake (Merluccius productus)fisheries.3.2 Methods3.2.1 Stock reduction analysis with catch at length data - length-SRAThe stock reduction analysis (SRA) described here proceeds through the following steps: (1) com-pute numbers at age (based on recruitment estimates and mortality in the previous year); (2) convertnumbers at age into numbers at length using the proportions of individuals at length given each ageclass; (3) calculate the exploitation rate at length using numbers at length and observed catch atlength; (4) convert the exploitation rate at length to exploitation rate at age; (5) compute numbersin the following year using the exploitation rate at age, natural mortality, and recruitment estimates.The model requires data on length composition of catch in numbers (used in step 3), a priordistribution for the recruitment compensation ratio, and a survey index of abundance that is used totune the model parameters to the most likely stock abundance trajectory. The model also requiresgood estimates of growth parameters, variability around mean length at age, and natural mortality.The stock assessment and simulation routines were written in ADMB (Fournier et al., 2012) andare available on crucial component of the length-SRA is the calculation of proportions of individuals at lengthgiven each age class (Pl|a - eqs. T3.1-T3.5). The calculation of such proportions (eq. T3.1) relieson four main assumptions regarding the distribution of length at age: (1) The mean length at agefollows a von Bertalanffy growth curve (eq.T3.4), (2) The length at age is normally distributed (eqs.T3.1 -T3.3), (3) The standard deviation of the length at age is defined (e.g., eq.T3.5), and (4) PL|ais constant for all lengths equal or greater than a maximum length L (eq.T3.3).39The proportions of length at age are used to convert the length-based quantities into age-basedquantities, which are used to propagate the age structured population dynamics forward (Table 3.3).I assume that recruitment follows a Beverton-Holt type recruitment curve (eq. T3.6), that harvestingoccurs over a short, discrete season in each time step (year or shorter), and that natural survival rateis known and constant over time (Equations T3.6-T3.10). The computation of numbers at age inthe initial year (i.e., first year in which data is reported - t = init) is different from that in theremaining years (Equation T3.13). Recruitment in the initial year is set to the unfished recruitmentlevel Ro times random recruitment deviates, which are used to indicate that the population was notat equilibrium at the start of the time series.I used equilibrium spawner per recruit (SPR) quantities to calculate management targets, forillustration purposes I use 40% as a SPR target and use YieldSPR=40% and USPR=40% as targetmanagement benchmarks (Table 3.4 - Equations T4.6 to T4.14). As in all spawner per recruitcalculations, the Yieldtarget and Utarget estimates depend on the selectivity curves calculated foreach year (Equation T4.9).To assess how well the model tracked changes in selectivity over time, I calculated the resultingselectivity estimates by normalizing the yearly vectors of exploitation rate at length (Ul,t) by theyearly average exploitation rate at length (U¯l) (Equation T3.11), which is more stable than the max-imum yearly exploitation rate (maxUl). This happens because observation errors tend to averageout over the length classes, diminishing variability of U¯l in relation to maxUl . When calculatingthe management targets, I used the same method to calculate the mean selectivity at age (EquationT4.7), however I also averaged the selectivity at age curves over the past two years (Equation T4.7)in order to further smooth the curves.The Length-SRA model estimates two main parameters: average unexploited recruitment R0and the recruitment compensation ratio κ . In addition, the annual recruitment deviations wt areestimated for all cohorts observed in the model. That is, the number of recruitment deviationsis equal to the number of years in the time series plus the number of age classes greater thanrecruitment age.40Table 3.1: Indexes, variable definition, and values used in simulation-evaluationSymbol Value Descriptionl {lo, ...,L} Central point of length bin, L = 50 cma {ao, ...,A} Age-class, A = 20 yearst {1, ...,T} Annual time step, T = 50 yearsao 1 First age or age of recruitmentlbin 2 cm Size of length binlo 8 cm Central point of first length bininit 21 Annual time step in which data starts to be re-portedDistribution of length given ageL∞ 50 cm Maximum average lengthK 0.3 Rate of approach to L∞to -0.1 Theoretical time in which length of individuals iszerocvl 0.08 Coefficient of variation for length at age curvePl|a Matrix of proportions of length at ageΦ Standard normal distributionzla,l Normalized z score for lower limit length binszua,l Normalized z score for upper limit length binsbll Lower limit of length binsbul Upper limit of length binsL¯a Mean length at ageσL Standard deviation of length at agePopulation dynamicsRo 100 Average unfished recruitmentκ 10 Goodyear recruitment compensation ratioS 0.7 Natural annual survivalσR 0.6 standard deviation for recruitment deviationswt N (0,σR) Recruitment deviations for years {init-A-ao,...,T}Na,t Numbers of fish at age and timeSBt Spawning biomass at timemata Proportion of mature individuals at agearec,brec Beverton & Holt stock recruitment parametersV Bt Biomass that is vulnerable to the survey at time tva {0,0.5,1,...,1} Survey vulnerability at ageUa,t Exploitation rate at age and timeUl,t Exploitation rate at length and timeCl,t Catch at length and timeNl,t Numbers at length and timelxa Unfished survivorship at ageφe Unfished average spawning biomass per recruitŝell,t Selectivity estimates at length and time41Table 3.2: Variable definition for operating model, MSY quantities, and values used insimulation-evaluationSymbol Value DescriptionOperating modelsell,t Fishing selectivity at length and timeg,d,k vary by scenario Parameters for selectivity functionUt vary by scenario Annual maximum exploitation rateIt Index of abundance at timeσIt 0.1 Standard deviation for index of abundance devi-atesq 1.0 Catchability coefficientτ multivariate logistic error term with στ = 0.1Management quantitieslza Fished survivorship at ageFz seq(0.0,1.0,by=0.001) Hypothetical average fishing mortality to calcu-late management targetsφz Average spawning biomass per recruitφeq Average exploited biomass per recruit under Uzŝela,t Selectivity at age and time tReq Average equilibrium recruitment under UzYieldz Equilibrium yield under UzYieldtarget Yield that would reduce spawner per recruit to40% of unfished levelsUtarget Exploitation rate that reduce spawner per recruitto 40% of unfished levelsThe objective function (Equation T5.8) is composed of a negative log-likelihood component,one penalty, and a prior component for the recruitment compensation ratio κ . The negative log-likelihood component minimizes the differences between the predicted and observed index of abun-dance (Equation T5.1). I assume that such differences are lognormally distributed (Equations T5.3-T5.4) and use the conditional maximum likelihood estimator described by Walters and Ludwig(1994) to estimate the survey catchability coefficient q (Equation T5.2). A lognormal penalty isadded to the negative log-likelihood function to constrain annual recruitment residuals so estimateshave mean of zero and fixed standard deviation σR (Maunder and Deriso, 2003) (Equation T5.5).Lastly, informative normal prior for log(κ) and log(q) were included in the objective function(Equations T5.6 and T5.7).42Table 3.3: population dynamics for Length-SRA and operating modelDistribution of length given agePl|a =∫ zua,lzla,lΦ(z)dz, l < L∫ ∞zla,l Φ(z)dz, l = L∫ zua,l−∞ Φ(z)dz, l = lo(T3.1)zla,l =bll− L¯aσLa(T3.2)zua,l ={bul−L¯aσLal < L1.0 l = L(T3.3)L¯a = L∞ · (1− e(−K·(a−to))) (T3.4)σLa = L¯a · cvl (T3.5)Population dynamicsNa,t>init =arec·SBt−11+brec·SBt−1 · ewt , a = aoNa−1,t−1 ·S · (1−Ua−1,t−1), ao < a < ANa−1,t−1·S·(1−Ua−1,t−1)1−S·(1−Ua,t) , a = A(T3.6)Ua,t =∑l(Pl|a ·Ul,t) (T3.7)Ul,t =Cl,tNl,t(T3.8)Nl,t =∑a(Pl|a ·Na,t) (T3.9)SBt =∑a(mata ·wa ·Na,t) (T3.10)ŝell,t =Ul,tU¯t(T3.11)V Bt =a∑Na,t ·wa (T3.12)Initial year and incidence functionsNa,t=init = lxa ·Ro · e(wt=init ...wt=(init−A+ao)) (T3.13)arec =κφe(T3.14)brec =κ−1Ro ·φe (T3.15)φe =∑alxa ·mata ·wa (T3.16)lxa =1, a = 1lxa−1 ·S, 1 < a < Alxa−1·S1−S , a = A(T3.17)43Table 3.4: Management quantities and operating modelOperating modelNa,t=1 = lxa ·Ro (T4.1)Ul,t =Ut · selOMl,t (T4.2)Cl,t = Nl,t ·Ul,t ·Pl|a · τ (T4.3)sell,t =11−g ·(1−gg)g· ed·g·(k−l)1+ ed·(k−l)(T4.4)It = q ·V Bt · e(N (0,σIt )) (T4.5)Management quantitieslza =lza = 1 a = aolza−1 ·S · exp(−Fz · ŝela−1,t) ao < a < Alza−1·S·exp(−Fz·ŝela−1,t)1−S·exp(−Fz·ŝelA,t)a = A(T4.6)ŝela,t =Ua,t−1U¯t−1+Ua,tU¯t2(T4.7)φz =∑alza ·mata ·wa (T4.8)Targetφ = |φzφe −0.4| (T4.9)φeq =∑alza · (1− exp(−Fz ∗ ŝela,t)) ·wa (T4.10)Req = Ro · κ−φe/φzκ−1 (T4.11)Yieldz = Req ·φeq (T4.12)Yieldtarget = Yieldz→ min(Targetφ ) (T4.13)Utarget = 1− exp(−Fz)→ min(Targetφ ) (T4.14)3.2.2 Simulation-evaluationModel performance was evaluated using a simulation-evaluation with the biological parameters ofan hypothetical fish species. I used the same model structure described in Table 3.3 for both thesimulation and estimation models. However, the operating model was modified to control annualexploitation rate (Equation T4.2), time-varying selectivity (Equation T4.4), and observation andprocess errors.The simulation model was initialized at unfished conditions (Equation T4.1) but only startedreporting data for the simulation-evaluation procedure after the tinit year. Selectivity in the operatingmodel was computed with the three parameter selectivity function described by Thompson (1994)44Table 3.5: Likelihood functions and penaltiesConditional LikelihoodZt = log(It)− log(V Bt) (T5.1)q = eZ¯ (T5.2)Zstatt = Zt − Z¯ (T5.3)LL1 ∼N (Zstat|µ = 0,σ = σIt ) (T5.4)PenaltiesPwt ∼{N (wt |µ = 0,σ = σR) phase < last phaseN (wt |µ = 0,σ = σR ·2) phase = last phase(T5.5)Priorsprior(log(κ))∼N (log(κ),σ = 0.5) (T5.6)prior(log(q))∼N (log(q),σ = 0.5) (T5.7)Objective functionOb j =−log(LL1)+(−log(Pwt ))+ prior(log(κ))+ prior(log(q)) (T5.8)(Equation T4.4). I chose to use this three parameter selectivity curve because of its flexibility,which allowed us to switch between logistic and dome-shaped selectivity curves in the scenariosin which time-varying selectivity was considered. The observation error in the operating modelincluded lognormal error in the index of abundance and logistic multivariate error (Schnute andRichards, 1995) in the catch numbers at length (Table 3.2). Recruitment deviations were assumedto be lognormally distributed with constant σR (Table 3.1).I considered a total of six different scenarios in simulation-evaluation trials, including threehistorical exploitation rate trajectories (contrast, one-way trip and U-ramp) and two selectivitypatterns (constant and time-varying). In the contrast scenario the exploitation rate (Ut) starts lowand increases beyond Utarget and then decreases until Ut ≈ Utarget . In the one-way trip scenarioU increased through time until U ≈ 2 ·Utarget . In the U-ramp scenario, Ut increases steadily untilUt ≈ Utarget and remains constant thereafter. In the constant selectivity scenario, selectivity wasassumed to follow a sigmoid shape. In the time-varying selectivity scenario, the selectivity curvewas assumed to vary every year, progressively changing from a dome shaped curve to sigmoid andback to dome shaped. The complete list of scenarios and the acronyms used are presented in inTable 3.6.45All simulations had 30 years of data, and 200 simulation trials were performed for each sce-nario. I evaluated the distribution of the relative proportional error ( esimated−simulatedsimulated ) for the mainparameter estimates (R0 and κ) and for four derived quantities (Depletion: SBtSB0 , Yieldtarget , Utarget ,and q).Table 3.6: Simulation-estimation scenariosScenario Code Selectivity U trajectoryCC constant contrastCO constant one-way tripCR constant U-rampVC time-varying contrastVO time-varying one-way tripVR time-varying U-ramp3.2.3 Misspecification of growth parametersOne important feature of the Length-SRA is that it assumes that growth follows a von Bertalanffycurve and that the growth parameters are known and constant over time. If this assumption is vi-olated, the model outcomes will be impacted as the model will try to explain the deviations fromthe true growth curve with changes in the selectivity pattern. Here I illustrate how the model out-comes are impacted by the misspecification of the growth parameters by purposefully misreportingthe values of L∞ (Table 3.7). I assumed a simple logistic selectivity curve for this exercise andtherefore expect the model to produce logistic patterns in the exploitation rate at length Ul,t .Table 3.7: Scenarios for testing misspecification of L∞Scenario name version L∞ valuetrue true 68plus10 10% overestimated 74.8minus10 10% underestimated Real data examplesTwo case studies were chosen to illustrate the application of the Length-SRA to real datasets:Pacific hake and Peruvian jack mackerel. Both species are believed to be subject to time-varyingselectivity.46The Pacific hake fishery is believed to exhibit time-varying selectivity due to cohort targetingand annual changes, fleet spatial distribution (Ruttan, 2003). The population is know to have spas-modic recruitment, with high recruitment events occurring once or twice every decade (Ressleret al., 2007). Pacific hake tends to segregate by size during their annual migration (Ressler et al.,2007), allowing the fishing fleet to target strong cohorts by changing the spatial distribution of fish-ing effort as the cohort ages. Pacific hake catch at length data was available for the period between1975 and 2013. The survey index of abundance was available intermittently from 1995 to 2013.The movement pattern of jack mackerel is not as well known, although fish appear to movebetween spawning and feeding areas (Gerlotto et al., 2012). Variability in selectivity patterns forthe jack mackerel fishery are believed to be associated both with evolution of fleet capacity and gearutilization and with compression and expansion of the species range associated with abundancechanges (Gerlotto et al., 2012). Jack mackerel catch at length data was available from 1980 to2013, and the survey index was available between 1986 and 2013, with the exception of 2010.3.3 Results3.3.1 Simulation-evaluationI evaluated the performance of the model in relation to the main parameters, and derived manage-ment quantities with boxplots of the relative proportional error. Throughout, I use the terms positiveand negative median bias to indicate that the median relative proportional error is above or belowzero. The median relative proportional error sign indicate if a parameter has been underestimatedor overestimated the majority of the time.The simulation-evaluation of the Length-SRA model resulted in a small positive median biasfor the κ parameter in all simulations scenarios. Negative median bias was seen in the R0 estimatesin all but one scenario (Constant selectivity and one-way trip exploitation history). The relativeerror estimates for κ indicate that this parameter was nearly unbiased, an effect of the informativeprior considered for that parameter (Figure 3.1 - bottom panel). The estimates for κ are also more47precise than the estimates for R0, this is a result of the use of the informative prior as well as thelikelihood function which lets σR be higher in the last phase of the estimation (Equation T5.5).κRo-1.0 -0.5 0.0 0.5 1.0VRVOVCCRCOCCVRVOVCCRCOCCRelative Proportional ErrorScenarioFigure 3.1: Relative proportional error for main parameters for all scenarios considered in thesimulation-evaluation. Boxplots center lines indicate the median estimate. Lower andupper hinges indicate first and third quartiles. Upper and lower whiskers are given bythe maximum and minimum values within the intervals given by the hinge value +/- 1.5· inter-quartile range (distance between the first and third quartiles).The depletion in the last year of data (SBT/SBo) estimates resulted in negative median relativeerror for all scenarios (Figure 3.2 - top panel). Yieldtarget median relative error were variable beingoverestimates for the VO and CO scenarios and underestimated for the CR, VC and VR scenarios.The absolute median relative error for Yieldtargetwas relatively low (< 7.5%) (Figure 3.2 - secondpanel). The Utarget and q relative error estimates were positively biased with higher median biasesseen for the CR, VO and VR scenarios (Figure 3.2 - third and fourth panels).48qUtargetYieldtargetDepletion(SBT SB0)-1.0 -0.5 0.0 0.5 1.0VRVOVCCRCOCCVRVOVCCRCOCCVRVOVCCRCOCCVRVOVCCRCOCCRelative Proportional ErrorScenarioFigure 3.2: Relative proportional error for main parameters for all scenarios considered in thesimulation-evaluation. Boxplots center lines indicate the median estimate. Lower andupper hinges indicate first and third quartiles. Upper and lower whiskers are given bythe maximum and minimum values within the intervals given by the hinge value +/- 1.5· inter-quartile range (distance between the first and third quartiles).The simulation-evaluation exercise showed that the Length-SRA model is able to track selec-tivity changes through time relatively well (Figure 3.3). However, the selectivity estimates arequite variable, which is likely to be a associated with the observation error in the catch at lengthcomposition.49year: 21 year: 30 year: 40 year: 50scenario: CCscenario: COscenario: CRscenario: VCscenario: VOscenario: VR20 30 40 50 60 20 30 40 50 60 20 30 40 50 60 20 30 40 50 60012301230123012301230123LengthSelectivitysimulatedestimatedFigure 3.3: Simulated and realized selectivity estimates for a set of years within simulation-evaluation time series. The estimatedsolid lines indicate median, 2.5% and 97.5% quantiles for the derived selectivities.503.3.2 Misspecification of growth parametersI found that misspecification of L∞ has severe implications in the capability of the model to estimateexploitation rate at length Ul,t(Figure 3.4). If the value of L∞ was reported to be lower than true,the estimates of Ul,t were lower than true for most length and extremely high for high lengths(approaching the true L∞ ). In the scenario where L∞ was reported to be higher than true, Ul,t wasestimated to follow a dome shaped pattern, with very low exploitation rates for the higher lengths.This patterns occur because the model is trying to adjust the mismatch between proportions ofcatch at length and the Pl|a matrix by changing the predicted selectivity pattern. As a result, failureto adequately specify L∞ leads to erroneous estimation of selectivity patterns and, consequently, tofailure in estimating management quantities.51plus1047plus1048plus1049plus1050minus1047minus1048minus1049minus1050true47true48true49true50182226303438424650545862 182226303438424650545862 182226303438424650545862 182226303438424650545862012340123401234LengthUlengthtypesimulatedestimatedFigure 3.4: Simulated and realized exploitation rate at length Ul,t when L∞ is misspecified. Results shown for the last four years ofsimulation-evaluation time series. Boxplots center lines indicate the median estimate. Lower and upper hinges indicate firstand third quartiles. Upper and lower whiskers are given by the maximum and minimum values within the intervals given bythe hinge value +/- 1.5 · inter-quartile range (distance between the first and third quartiles).523.3.3 Real data examplesThe model fit the Pacific hake and jack mackerel indexes of abundance relatively well (Figure 3.5),despite some limitation in the available data. The Pacific hake index of abundance time-series isrelatively short and intermittent (survey happens every two or three years). The index of abundancetime series for jack mackerel was longer, but it indicates a downward trend in abundance with lowcontrast in the last ten years of data.The model fit for both species resulted in time-varying selectivities that lead to variation inYieldtarget and consequent changes in Utarget (Figure 3.5). This is because changes in selectivity re-sult in changes to the vulnerable biomass even if total biomass is constant. Variability in selectivityand, consequently in Utarget , are more pronounced for Pacific hake, if compared to jack mackerel.The relationship between Yieldtarget and Utarget is also more variable for the Pacific hake case, againan indication of temporal changes in selectivity. Whereas for jack mackerel, both Yieldtarget andUtarget seem to follow the same trend.The selectivity curves estimated for Pacific hake and jack mackerel are quite variable and fre-quently estimated to be dome shaped (Figure 3.6). The resulting selectivity curves presented forthe Pacific hake case, differ from those presented in the 2014 Pacific hake stock assessment (Tay-lor et al., 2014). The selectivities estimates presented in the Taylor et al. (2014) assessment alsovary through time but tend to follow an asymptotic shape. In the present study, the selectivitypatterns tends to alternate between dome-shaped and asymptotic. For the Peruvian jack mackerelcase, the resulting selectivity shapes match those presented in the assessment closely. Both theresults presented here and the results presented in the assessment indicate that the Peruvian fleetselectivity for the Peruvian jack mackerel is dome shaped and with peak selectivity at young ages(Figure 3.6 and Anonymous (2013)). It is important to note, however, that the observed variabilityin selectivity estimates for both examples might indicate real changes in selectivity (e.g., cohorttargeting) or might also be caused by misspecification of the growth parameters (see Figure 3.4).At this point it impossible to determine what are the causes for the resulting patterns in selectivity53observed with the Length-SRA fit. Further investigation would be needed if this model is to beused for management purposes.1231980 1990 2000 2010 Index of abundancePacific hake050001000015000200001980 1990 2000 2010 valuejack mackerel0.5550.5600.5650.5701980 1990 2000 2010 YieldtargetPacific hake1701801902002101980 1990 2000 2010  jack mackerel1002003001980 1990 2000 2010YearCatches in kt02505007501980 1990 2000 2010Year 0.20.31980 1990 2000 2010Year Utarget0.080.090.101980 1990 2000 2010year Figure 3.5: Fit to index of abundance, historical catches, and Yieldtarget and Utarget estimatesfor Pacific hake and jack mackerel. Observed indexes of abundance are shown in opencircles, closed dots in Yieldtarget and Utarget panels indicate model estimates.3.4 DiscussionI present a length-based stock reduction analysis (Length-SRA) that allows monitoring of time-varying selectivity. In the Length-SRA model, catch at length is assumed to be known withouterror, and exploitation rate at length is calculated directly from estimates of numbers at length. Inturn, numbers at length are produced based on numbers at age and on probabilities derived fromgrowth curve parameters and the assumed variability (standard deviation) around mean length atage. This fact is important because it allows the model to bypass the requirement for the estimationof a selectivity ogive, as is required in more traditional age- and length-based models (e.g., Sullivanet al., 1990; Mesnil and Shepherd, 1990) and in more recent length based state-space modelingapproaches (White et al., 2016). Estimation of selectivity ogives can be very difficult, especially541978197919801981198219831984198519861987198819891990199119921993199419951996199719981999200020012002200320042005200620072008200920102011201220132.5 5.0 7.5 10.0 12.5Age (years)YearPacific hake19801981198219831984198519861987198819891990199119921993199419951996199719981999200020012002200320042005200620072008200920102011201220132.5 5.0 7.5 10.0 12.5Age (years)Yearjack mackerelFigure 3.6: Realized selectivity at age patterns across years for Pacific hake and jack mack-erel.55if selectivity is believed to vary over time unpredictably (Martell and Stewart, 2014; Linton andBence, 2011).Nielsen and Berg (2014) presents a stock assessment approach that accounts for time-varyingselectivity by treating fishing mortality at age as stochastic processes that are correlated over ageand time. The accuracy in the estimates of selectivity obtained with the Length-SRA are compara-ble with those presented by Nielsen and Berg (2014), especially for the one-way trip scenarios. TheLength-SRA selectivity estimates are less precise than those shown by Nielsen and Berg (2014),likely because the Length-SRA estimates incorporate observation error. Their model seems to per-form extremely well, however they only considered one exploitation rate trajectory, with significantcontrast in the data. In addition, the changes in selectivity considered in their study are more subtlethan the ones considered here.An important advancement of Length-SRA over conventional stock-assessment models is theindirect calculation of time-varying selectivity. This information alone can be used to characterizethe complexity of the fishery system. Length-SRA on its own is reasonably accurate in deriving im-portant management-oriented parameters (depletion and Yieldtarget), however another option maybe to combine findings from this model with another assessment model, such as a statistical catchat age (SCA) model. In this framework, Length-SRA can be used to calculate annual selectivitypatterns and provide an indication of possible changes over time. These selectivity estimates canthen become an input into an SCA to calculate other important variables and produce managementadvice. This combination of models has been used in the past (Walters and Punt, 1994); I suggestLength-SRA may be a useful tool in this context.Accurate estimates of selectivity are particularly important if the fishery management is basedon yield per recruit reference points. Fishery yield per recruit depends on the selectivity curve(Beverton and Holt, 1957) and for this reason, changes in selectivity over time will directly affectreference points (Beverton and Holt, 1957; Hilborn and Walters, 1992). I observed selectivitychanges for both Pacific hake and jack mackerel, and show how this variability can lead to a notinsignificant difference between the maximum and minimum estimates of Yieldtarget and Utarget56calculated along the time series. I believe that tracking these changes is important not only toensure appropriate management recommendations, but also to illustrate the relationship betweenselectivity patterns and management targets (Vasilakopoulos et al., 2016).One potential point of concern that should be considered when using the Length-SRA is thatit assumes that the biological parameters used in the growth curve and catch at age relationshipare known without error and constant over time. I have tested the Length-SRA under misspeci-fication of the von Bertalanffy growth parameters, and I observed additional bias in the estimatesof parameter and management quantities as well as strong distortions in the resulting selectivityparameters. Similarly, Minte-Vera et al. (2017) showed that misspecification in biological param-eters, especially in asymptotic length, can have a significant impact in assessment results. Otherlength models, e.g., MULTIFAN-CL (Fournier et al., 1998), overcome the assumption of knowngrowth parameters by estimating the von Bertalanffy parameters alongside the assessment param-eters. Once a selectivity curve is assumed, additional deviations in observed catch at length areexplained by adjusting the growth parameters. This assumption can also lead to bias in parameterestimates, as other studies show that variability in selectivity and non-asymptotic patterns are com-mon (Waterhouse et al., 2014). In reality, in most cases it is difficult to know if patterns observedin catch at length are caused by fisheries targeting (i.e., selectivity) or if they would be more appro-priately explained by adjusting the growth parameters. Therefore, I recommend that, when usingthe Length-SRA, the user should perform extensive sensitivity analyses over the possible rangeof values for the growth parameters, particularly if the predicted selectivity patterns are highlyvariable.As mentioned previously, the model and simulation exercise presented here assumes that thegrowth parameters are known and constant through time. Consequently, time variability in growthpatterns could also impacts the results produced by the model. I would not recommend attemptingto estimate time-varying growth parameters within the Length-SRA because growth and exploita-tion rates at length are confounded. However, if estimates of time-varying growth are available,preferably from fishery independent data (empirical age-length keys), those could be used as an57input to the Length-SRA model. Non-stationarity in growth is a relatively easier phenomenon tostudy, particularly if a fishery independent survey is already established. Changes in growth canbe investigated independently by developing empirical age-length keys. In contrast, measuringchanges in selectivity directly is much more difficult and data intensive, requiring independenttagging programs.The approach used in the Length-SRA is analogous to that used in virtual population analysisin that the length composition data is assumed to be known without error. For this reason, theselectivity estimates include extra variability due to observation and sampling error. I attempted tominimize this effect by smoothing the predicted selectivity over two years, however this methodis not capable of completely removing the observation error effect from the selectivity estimates.Because of the assumption of known catch at length, it is important that the catch sampling isrepresentative of the total removals from the population (Pope, 1972). As in any other fisheriesmodel, biased sampling and/or low sampling effort will result in bias in parameter and fisheryreference point estimates (Coggins and Quinn, 1998; Bunch et al., 2013).Some management parameters are consistently overestimated (Yieldtarget) and underestimated(depletion), which may be cause for concern. However, it is important to note that both parametershave low absolute median relative error (<7%). The magnitude of the bias in the estimates ofYieldtarget and Utarget observed in this study are comparable (in magnitude) to the results obtainedby Martell and Stewart (2014) for MSY and FMSY in a simulation study on the impacts of time-varying selectivity on the estimates generated by a statistical catch at age model. Other studies showeven higher biases in face of time-varying selectivity (e.g., Linton and Bence, 2011; Henríquezet al., 2016). The estimates of depletion are also comparable to those produced with other SRA typeassessments evaluated by Thorson and Cope (2015). Overall, parameter and derived parametersestimates are generally within the range of many other stock assessment models.The Length-SRA approach presented in this study can be a useful tool for fisheries stock as-sessment. I believe that this is particularly true when time-varying selectivity is thought to occur,especially if the variability is not easily predictable from historical changes in gear use/fleet com-58position. However, I would like to acknowledge that the selectivity estimates will only be reliable ifthe growth parameters for the population being assessed are known. In addition, the simple natureof the Length-SRA model makes it a good candidate model for inclusion in closed-loop simulationstudies. Further testing of this model in a closed-loop simulation set up would provide more in-sight on the model performance on achieving management outcomes (Punt et al., 2016). I foreseethe application of this model as an investigative tool to evaluate potential time-varying selectivitypatterns, as a stock assessment tool and as part of closed loop simulation studies.59Chapter 4Evaluation of harvest control rules fortransboundary stocks4.1 IntroductionIn fisheries management, harvest control rules are previously agreed upon management actionsthat should be taken in response to stock status indicators (Deroba and Bence, 2008). Despite thelong history of harvest control rule evaluations in the literature (Hilborn, 1986; Walters and Parma,1996), the practice of formally adopting harvest control rules in fisheries management is morerecent (Deroba and Bence, 2008). Over the past decade, harvest control rules have been recognizedas a mechanism to increase consistency and transparency in fisheries management (Punt, 2010;Kvamsdal et al., 2016). In addition a variety of policy documents have fostered the adoption ofharvest control rules in many fisheries around the world (Government of Canada, 2009; Food andAgriculture Organization of the United Nations, 1996; Magnuson-Stevens act, 2007).Harvest control rules can be generally grouped into three categories: fixed escapement, fixedexploitation rates and threshold harvest control rules (Punt, 2010). Fixed escapement rules im-ply that fishing should only occur if the biomass is higher than a given reference biomass, i.e.,a given amount of exploitable biomass is allowed to “escape” harvest. This kind of harvest con-trol rule is more commonly applied to salmonids (Hawkshaw and Walters, 2015), and it has been60shown to lead to high variability in yields and increased frequency of closures (Deroba and Bence,2008). The fixed exploitation rate harvest control rules works by adjusting the yield in propor-tion to the population size. Fixed exploitation rate harvest control rules reduce inter annual fisheryvariability when compared to fixed escapement harvest control rules (Hilborn, 1986). Finally thethreshold harvest control rules are usually designed to produce stepwise changes in exploitationrate as biomass decreases below some threshold; many such control rules also include a minimumescapement threshold i.e., the harvest rate is set to zero if biomass is below a limit level. Thresholdharvest control rules are popular among many fisheries management agencies (Punt, 2010) becausethey tend to maintain higher biomass and allow for faster rebuilding of depleted stocks (Quinn IIet al., 1990).Ideally, the choice of harvest control rule is based on an evaluation process that optimizes a setof performance metrics associated with the objectives of a fishery. Some examples of commonlyused performance metrics include average expected yield, annual average variability in yield andminimum biomass threshold (Punt et al., 2016). This process usually also exposes the trade-offsbetween potentially conflicting objectives (Hall et al., 1988; Punt and Donovan, 2007). The evalu-ation of harvest control rules is usually done through closed-loop simulations (e.g., Walters, 1998;Ishimura et al., 2005) which are computer models used to simulate the fisheries management sys-tem and evaluate the performance of management options given the best available understandingof dynamic of the resource as well as observation and implementation error models.The performance of a harvest control rule is often evaluated for an entire stock (e.g., Ishimuraet al., 2005; Tong et al., 2014; Hawkshaw and Walters, 2015). However, in the case of transbound-ary stocks the aggregate evaluation of performance metrics may not reflect the outcomes experi-enced by each nation separately. This is especially true if the distribution of the fished stock variessystematically with abundance. For example, many stocks exhibit range contraction as abundancedecreases (e.g., Brodie et al., 1998). This can affect the availability of the resource in a given areaeven when the biomass is considered to be above the aggregate stock reference points. In otherstocks, spatial distribution may be a function of ontogeny, usually with larger individuals perform-61ing more extensive migrations (e.g., Ressler et al., 2007). These populations too can become lessavailable to a given fishing nation if the population age or size structure becomes truncated due tofishing mortality, even if such fishing mortality is equal to target exploitation goals for the aggregatestock. These effects are usually not accounted for in international fisheries management treaties de-spite the fact that treaties are often designed with the intention of securing equitable benefits of theresource to all the parties. When the application of a harvest control rule results in a change in thedistribution of the stock, then it is possible that the benefits of the treaty will not be realized by oneor more parties. In such situations, it may be useful to consider the relationship between biomass,age/size composition and spatial distribution of the stock when managing transboundary stocks.This would enable us to address the following questions: Is it possible to optimize a set of perfor-mance metrics across all nations that share the resource? If not, what are the trade-offs betweenthe performance metrics for each of the nations sharing the resource? The explicit consideration ofthese trade-offs can help to design effective management strategies for shared stocks. This couldbe achieved by using spatially explicit models to identify potential differences in management out-come experienced by nations separately, thus, exposing trade-offs between two nations that wouldnot otherwise have been visible using a spatially aggregate approach.Pacific hake (Merluccius productus) is an example of a transboundary stock whose spatial dis-tribution is thought to be affected by changes in age structure (Bailey et al., 1982). It exhibitsseasonal migratory behavior with spawning occurring off southern California during the winterand fish migrating north between spring and fall to feed (Ressler et al., 2007). Larger fish, typicallyolder than age-4, migrate longer distances and are found to be more abundant in Canadian waters(Methot and Dorn, 1995). Fish age-3 years and younger tend to remain in U.S. waters off the coastof California and Oregon (Methot and Dorn, 1995; Ressler et al., 2007).Management of the Pacific hake stock follows the regulations determined in an internationaltreaty between Canada and the U.S.A. (United States State Department, 2004). The treaty estab-lishes that the coastwide Total Allowable Catches (TAC) should be calculated following a 40:10threshold harvest control rule. The treaty also determines that the TAC be split between the two62countries following a fixed allocation, 73.88% to the U.S.A. and 26.12% to Canada. Given thespatial dynamics of the resource and location of the spawning grounds, the American fishing fleethas first access to the incoming cohorts. This has the potential to cause conflicts between the twonations because harvest will lead to age truncation, which in turn could lead to limited availabilityof the resource in Canadian waters. Despite the potential for conflict associated with harvest levelsand spatial distribution, the impacts of the current treaty based management procedures have onlybeen evaluated for the aggregate stock (Ishimura et al., 2005; Punt et al., 2008; Taylor et al., 2014;Hicks et al., 2016). However, given the life history of the stock and the fishing practices (Baileyet al., 1982; Ressler et al., 2007), it is important to consider spatial effect of management strategiesthat are otherwise masked by non-spatial models.In this chapter, I aim at addressing how the application of a harvest control rule to an aggregatestock affects spatial fishing opportunities. I evaluate the performance of a large set of harvestcontrol rules for the Pacific hake stock using a spatially explicit model in a closed-loop simulationroutine. I illustrate some differences in performance between harvest control rules using yield andconservation related metrics, for the whole stock and relative to each individual nation. I map thetradeoffs between alternative harvest policies for the two nations sharing the Pacific hake resource.4.2 MethodsIn order to evaluate the performance of various harvest control rules, I performed a series of closed-loop simulations. The spatial operating model was used to simulate fishery and scientific data everyyear. Observations on the stock status were generated with process and observation error. Annualtotal allowable catch was set each year using alternative harvest control rules and the data generatedfrom the spatial operating model. In the following sections, I describe the simulation model, theharvest control rules, and performance metrics in detail.634.2.1 The simulation modelI used the spatial model described in Chapter 2 as an operating model to describe the populationdynamics of the Pacific hake resource. The Lagrangian model used in that study allows the fish tomove in a seasonally cyclic manner that is characteristic of the Pacific hake offshore stock (Ressleret al., 2007). This model applies a sine function to model the cyclic the movement of each indi-vidual cohort along the coast of U.S.A. and Canada. The mean position of an individual cohort ata given time step t is given by X¯a,t , which is a function of the mean minimum position X¯min, themean maximum cohort specific position X¯max,a, the number of time steps within a migration cycletmax and the time step at which the migration cycle starts t0 (Equation 4.1).X¯a,t = X¯min +(X¯max,a− X¯min) ·(0.5+0.5∗ sin(t · 2pitmax− t0 · 2pitmax −pi2))(4.1)As noted above, the mean maximum position X¯max,a is specific for each separate cohort. Theextent of the migrations is given by a logistic function of age, similar to the one used in Methot andDorn (1995) (Equation 4.2), which allows older (and larger) fish to move further away from thespawning grounds. In the following function, the parameters a50 and σXmax are the logistic functionparameters, and σvt is normally distributed random error component.X¯max,a =11+ exp(−(a−a50)/σXmax)· e(vt∼N (0,σvt)) (4.2)For a more detailed description of the movement model, as well as the population dynamicscomponents and effort dynamics components of the operating model used in this study, I recom-mend that the reader refer to Chapter 2. In this study, I used the multiple groups version of the modeldescribed in Chapter 2. The model parameterization was extracted from the 2017 stock assessment(Berger et al., 2017b), and the movement parameters were set to approximate the movement dy-namics of Pacific hake described in the literature (Methot and Dorn, 1995; Ressler et al., 2007). All64parameter values are given in Table 4.2. In order to mimic the historical trends in abundance, I setthe historical catch limits for the historical data equal to the realized catches by each nation.The operating model, therefore has the ability to model Pacific hake movement and effort dis-tribution allowing the characterization of the relationship between fishing mortality, strong recruit-ment events and cohort targeting on the distribution of the stock, and hence the realized perfor-mance in each nation. As the fishing mortality increases, the age structure of the population willtend to become truncated (due to the cumulative effects of fishing mortality), and therefore, in theoperating model, the mean biomass distribution of the stock will shift southwards, reducing theavailability of the resource in northern waters.The Pacific hake stock distribution is also affected by the recruitment dynamics. When a strongrecruitment event occurs, the mean biomass distribution of the resource will be strongly linked tothe age of that strong cohort, i.e., the mean biomass distribution of the stock will tend to shift northas that cohort grows. Strong recruitment events occur recurrently for the Pacific hake offshorestock (Berger et al., 2017b; Ressler et al., 2007). These strong recruitment events usually increasethe overall biomass of the stock significantly, and individuals from strong cohorts dominate thefisheries catch for a few years. Despite the intermittent strong recruitment events in the Pacific haketime series, the causes for such strong recruitment events are unknown so predictions of when astrong recruitment will occur in the future are unreliable. In order to model this uncertainty, I optedto simulate three possible recruitment scenarios based on the historical recruitment time seriesreported by Berger et al. (2017b). I repeated the last 30 years of recruitment deviations reportedby Berger et al. (2017b) twice to generate the recruitment for the next 60 projection years. Theseprojection recruitment time series were modified to generate three scenarios: no strong recruitmentevents, one strong recruitment event per decade, and two strong recruitment events per decade. Iopted for the repetition of the historical recruitment time series in order to maintain the historicalpatterns in recruitment autocorrelation, as well as the cycles in strong recruitment occurrences, i.e.,a minimum interval of thee years between strong recruitment events.65llllllllllll ll l l l l l l l5 10 15 20354045505560AgeX max, a (Latitude)Basea50 = 4, σXmax = 2lllllll l l l l l l lEarly movementa50 = 2.5, σXmax = 1lllllllllllllll ll l llLate movementa50 = 5.5, σXmax = 3CanadaU.S.AFigure 4.1: Logistic functions used in the three movement scenarios considered in this study.In addition to scenarios that explore the impacts of strong recruitment events, I also exploredthe sensitivity to the model in relation to the movement parameters. More specifically, I exploredchanges in the logistic function that describes the maximum average position reached by each co-hort (Equation 4.2). The modification to the logistic curve are shown in Figure 4.1. The movementparameter sensitivity analysis was only performed for the no strong recruitment scenario. A list ofscenarios evaluated is presented on Table 4.1.The historical population was reconstructed for the years 1966 to 2016 by simulating the popu-lation dynamics using the stock assessment parameters reported in Berger et al. (2017b), and settingthe historical catches equal to those extracted by the U.S.A. and Canada fleets during that period.Catches for both nations for the historical period were also reported by Berger et al. (2017b). Then66Table 4.1: List of scenarios for closed loop simulationsNumber Recruitment Movement1 no strong recruitment Base2 one strong recruitment per decade Base3 two strong recruitment per decade Basem2 no strong recruitment Early movementm3 no strong recruitment Late movementthe closed loop simulations were carried on for 60 years into the future. For the closed loop simu-lations, I assumed that the TAC allocation between nations remains constant as determined in thePacific hake treaty, 73.88% to the U.S.A. and 26.12% to Canada.In the simulations, the observation model was represented by adding uncertainty around themodel predictions of relative spawning biomass, i.e., spawning biomass levels in relation to the un-fished average. This study focuses on exploring and comparing the performance of harvest controlrules, not on the stock assessment components of the management process. Therefore, I chose tonot implement an assessment model in the closed loop simulations. Instead, the observation uncer-tainty was represented by adding autocorrelated and normally distributed error around the modelpredictions of true biomass and spawning biomass (Equation 4.3 and 4.4). Walters (2004) showsthat estimates of biomass from stock assessments are usually autocorrelated over time; I assumedthat the autocorrelation coefficient, ρ , was 0.5. This approach was chosen because of computa-tional convenience, i.e., faster simulation running time than the simulation of the stock assessmentmethodology. I believe that the high autocorrelation coefficient I use provides a comparable effectto that of using more complicated assessment models.B̂t = Bt · eεt (4.3)ŜBtSB0=SBtSB0· eεt (4.4)εt = ρ · εt−1 +νt (4.5)νt ∼N (0,σ = 0.3) (4.6)67Table 4.2: Pacific hake operating model dimensions and parameter valuesSymbol Value or Range DescriptionModel Dimensionst 1−12 Time steps within a migration cycley 50−160 YearsY 60 Total Number of projection yearshy 50−100 Historical yearsa 1−20 Ager 30−60 Areak 5 Fishing groundskb 42, 46, 48.5 and 51 Fishing ground boundaries in latitude de-greesn 2 Number of nations48.5 Nation boundary in latitude degreesg 1−20 Groupsdr 1 Interval between two adjacent areasPopulation dynamics parametersM 0.223 Annual natural mortalityR0 2923 thousand tons Average unfished recruitmenth 0.814 Beverton & Holt recruitment steepnessσR 1.4 Standard deviation for recruitment devia-tions - used in bias correction onlyMovement parameterst0 1 Time step at which individuals are at theirminimum average positionCV 0.07 Coefficient of variation for position at atagea50 4.0 Inflection point for maximum average posi-tion logistic functionσXmax 2.0 Standard deviation for maximum averageposition logistic functionerror levelsσwx 0.08 Standard deviation for lognormal variationaround the maximum average positionσvt 0.1 Standard deviation for lognormal variationaround the effort scalerEffort parametersEy,n 1 for nation 1 and 0.35 for nation 2 Yearly effort scaler - constant for all yearsEt,k (0,0,0,0,0.5,1.0,1.0,1.0,0.5,0.1,0.0,0.0)for k = 1,2,3 and(0,0,0,0,0,1.0,1.0,1.0,0.5,0.3,0,0)for k = 4,5Monthly effort scalerq 3 Effort scaler684.2.2 The 40:10 harvest control ruleThe Pacific hake treaty determines that the annual quota for the stock should be determined accord-ing to a 40:10 threshold harvest control rule (United States State Department, 2004; Hicks et al.,2016). The 40:10 harvest control rule used in the Pacific hake agreement determines that the coastwide total allowable catches are calculated using a proxy for the fishing mortality that will pro-duce maximum sustainable yield (FMSY - in practice replaced by the proxy FSPR=40%) whenever thespawning biomass is above 40% of average unfished levels. The 40:10 adjustment refers to the re-duction in harvest rate when the spawning biomass falls below 40% of unfished average equilibriumlevel. The harvest rate adjustment corresponds to a linear decrease in TAC if the spawning biomassis between 40% and 10% of unfished spawning biomass and set to zero if spawning biomass isbelow the 10% threshold. I used the same formulation as the one described by Hicks et al. (2016).The 40:10 harvest control rule results in very high quota recommendations when the stock abun-dance is high, which happens whenever a strong cohort reaches maturity (around age 3). However,in reality the actual quotas recommended by the Pacific hake Joint Management Committee areusually capped and tend not to exceed 600 thousand metric tonnes (Hicks et al., 2016). For thisreason, I implemented a cap on the 40:10 harvest control rules recommendations to 600 thousandtonnes.4.2.3 Linear harvest control-rulesI evaluate a series of harvest control rules that are given by a linear relationship between relativespawning biomass (SBt/SB0) and TAC. I opted for evaluating linear harvest control rules becausethey encompass a broad range of harvest control rules commonly used in fisheries management(i.e., fixed exploitation rate and fixed escapement). Similar rules are frequently considered in policyoptimization studies, (e.g., Reed, 1979; Moxnes, 2003; Hawkshaw and Walters, 2015). In addition,this type of control rule has also been previously evaluated for Pacific hake (Ishimura et al., 2005)but not within a spatially explicit model. The re-evaluation of the linear harvest control rules allowfor direct comparison with the Ishimura et al. (2005) study, building up on previous knowledge.69The linear harvest control rule functional form is given in Equation 4.7 and illustrated on Fig-ure 4.2. The slope is the harvest rate at which the stock biomass is harvested. The intercept isthe relative spawning biomass threshold (biomass threshold, for short). The biomass threshold isthe minimum relative spawning stock biomass level required for harvest to occur. B̂t is the ob-served total biomass in the last year of data and SB0 is the unfished average spawning biomass.Analogously to what was implemented for the 40:10 harvest control rule described in the previoussection, I impose a cap on TAC so that it does not exceed 600 thousand tonnes. Figure 4.2 shows anillustration of the resulting TACs as a function of relative spawning biomass for two linear harvestcontrol rules examples and the implementation of 40:10 rule.TAC = slope · B̂t · (ŜBt−SB0 · intercept)SB0 (4.7)To evaluate and compare the performance of linear harvest control rules, I systematically cal-culated a series of performance metrics over a range of 10 slopes (to capture alternative harvestrated from 0.05 to 0.5 in 0.05 intervals) and 6 intercepts (representing different relative biomassthreshold from 0.0 to 0.5 in 0.1 intervals). I then mapped the performance metrics values over theslope-intercept surface. I computed the performance metrics for 54 combinations of slope and in-tercept and the 40:10 harvest control rule. Each combination is evaluated with 100 simulation runswith the same set of random number seeds for each harvest control rule, which was sufficient toaccurately and precisely characterize the distribution for each performance metric for each harvestcontrol rules.4.2.4 Performance metricsIn order to evaluate the harvest control rules, I calculate a set of performance metrics for the aggre-gate fisheries and for each nation’s fleet separately. The performance metrics that were calculatedfor each nation separately include log utility (Average of annual log of yield plus a small value),total yield and annual average variability in yield (AAV) (Table 4.3). Total yield and AAV were700.0 0.2 0.4 0.6 0.8 1.00100200300400500600SBt SB0TACslope = 0.2intercept = 0.0slope = 0.5intercept = 0.440:10 ruleFigure 4.2: Illustration of the relationship between relative spawning biomass and TAC fortwo example linear harvest control rules and the 40:10 harvest control rule used in thisstudy.chosen to illustrate potential impacts of the magnitude and variability of yield, respectively. Logutility was included to represent a composite view of these two quantities as it increases as yieldincreases, but it also strongly penalizes fisheries closures (i.e., yield equal to zero). In addition,I computed the following performance metrics for the stock as an aggregate (both nations com-bined), % of times the fishery closed, % of times the total biomass was below 40% of averageunfished levels (Table 4.3).This set of performance metrics was chosen based of potential objectives that have been ana-lyzed in other closed loop simulations for hake (Taylor et al., 2014; Hicks et al., 2016) and other71Table 4.3: Equations used to calculate performance metrics for evaluation of harvest controlrules. All quantities were averaged across simulation runs.Average log utilitylog(U) =∑y log(Yieldy +1)Y(T3.1)Average annual yieldYield =∑yYieldyY(T3.2)% of closures%Closure =y∑CL/Y (T3.3)Mean Annual Average Variability in YieldAAV = mean(|Yieldi−Yieldi−1|∑i=yi=y−1Yieldi) (T3.4)% of years with of stock spawning biomass below 40% of SBo%B < 40% =∑y B < 40%SBoY(T3.5)fisheries (Cox et al., 2013; Punt et al., 2008): preference for higher average yield, stability (i.e.,low inter annual variation in TAC), avoidance of fisheries closures, and the conservation objectiveof maintaining the biomass at or above the target level.The analysis presented in this study include 305 unique simulation configurations (61 harvestcontrol rules, three recruitment scenarios and three movement scenarios). In order to summarize theoutputs I chose to compute the average performance over the projection years and the simulationruns for each harvest control rule, recruitment scenario and movement scenario configuration.For both yield and log utility, I compare the relative performance by nation by using trade-offplots, i.e., by contrasting the relative performance of each metric by nation. The average perfor-mance metrics for each nation are rescaled so that the maximum value is set to one. For these trade-offs plots, I limit the linear harvest control rules to those with intercept (biomass threshold) equalor below 0.3. This threshold was chosen so that the figures would be more easily interpretable.To compare the performance of the linear and 40:10 harvest control rules in terms of AAV,percentage of times the fishery closed, and percentage of times the total biomass was below 40%, Icreated surface maps where the x and y axes represent the values for biomass threshold and harvest72rate for the linear harvest control rules, and the z axis (color gradient) represents the values obtainedfor each performance metric. I also indicated in the surface, which of the linear harvest control rulesproduced the minimum or maximum (the optimum) performance metric. The performance metricvalue obtained with the 40:10 rule was also mapped to that surface, i.e., the 40:10 indicator wasplaced on the point that more closely matched the performance metric value obtained with the40:10 rule.4.3 ResultsEach combination of harvest control rule and recruitment scenario gives rise to a biomass and yieldtrajectory. As an example, I show three of these catch trajectories under the scenario with no strongrecruitment events (Figure 4.3 ). These trajectories represent median and 95% intervals for catchtrajectories simulated under the 40:10 rule and for two linear harvest control rules with 0.1 biomassthreshold and 0.1 and 0.5 harvest rates. The linear harvest control rule with 0.1 harvest rate (redline, Figure 4.3 ) produces lower catches than the 40:10 rule for both nations. However, the linearharvest control rule with higher harvest rate (0.5 - blue line, Figure 4.3) produces higher medianyields than the 40:10 rule for nation 1 (U.S.A.) and lower median yield than the 40:10 rule fornation 2 (Canada). This result is an indication of how high harvest rates decrease the availabilityof the resource in northern areas, leading to lower yield to Canada. One interesting fact to notice isthat, in Figure 4.3, the yield variability increases as the exploitation rate increases. This result couldbe considered counterintuitive as one would expect less variability in catches if a stock is lightlyfished. This result is likely a byproduct of the simulation design. Here, I assume that the populationhas exactly the same recruitment deviation trajectories regardless of stock size for all simulationruns. Therefore the variability around biomass and yield becomes directly proportional to thepopulation size. High biomass, will produce higher recruitment variability and higher observationerror, consequently producing higher yield variability, while the inverse will occur for low biomass.7301002003004002000 2040 2080Year Yield (1000t) exploitation rate0.10.540:10 ruleNation 10501002000 2040 2080 Yield (1000t) exploitation rate0.10.540:10 ruleNation 2Figure 4.3: Historical and median and 95% intervals for projected catches for three exam-ple harvest control rules under the "no strong recruitment" scenario. Harvest controlrules include the 40:10 rule with cap, and two linear harvest control rules with biomassthreshold (SBt/SBo) of 0.1 and exploitation rates of 0.1 and 0.5.Reproducing Figure 4.3 for all harvest control rules evaluated in this study would be a formidabletask. For this reason, I chose to summarize the results for all scenarios and harvest control rules bycomparing the mean performance across projection years and across simulation runs.In Figure 4.4, I show the differences in mean performance when it comes to average log util-ity and average yield for the two nations sharing the resource. If no interaction between harvestlevels and spatial distribution of the stock existed the lines shown in Figure 4.4 would be straight,following a 1:1 ratio. However, the lines tend to bend for low biomass threshold and higher harvestrates, indicating that Canada (Nation 2) experiences lower than expected catches. This is a result ofdecreased availability of the Pacific hake resource in Canadian waters. As the harvest rate (slope)increases, the population age structure will become more truncated, and it becomes harder for the74Canadian fleet to catch its full quota. In general the difference in performance between nationsand harvest control rules become more important as the incidence of strong recruitment decreases(4.4). In terms of log utility, the 40:10 harvest control rule tends to perform near optimum for bothnations across all scenarios (Figure 4.4). However, when it comes to average yield, the 40:10 har-vest control rule performs relatively better for Canada than the U.S.A. (Figure 4.4). For the U.S.A(Nation 1 - x axis on Figure 4.4), both yield and log utility tend to increase as harvest rate increasesand biomass threshold decreases. For Canada (Nation 2 - y axis on Figure 4.4) results for log utilityand yield indicate that better performance is obtained at harvest rates between 0.1 and 0.25 for theno strong recruitment scenario and higher harvest rates (> 0.35) for the other recruitment scenarios(Figure 4.4). In addition, for the no strong recruitment scenario, the relative yield and log utilityfor Canada declines when harvest rates are greater than 0.35, while the relative yield and log utilityresults for the U.S.A. tend to remain constant (Figure 4.4).7540:10 40:10 40:10no strong recruitment one strong recruitment two strong recruitments0.7 0.8 0.9 1.0 0.6 0.7 0.8 0.9 1.0 0.6 0.7 0.8 0.9 log utility − Nation 1average log utility − Nation 2no strong recruitment one strong recruitment two strong recruitments0.4 0.6 0.8 1.0 0.25 0.50 0.75 1.00 0.25 0.50 0.75 average yield − Nation 1 average yield − Nation 2harvest rate0. threshold 0 0.1 0.2 0.3 40:10Figure 4.4: Comparison of relative performance (normalized to the maximum value by nation) of average log yield (log utility) andaverage yield between nations (trade-offs). Quantities were normalized by the maximum observation per nation and scenariofor comparison. Nation 1 = U.S.A and Nation 2 = Canada. Colors indicate harvest rates for linear harvest control rules andpoint shape indicate biomass threshold relative to unfished spawning biomass (SBt/SBo). Text indicates performance of the40:10 harvest control rule.76The AAV values increased as the incidence of strong recruitment events decreased (see colorscale legend on Figure 4.5). Across all recruitment scenarios, the minimum variability occursat lower values of biomass threshold (≤ 0.1), and high harvest rates (> 0.4) (Figure 4.5). The40:10 rule AAV gets closer to the point where minimum AAV is found as the incidence of highrecruitment events increase (Figure 4.5). This results indicate that, for the scenarios consideredin this study, minimum variability in catch is obtained when the stock is fished to very low levelsand kept near the origin of stock-recruitment relationship. At some point, in reality, the biomasswould get so low that it would no longer be economically viable to operate a fishery. Similarlyto what happened in Figure 4.3, the results observed in Figure 4.5 could also be influenced bythe simulation design. Lower variability in both biomass and yield occur at low population size,because recruitment deviations here were assumed to follow the same trajectory irrespective ofpopulation size. This results would likely change, if, for example, strong recruitment events becamemore likely when abundance is either high or low. To date, there is no evidence of such phenomenonfor the Pacific hake stock, as historical strong recruitment events seem to happen at both high andlow spawning biomass levels (Berger et al., 2017b).Fisheries closures were only computed for the aggregate stock, as closures are determined basedon the total biomass threshold. As expected, closure rates decrease as biomass threshold decreases(Figure 4.6). Also, the area with very low closure percentage (blue area in the colored surface)increases as the occurrence of strong recruitment events increases. This happens because whenstrong recruitment events occur the population biomass tends to increase considerably, rendering itunlikely that the biomass levels will fall below the threshold. The 40:10 also performed very wellin avoiding fisheries closures, with performance comparable to the minimum closure percentageacross all scenarios (Figure 4.6).With reference to the conservation performance metric, the % of time that biomass is below40% of unfished levels, minimum values occurred for high biomass threshold values (Figure 4.7 -blue area). For the scenario with two strong recruitment events per decade low harvest rates alsoproduced low % of biomass below 40%. In relation to the 40:10 rule performance, the difference77between the 40:10 rule and the minimum % of time that biomass is below 40% of unfished levelstended to increase as the number of strong recruitments per decade decreased (Figure 4.7). For thescenarios where no strong recruitments occur, biomass was predicted to be below 40% of unfishedlevels about 86% of the time.The sensitivity analysis for the movement parameters of the migration model indicate that themovement parameter values have less impact on the difference in performance between nationsthan the occurrence of strong recruitments (Figure 4.8). However, the movement parameters dohave an impact. For the Early movement scenario, i.e., where fish migrate farther at younger ages,the differences in performance between the nations are less prominent than in the Base case andLate movement scenarios. In addition, in the Early movement scenario, the 40:10 rule performanceis closer to that of linear harvest control rules with low harvest rate and low biomass threshold forboth nations (Figure 4.8). This result indicates that if a higher percentage of younger fish migratefurther north, there will be less of an impact on the long term yield due to reduced availability ofthe resource in the northern range of the distribution.7840: 0.2 0.4Biomass thresholdHarvest rateNation 140: 0.2 0.4Biomass ThresholdNation 240: 0.2 0.4Biomass threshold0. AAVTotalno strong recruitment40: 0.2 0.4Biomass thresholdHarvest rateNation 140: 0.2 0.4Biomass ThresholdNation 240: 0.2 0.4Biomass threshold0.20.30.4Median AAVTotalone strong recruitment40: 0.2 0.4Biomass thresholdHarvest rateNation 140: 0.2 0.4Biomass ThresholdNation 240: 0.2 0.4Biomass threshold0. AAVTotaltwo strong recruitmentsFigure 4.5: Median Average Annual Variability in yield (AAV). Surfaces indicate perfor-mance for linear harvest control rules. White circle indicates performance level compa-rable to the 40:10 harvest control rule. Minimum value along the surface is indicatedwith black square. Recruitment scenarios are indicated on the rows. Nation 1 = U.S.Aand Nation 2 = Canada.Total is aggregate result.7940:100.250.500.750.0 0.2 0.4Biomass thresholdHarvest rate0 0.6525Closure rateTotalA40:100.250.500.750.0 0.2 0.4Biomass thresholdHarvest rate0 0.3622Closure rateTotalB40:100.250.500.750.0 0.2 0.4Biomass thresholdHarvest rate0 0.1852Closure rateTotalCFigure 4.6: % closures for total fisheries over the 60 years of simulation. Surfaces indicateperformance for linear harvest control rules. White circle indicates performance levelcomparable to the 40:10 harvest control rule. Minimum value along the surface is indi-cated with a black square. Recruitment scenarios are indicated on the columns: A - nostrong recruitments, B - one strong recruitment per decade,C - two strong recruitmentsper decade.40: 0.2 0.4Biomass thresholdHarvest rate0.246 0.872% time Bt < 0.4 BoTotalA40: 0.2 0.4Biomass thresholdHarvest rate0 0.633% time Bt < 0.4 BoTotalB40: 0.2 0.4Biomass thresholdHarvest rate0 0.051% time Bt < 0.4 BoTotalCFigure 4.7: % of time that biomass is above 40% of average unfished levels. Surfaces indicateperformance for linear harvest control rules. White circle indicates performance levelcomparable to the 40:10 harvest control rule. Minimum value along the surface is indi-cated with a black square. Recruitment scenarios are indicated on the columns: A - nostrong recruitments, B - one strong recruitment per decade,C - two strong recruitmentsper decade.8040:10 40:10 40:10Base case movement Early movement Late movement0.7 0.8 0.9 1.0 0.7 0.8 0.9 1.0 0.7 0.8 0.9 log utility − Nation 1average log utility − Nation 2Base case movement Early movement Late movement0.4 0.6 0.8 1.0 0.4 0.6 0.8 1.0 0.4 0.6 0.8 average yield − Nation 1 average yield − Nation 2harvest rate0. threshold 0 0.1 0.2 0.3 40:10Figure 4.8: Comparison of relative performance of average log yield (log utility) and average yield between nations (trade-offs) forthe three alternative movement scenarios. Quantities were normalized by the maximum observation per nation and scenariofor comparison. Nation 1 = U.S.A and Nation 2 = Canada. Colors indicate harvest rates for linear harvest control rules andpoint shape indicate biomass threshold relative to unfished spawning biomass (SBt/SBo). Text indicates performance of the40:10 harvest control rule.814.4 DiscussionThe results shown in this study corroborate those presented by Ishimura et al. (2005) in that lowervalues of both biomass thresholds and harvest rates tend to produce higher yields and lower vari-ability in yield for the aggregate stock. The present results are also in agreement with Ishimura et al.(2005) in pointing that the 40:10 harvest control rule performs similarly to the low harvest rate andlow biomass threshold alternatives in terms of maximizing yield and minimizing yield variability(bottom left are of the heat maps). The results obtained for the 40:10 harvest control rule withcap also corroborate with those presented Taylor et al. (2014) and Hicks et al. (2016), indicatingthat the 40:10 harvest control rule, with cap, performs well in the long term; i.e., it secures yieldwhile minimizing fisheries closures and maintaining the biomass above the 10% of unfished levelsthreshold.However, in addition to evaluating the performance of harvest control rules for the stock as anaggregate, I also evaluated the impacts experienced by each nation individually. I show that there isa difference in performance between the two fishing nations when the log utility and sum of yieldsmetrics are considered. In general, after rescaling by the maximum observed values, the U.S.A. hashigher relative yield and log utility when fishing at higher harvest rates when compared to Canada.These differences diminish as the frequency of strong recruitment events increase. This result isassociated with the migration of the stock because the Pacific hake migration is associated with ageand size of the fish (Ressler et al., 2007). As harvest rates increase, the overall age structure of thepopulation will become more truncated, i.e., older and larger fish become less abundant; thereforethe resource becomes less abundant in the northern range of the distribution which correspondsto waters off Washington state and Canada. This effect is minimized as the frequency of strongrecruitment events increase because of a combination of higher overall abundance and reducedharvest rates induced by the 600 thousand tons cap imposed on all harvest control rules. Thedifferences in long term yield between U.S.A and Canada are also minimized when if fish startmigrating further north at younger ages, as was demonstrated on the early movement scenario.82Overall the 40:10 harvest control rule seems to be a good compromise for the two nations, andit performs relatively well across all performance metrics considered in this study. In the scenariowhere no strong recruitment events occurred, the 40:10 rule tended to reduce the biomass below40% of unfished biomass about 86% of the time. However, even when no strong recruitment eventswere considered, the biomass rarely went below 10% of unfished levels. Our results are similar tothose obtained by Taylor et al. (2014) and Hicks et al. (2016), who considered even lower TAC capsfor the aggregate stock; they also found that the value chosen for the TAC cap is inversely relatedto the conservation metrics, with lower caps resulting in more conservative outcomes. In alignmentwith other studies, the adoption of a formal TAC cap for the Pacific hake fishery is likely to resultin prevention of fishery closures, prevention of severe depletion of the stock (i.e., stock biomassfalling below 10-20% of average unfished levels), and the provision of a more stable distribution ofthe resource (given limited exploitation rates on strong cohorts).The present study has two major limitations. These are the limited sensitivity analyses for thekey biological parameters, and the simplified representation of measurement and observation errorin the analyses. I assume that the biological population parameters for the Pacific hake stock arestable over time (i.e., fixed and subject to random variability only) and equal to those reportedin the 2017 stock assessment (Berger et al., 2017b). Punt et al. (2008) shows that uncertaintyin stock recruitment parameters, more specifically, productivity and recruitment variability, canhave an impact on various performance metrics, including the percentage of time that biomassis above some threshold quantity. I observed similar results regarding changes in productivitywhen comparing the results across the three recruitment scenarios considered in this study. I alsoparameterized the movement dynamics of the resource based on literature documentation insteadof using a statistical model fit to data. In Chapter 2, I showed that the movement parameterscan be estimated if spatial catch at age composition data is available. An expansion of this studyshould include an evaluation of uncertainties in other key population dynamics parameters, suchas the growth rates, (and how they might vary for strong cohorts), as well as an estimation of the83movement parameters, which will affect the difference in performance between the two fishingnations.In this study, I used an auto-correlated error time series to represent the measurement error, i.e.,the error associated with stock assessments. This methodology was suggested by Walters (2004),and employed in other studies of evaluation of harvest control rules (e.g., Punt et al., 2008). How-ever, even though I believe that the results shown here are representative of the overall trends andtrade-offs between performance metrics experienced by each nation, the results might change ifother levels of autocorrelation and variance are used in the observation model. In a sensitivityanalysis (results not shown), I found that if lower observation error variances are considered, thelong term yields become higher and the difference in performance of long term yield betweenthe U.S.A. and Canada become more pronounced. (Moxnes, 2003) reviews evidence of this phe-nomenon when using a policy optimization model to investigate the impacts of measurement errorand stock uncertainty in policy outcomes. (Moxnes, 2003) found that if stock measurements aremore uncertain, the revenue from a fishery tends to decrease and the optimum policies becomemore conservative.In addition to the weaknesses described above, I would like to also emphasize that the perfor-mance metrics chosen for this study are somewhat arbitrary and chosen as a potential representationof the fishery objectives. The performance metrics used in this study are equal or similar to the onesused in the management strategy evaluation process carried out by the Pacific hake Joint TechnicalCommittee (Taylor et al., 2014; Hicks et al., 2016). However, if a similar exercise is to be usedto evaluate management options in the Pacific hake fishery, the preferred policies in practice willdepend on performance measures and objectives defined in the Pacific hake management process.Other studies have evaluated management options for migratory and transboundary fish stocks,usually accounting for the spatial effect indirectly. For example, Jones et al. (2016) implementedclosed-loop simulations for the Lake Erie walleye fisheries with the objective of comparing the per-formance of alternative harvest control rules, and report trade-offs between the Canadian commer-cial fisheries and the U.S. recreational fisheries. To account for changes in availability associated84with fish movement, they included time-varying catchability in their non-spatial model. In addi-tion, in recent years there has been increased interest in investigating the effects of spatial structureon stock assessment outcomes and management benchmarks (Berger et al., 2017a). Many simu-lation evaluation studies have been done to assess the impacts of spatial structure and movementdynamics on stock assessment performance and reference points (e.g., Lee et al., 2017; Goetheland Berger, 2017; Kerr et al., 2017; Carruthers et al., 2015). However the use of spatially explicitmodels in closed loop simulations remain relatively scarce, in part due to the computational burdenand the difficulties in fitting spatially structured models to data (Goethel et al., 2016).This study focused on Pacific hake as a case study and, for this reason, the model parameteriza-tion and structure set to mimic the dynamic of that resource. However, I believe that the frameworkdesigned here can be applied to many other transboundary resources that are subject to size seg-regation and/or migration range variability. The issue of management of transboundary stockssusceptible to changes in migration range has been explored in the literature under a game theoryapproach (e.g., Bailey et al., 2013; Hannesson, 2013; Liu et al., 2016). These studies generallypoint out that cooperative approaches to management perform better over the long term. However,Bjørndal and Ekerhovd (2014) point out that changes in migration range and spatial distributionare likely to impact international management agreements. Here, I present a framework that canaid cooperative management agreements to evaluate the impacts of harvest control rules and othermanagement procedures on the distribution of the stock. This tool can be used to detect changes inspatial distribution over time, and search for management procedures that minimize change in stockdistribution, helping to ensure equitable access to the resource by the parties sharing the stock.In summary, the framework and evaluation exercise presented in this study are valuable forthe management of the Pacific hake resource and for the management of migratory transbound-ary species in general. For the Pacific hake resource, I demonstrate the potential implications ofpopulation size truncation and change in migration range that can arise from using different harvestcontrol rules. I also demonstrate the effectiveness of the 40:10 harvest control rule and the potentialbenefits of imposing a TAC cap. For migratory transboundary species in general, I present a frame-85work to aid the identification of effective management procedures, that promote the sustainabilityof the stock and stable spatial distribution of the resource.86Chapter 5ConclusionThe overarching objective of this dissertation work was to explore questions that relate to the man-agement of migratory transboundary species. I focused on two main topics: (1) the interactionbetween age/size based migratory movement and the spatial availability of the resource and (2)time-varying fisheries selectivity associated with size segregation, migratory movement and cohorttargeting. I have developed two new modeling tools to address research questions related to thesetwo topics. In chapter 2, I present a continuous migration model capable of modeling cyclic mi-grations that are commonly found in iteroparous fish species. This model is suitable for exploring,testing and demonstrating potential issues of resource availability that arise from migratory move-ment and age/size segregation. In chapter 3, I present a length based stock assessment tool thatattempts to provide better estimates of time-varying fisheries selectivity. Finally, in chapter 4, I usethe model developed in chapter 2 in a closed-loop simulation framework to explore the impacts ofa large set of harvest control rules on management outcomes experienced by two nations sharinga transboundary resource. In the following section, I summarize the work in chapters 2 through 4,and discuss how these models and findings may contribute to the management of Pacific Hake andother migratory transboundary species around the globe.875.1 Research summaryIn chapter 2, I introduced a migration model that characterizes the cyclic migrations between feed-ing and spawning grounds that are common in iteroparous migratory species. A Lagrangian ap-proach to model movement is used to track individual cohorts (or sub groups within cohorts)through space and time. I demonstrate how the movement parameters in the model can be esti-mated from spatial catch at age data, a commonly available data type for many temperate exploitedstocks. The Lagrangian movement model I present is continuous in space and time. This is im-portant because in a continuous model it becomes unnecessary to delimit spatial areas from whichmovement rates are measured, and time steps can vary in size as needed. On the other hand, theLagrangian model requires an explicit migration hypothesis to generate the movement trajecto-ries. Such hypotheses, however, exist for many exploited species. For example, many migrationhypotheses for the ocean phase of Pacific salmon species are described by Groot and Margolis(1991). These generally include northward migrations following the North American coast thenreturning to re-enter their natal streams as fish mature. Tunas are another example of a groupthat in which some species perform trans oceanic cyclic migrations between spawning and feedinggrounds (Nikolic et al., 2017; Nakamura, 1969). Other species, like pelagic sharks are also believedto perform cyclic migrations, usually between inshore and offshore waters (Campana et al., 2011;Jorgensen et al., 2010). Similar behavior has also been demonstrated for flatfish like plaice (Hunteret al., 2003) and some cod populations (Robichaud and Rose, 2004).I illustrated the model performance using Pacific hake as a case study, which perform cyclicmigrations between the spawning grounds off southern California in the winter and the feedingareas along the North American coast all the way to northern British Columbia, Canada (Ressleret al., 2007). The migration extent is associated with fish age/size with fish moving further awayfrom spawning grounds as they age/grow. Because I used Pacific hake as an example, I made themigration range a function of age. However, the model can be extended to incorporate covariatesrepresenting biological and environmental forces that alter the distribution and migration rangeof exploited populations. I expect that this movement model will be a useful tool to model fish88migration and to illustrate how fisheries dynamics are affected by fish migration. The model couldalso be used as the basis for an operating model in closed loop simulation exercises to test therobustness of management frameworks that apply to populations and fisheries that are subject tospatial structure.In chapter 3, I introduce a new length-based stock assessment model: the length-SRA. Thismethod bypasses the requirement of estimating selectivity by calculating exploitation rate at lengthdirectly from observed catch at length data. The objective was to come up with a method thatwould be robust to time-varying selectivity, a phenomenon commonly associated with resourcessubject to spatial complexities and cohort targeting. I tested the performance of the Length-SRAwith a simulation-evaluation framework under three exploitation rate trajectories and under fixedand time-varying selectivity scenarios. The model produced parameter and derived managementquantities estimates with precision and accuracy that are comparable to that of other assessmentmodels, especially when considering time-varying selectivity (e.g., Martell and Stewart, 2014; Lin-ton and Bence, 2011). The selectivity estimates produced by the model were accurate over mostof the simulation scenarios, except when the exploitation rate time series showed no contrast, i.e.,exploitation rate was kept at values near the management target for most of the time series. Ingeneral the selectivity estimates were not very precise. The imprecision in selectivity estimates isprobably associated with the fact that the model assumes no error in the catch at length data, andtherefore incorporates all the observation error in the selectivity estimates.In addition, I explored the effects of misspecification of growth parameters on the length-SRAresults. The model was found to be extremely sensitive to the growth parameters input, particularlywhen it comes to selectivity estimates. I used the length-SRA model to assess two species: Pacifichake and Peruvian jack mackerel. Both species are believed to be subject to time-varying selectivityassociated with fisheries targeting of areas of high abundance and changes in population distributionover time. I found that both species presented time-varying and mainly dome shaped selectivity,which corroborated with the findings by Waterhouse et al. (2014) and Butterworth et al. (2014).However, as I pointed out earlier, the model is extremely sensitive to growth parameters, and it89remains uncertain whether the growth parameter estimates used for both Pacific hake and Peruvianjack mackerel are reliable. I recommend that whenever this model is used, extra caution is exercisedwhen choosing growth parameter estimates and that extensive sensitive analyses are performed onthe growth parameter values. Another alternative might be to use a Bayesian approach to integrateover the uncertainty in growth parameters by using distributions for the growth parameters values.This approach likely will not address the potential biases, but it will better capture the uncertaintycaused by misspecification of growth parameters.In chapter 4, I used the movement model described in chapter 2 in a closed loop simulationapproach to evaluate the performance of a large set of harvest control rules for the Pacific hakepopulation. I took advantage of the spatial capability of the movement model to explore the differ-ences in performance of each harvest control rule in relation to the outcomes experienced by theU.S.A. and Canada, the nations sharing the Pacific hake resource. I found that when the harvestcontrol rules allow for higher exploitation rates, issues of availability of the resource in Canadianwaters become more prominent, resulting in lower average yields for Canada when compared tothe lower exploitation rates. In order to assess the impacts of the occurrence of strong recruit-ment events in the population, I ran the evaluation under three distinct recruitment scenarios: nostrong recruitment, one strong recruitment event per decade, and two strong recruitment events perdecade. These scenarios were devised based on the historical recruitment of Pacific hake; strongrecruitment events are believed to have occurred in 1980, 1985, 1999, 2010 and 2014, and seemto have no apparent relationship with stock size (Berger et al., 2017b). I found that the availabilityissue became more acute, i.e., larger decreases of average yield for Canada when strong recruit-ment events are less frequent or absent. I also tested the sensitivity of the results to movementparameter assumptions, and found that the differences in relative average long term yield betweennations become less prominent if fish are assumed to migrate farther at younger ages. However,the impact of the movement parameter assumptions over the range of values explored tend to besmaller than that of the occurrence of strong recruitment events. This is because the occurrence of90strong recruitment events tends to inflate the biomass of the stock, eliminating problems related toavailability throughout the species range, even if movement rates are diminished.Among the harvest control rules being tested, I included the 40:10 harvest control rule with amaximum TAC cap. This rule is currently used for the Pacific hake management. The TAC cap isnot officially part of the Pacific hake treaty, but a maximum cap has effectively been implementedconsistently by the Pacific hake JMC over the past few years (Hicks et al., 2016). I found that thisrule performed well in terms of maintaining the resource exploitation at sustainable levels, and interms of mitigating the potential losses experienced by the Canadian fleet due to availability issues.In addition to demonstrating the performance of a set of harvest control rules for the man-agement of Pacific hake, the exercise presented in chapter 4 also has broader implications. Theapproach presented can be used to evaluate harvest control rules for other transboundary stocksthat are subject to changes in distribution and migration range.5.2 Future research directionsThe material presented in this dissertation has value for future research and for management ap-plications. The modeling tools are ready for management use provided that the necessary datais available. Future research could include further development of the modeling tools, applica-tion of the methods to other species, and further use of the tools in closed loop simulations andmanagement strategy evaluations. There are several ways in which the research presented in thisdissertation can be continued and enriched; I discuss some of those in the following paragraphs.Most of the methodology, particularly chapters 2 and 4, were produced with the Pacific hakestudy case in mind, and therefore, will need to be modified if other fisheries and resources areconsidered. The Lagrangian movement model can be modified to include alternative (non cyclic)movement functions, to model more than one spatial dimension simultaneously (e.g., centroids ofa distribution, longitude, latitude, and depth), and to incorporate other covariates, such as environ-mental variables. The explicit inclusion of environmental covariates in the Lagrangian movement91model could be used to test hypotheses derived from correlation between environmental variablesand fish distribution (e.g., Chen et al., 2005; Agostini et al., 2006; Mourato et al., 2014).In relation to the length-SRA model presented in chapter 3, an interesting expansion of themodel would be to produce a stochastic version of the stock reduction analysis model, similarto what was done by Walters et al. (2006). In that approach the SRA projections are producedbased on a range of hypothetical values for the main model parameters: unfished recruitment (R0)and recruitment compensation ratio (κ). Then, these projections are “filtered” either by removingthose in which the stock was driven into extinction or by using a sampling-importance resamplingroutine. This expansion would likely lead to a better characterization of uncertainties associatedwith the model estimates. In addition, further testing of the model could be done using a closed-loop simulation approach. Such testing would enhance the understanding of the model, not onlyin terms of accuracy and precision, but also in terms of achieving management goals, such assustainable and stable catches as well as conservation of the stock in the long term.The closed-loop simulation approach described in chapter 4 could be expanded in a variety ofways, for example, to test different stock assessment methods. Testing of alternative allocationsmethods, such as the fisheries footprint approach, suggested by Martell et al. (2015) could alsobe considered. The closed loop simulation framework included in this dissertation could also beincorporated into a full management strategy evaluation process, in which the questions explored,as well as the management objectives and performance metrics, would be the result of consultationwith managers, fisheries stakeholders, and scientists. In the case of Pacific hake, a managementstrategy evaluation process has been underway since 2013 (Taylor et al., 2014), and is likely tocontinue in the coming years with the expansion to a spatially structured operating model (Anony-mous, 2017). The work presented in this thesis, particularly the spatial migration model could be arelevant contribution to the future work carried by the Pacific hake Joint Technical Committee.In conclusion, the research presented here is comprised of two new modeling tools and a closedloop simulation framework. Each of these could be useful in research and management of trans-boundary fish stocks that are susceptible to age/size based migration. The use of these tools may92thus lead to improvement in understanding of the potential spatial consequences of managementactions and, thereby leading to improved resource management.93BibliographyAgostini, V. N., Francis, R. C., Hollowed, A. B., Pierce, S. D., Wilson, C., and Hendrix, A. N.(2006). The relationship between Pacific hake ( Merluccius productus ) distribution andpoleward subsurface flow in the California Current System. Canadian Journal of Fisheries andAquatic Sciences, 63(12):2648–2659.Alheit, J. and Pitcher, T. J. (1995). Hake: Biology, fisheries and markets. Springer Science &Business Media. Google-Books-ID: NLXvCAAAQBAJ.Anonymous (2013). 2013 jack mackerel stock assessment. Technical Report SC-01 Annex 05,South Pacific Regional Fisheries Management Organisation.Anonymous (2017). JOINT REPORT OF THE CANADIAN AND U.S. ADVISORY PANEL ofThe Agreement between the Government of the United States of America and the Governmentof Canada on Pacific Hake/Whiting (the Agreement). Technical report.Bailey, K. M., Francis, R. C., and Stevens, P. R. (1982). The life history and fishery of Pacificwhiting, Merluccius productus. Northwest and Alaska Fisheries Center, National MarineFisheries Service, US Department of commerce.Bailey, M., Rashid Sumaila, U., and Lindroos, M. (2010). Application of game theory to fisheriesover three decades. Fisheries Research, 102(1–2):1–8.Bailey, M., Sumaila, U. R., and Martell, S. (2013). Can Cooperative Management of TunaFisheries in the Western Pacific Solve the Growth Overfishing Problem? Strategic Behaviorand the Environment, 3(1–2):31–66.Barbaro, A., Einarsson, B., Birnir, B., Sigurðsson, S., Valdimarsson, H., Pálsson, Ó. K.,Sveinbjörnsson, S., and Sigurðsson, Þ. (2009). Modelling and simulations of the migration ofpelagic fish. ICES Journal of Marine Science: Journal du Conseil, page fsp067.Berger, A. M., Goethel, D. R., Lynch, P. D., Quinn, T., Mormede, S., McKenzie, J., and Dunn, A.(2017a). Space oddity: The mission for spatial integration. Canadian Journal of Fisheries andAquatic Sciences, 74(11):1698–1716.Berger, A. M., Grandin, C., Taylor, I. G., Edwards, A. M., and Cox, S. (2017b). Status of thePacific Hake (whiting) stock in U.S. and Canadian waters in 2017. Technical report.Berger, A. M., Jones, M. L., Zhao, Y., and Bence, J. R. (2012). Accounting for spatial populationstructure at scales relevant to life history improves stock assessment: The case for Lake Eriewalleye Sander vitreus. Fisheries Research, 115:44–59.94Beverton, R. J. H. and Holt, S. J. (1957). On the Dynamics of Exploited Fish Populations,volume 19 of Investment series. U.K. Ministry of Agriculture and Fisheries, London.Google-Books-ID: BqbnCAAAQBAJ.Bez, N., Oliveira, E. D., and Duhamel, G. (2006). Repetitive fishing, local depletion, and fishingefficiencies in the Kerguelen Islands fisheries. ICES Journal of Marine Science: Journal duConseil, 63(3):532–542.Bjørndal, T. and Ekerhovd, N.-A. (2014). Management of Pelagic Fisheries in the North EastAtlantic: Norwegian Spring Spawning Herring, Mackerel, and Blue Whiting. Marine ResourceEconomics; Chicago, 29(1):69–83.Branch, T. A. and Hilborn, R. (2010). A general model for reconstructing salmon runs. CanadianJournal of Fisheries and Aquatic Sciences, 67(5):886–904.Brodie, W. B., Walsh, S. J., and Atkinson, D. B. (1998). The effect of stock abundance on rangecontraction of yellowtail flounder (Pleuronectes ferruginea) on the Grand Bank ofNewfoundland in the Northwest Atlantic from 1975 to 1995. Journal of Sea Research,39(1–2):139–152.Brooks, E. N., Deroba, J. J., and Wilberg, M. (2015). When “data” are not data: the pitfalls of posthoc analyses that use stock assessment model output. Canadian Journal of Fisheries andAquatic Sciences, 72(4):634–641.Bunch, A. J., Walters, C. J., and Coggins, L. G. (2013). Measurement Error in Fish Lengths:Evaluation and Management Implications. Fisheries, 38(7):320–326.Butterworth, D. S., Rademeyer, R. A., Brandão, A., Geromont, H. F., and Johnston, S. J. (2014).Does selectivity matter? A fisheries management perspective. Fisheries Research,158:194–204.Caddy, J. F. (1975). Spatial model for an exploited shellfish population, and its application to theGeorges Bank scallop fishery. Journal of the Fisheries Board of Canada, 32(8):1305–1328.Camhi, M. D., Pikitch, E. K., and Babcock, E. A. (2009). Sharks of the Open Ocean: Biology,Fisheries and Conservation. John Wiley & Sons. Google-Books-ID: lc9MyMaXHgEC.Campana, S. E., Dorey, A., Fowler, M., Joyce, W., Wang, Z., Wright, D., and Yashayaev, I. (2011).Migration Pathways, Behavioural Thermoregulation and Overwintering Grounds of BlueSharks in the Northwest Atlantic. PLoS ONE, 6(2):e16854.Canadian Advisory Panel (2013). Canadian Advisory Panel TAC Position.Carlson, A. E., Hoffmayer, E. R., Tribuzio, C. A., and Sulikowski, J. A. (2014). The Use ofSatellite Tags to Redefine Movement Patterns of Spiny Dogfish (Squalus acanthias) along theU.S. East Coast: Implications for Fisheries Management. PLoS ONE, 9(7):e103384.Carruthers, T. R., McAllister, M. K., and Taylor, N. G. (2011). Spatial surplus productionmodeling of Atlantic tunas and billfish. Ecological Applications, 21(7):2734–2755.95Carruthers, T. R., Walter, J. F., McAllister, M. K., and Bryan, M. D. (2015). Modellingage-dependent movement: an application to red and gag groupers in the Gulf of Mexico.Canadian Journal of Fisheries and Aquatic Sciences, 72(8):1159–1176.WOS:000358608600004.Cave, J. D. and Gazey, W. J. (1994). A Preseason Simulation Model for Fisheries on Fraser RiverSockeye Salmon (Oncorhynchus nerka). Canadian Journal of Fisheries and Aquatic Sciences,51(7):1535–1549.Chen, I.-C., Lee, P.-F., and Tzeng, W.-N. (2005). Distribution of albacore (Thunnus alalunga) inthe Indian Ocean and its relation to environmental factors. Fisheries Oceanography,14(1):71–80.Chittaro, P. M., Zabel, R. W., Palsson, W., and Grandin, C. (2013). Population interconnectivityand implications for recovery of a species of concern, the Pacific hake of Georgia Basin.Marine Biology, 160(5):1157–1170.Ciannelli, L., Fauchald, P., Chan, K. S., Agostini, V. N., and Dingsør, G. E. (2008). Spatialfisheries ecology: Recent progress and future prospects. Journal of Marine Systems,71(3–4):223–236.Coggins, L. G. and Quinn, T. J. (1998). A simulation study of the effects of aging error andsample size on sustained yield estimates. Fishery stock assessment models, pages 955–975.Costa, D. P., Breed, G. A., and Robinson, P. W. (2012). New Insights into Pelagic Migrations:Implications for Ecology and Conservation. Annual Review of Ecology, Evolution, andSystematics, 43(1):73–96.Cox, S. P., Kronlund, A. R., and Benson, A. J. (2013). The roles of biological reference points andoperational control points in management procedures for the sablefish (<spanclass="italic">Anoplopoma fimbria</span>) fishery in British Columbia, Canada.Environmental Conservation, 40(4):318–328.Criddle, K. R. and Strong, J. W. (2014). Straddling the line: cooperative and non-cooperativestrategies for management of Bering Sea pollock. Fisheries Science, 80(2):193–203.Deroba, J. J. and Bence, J. R. (2008). A review of harvest policies: Understanding relativeperformance of control rules. Fisheries Research, 94(3):210–223.Food and Agriculture Organization of the United Nations (1996). Precautionary approach tocapture fisheries and species introductions. Technical report, FAO Rome (Italy).Forrester, C. R., Beardsley, A. J., and Takahashi, Y. (1978). Groundfish, Shrimp, and HerringFisheries in the Bering Sea and Northeast Pacific - Historical catch statistics through 1970.Fournier, D. A., Hampton, J., and Sibert, J. R. (1998). MULTIFAN-CL: a length-based,age-structured model for fisheries stock assessment, with application to South Pacific albacore,Thunnus alalunga. Canadian Journal of Fisheries and Aquatic Sciences, 55(9):2105–2116.96Fournier, D. A., Skaug, H. J., Ancheta, J., Ianelli, J., Magnusson, A., Maunder, M. N., Nielsen, A.,and Sibert, J. (2012). AD Model Builder: using automatic differentiation for statisticalinference of highly parameterized complex nonlinear models. Optimization Methods andSoftware, 27(2):233–249.Francis, R. I. C. C. and Clark, M. R. (1998). Inferring spawning migrations of orange roughy(Hoplostethus atlanticus) from spawning ogives. Marine and Freshwater Research,49(2):103–108.Gerlotto, F., Gutiérrez, M., and Bertrand, A. (2012). Insight on population structure of the Chileanjack mackerel ( Trachurus murphyi ). Aquatic Living Resources, 25(4):341–355.Giske, J., Huse, G., and Berntsen, J. (2001). Spatial modelling for marine resource management,with a focus on fish. Sarsia, 86(6):405–410.Goethel, D. R. and Berger, A. M. (2017). Accounting for spatial complexities in the calculation ofbiological reference points: effects of misdiagnosing population structure for stock statusindicators. Canadian Journal of Fisheries and Aquatic Sciences, 74(11):1878–1894.Goethel, D. R., Kerr, L. A., and Cadrin, S. X. (2016). Incorporating spatial population structureinto the assessment-management interface of marine resources. In Management Science inFisheries: An Introduction to Simulation-based Methods, page 319.Goethel, D. R., Quinn, T. J., and Cadrin, S. X. (2011). Incorporating Spatial Structure in StockAssessment: Movement Modeling in Marine Fish Population Dynamics. Reviews in FisheriesScience, 19(2):119–136.Government of Canada, F. a. O. C. (2009). A fishery decision-making framework incorporatingthe precautionary approach.Groot, C. and Margolis, L. (1991). Pacific salmon life histories. UBC Press.Grüss, A., Kaplan, D. M., Guénette, S., Roberts, C. M., and Botsford, L. W. (2011).Consequences of adult and juvenile movement for marine protected areas. BiologicalConservation, 144(2):692–702.Gudmundsson, G. and Gunnlaugsson, T. (2012). Selection and estimation of sequentialcatch-at-age models. Canadian Journal of Fisheries and Aquatic Sciences, 69(11):1760–1772.Hall, D. L., Hilborn, R., Stocker, M., and Walters, C. J. (1988). Alternative harvest strategies forPacific herring (Clupea harengus pallasi). Canadian Journal of Fisheries and Aquatic Sciences,45(5):888–897.Hammerschlag, N., Gallagher, A. J., and Lazarre, D. M. (2011). A review of shark satellitetagging studies. Journal of Experimental Marine Biology and Ecology, 398(1):1–8.Hannesson, R. (2013). Sharing the Northeast Atlantic mackerel. ICES Journal of Marine Science,70(2):259–269.Harden Jones, F. R. (1968). Fish migration. St. Martin’s Press, New York.97Hawkshaw, M. and Walters, C. (2015). Harvest control rules for mixed-stock fisheries coping withautocorrelated recruitment variation, conservation of weak stocks, and economic well-being.Canadian Journal of Fisheries and Aquatic Sciences, 72(5):759–766.Henríquez, V., Licandeo, R., Cubillos, L. A., and Cox, S. P. (2016). Interactions between ageingerror and selectivity in statistical catch-at-age models: simulations and implications forassessment of the Chilean Patagonian toothfish fishery. ICES Journal of Marine Science:Journal du Conseil, 73(4):1074–1090.Hicks, A. C., Cox, S. P., Taylor, N., Taylor, I. G., Grandin, C., and Ianelli, J. N. (2016).Conservation and yield performance of harvest control rules for the transboundary Pacific hakefishery in US and Canadian waters. Management Science in Fisheries: An Introduction toSimulation-based Methods, page 69.Hilborn, R. (1986). A comparison of alternative harvest tactics for invertebrate and managementof invertebrates. In North Pacific Workshop on Stock Assessment and Management ofInvertebrates: Nanaimo, British Columbia, May 7-10, 1984, number 92 in Canadian specialpublication of fisheries and aquatic sciences, pages 313–317.Hilborn, R. and Walters, C. J. (1992). Quantitative Fisheries Stock Assessment: Choice, Dynamicsand Uncertainty/Book and Disk. Springer Science & Business Media.Hunter, E., Metcalfe, J. D., and Reynolds, J. D. (2003). Migration route and spawning area fidelityby North Sea plaice. Proceedings of the Royal Society B: Biological Sciences,270(1529):2097–2103.Hurtado-Ferro, F., Punt, A. E., and Hill, K. T. (2014). Use of multiple selectivity patterns as aproxy for spatial structure. Fisheries Research, 158:102–115.Ishimura, G., Herrick, S., and Sumaila, U. R. (2013). Fishing games under climate variability:transboundary management of Pacific sardine in the California Current System. EnvironmentalEconomics and Policy Studies, 15(2):189–209.Ishimura, G., Punt, A. E., and Huppert, D. D. (2005). Management of fluctuating fish stocks: thecase of Pacific whiting. Fisheries Research, 73(1-2):201–216.Jones, M. L., Catalano, M. J., Peterson, L. K., and Berger, A. M. (2016). Stakeholder-centereddevelopment of a harvest control rule for Lake Erie walleye. In Management Science inFisheries: An Introduction to Simulation-based Methods, page 163.Jorgensen, S. J., Reeb, C. A., Chapple, T. K., Anderson, S., Perle, C., Van Sommeran, S. R.,Fritz-Cope, C., Brown, A. C., Klimley, A. P., and Block, B. A. (2010). Philopatry and migrationof Pacific white sharks. Proceedings of the Royal Society B: Biological Sciences,277(1682):679–688.Kerr, L. A., Cadrin, S. X., Secor, D. H., and Taylor, N. G. (2017). Modeling the implications ofstock mixing and life history uncertainty of Atlantic bluefin tuna. Canadian Journal ofFisheries and Aquatic Sciences, 74(11):1990–2004.98Kerr, L. A. and Goethel, D. R. (2014). Simulation Modeling as a Tool for Synthesis of StockIdentification Information. In Stock Identification Methods, pages 501–533. Elsevier.Kimura, D. K., Balsiger, J. W., and Ito, D. H. (1984). Generalized Stock Reduction Analysis.Canadian Journal of Fisheries and Aquatic Sciences, 41(9):1325–1333.Kimura, D. K. and Chikuni, S. (1987). Mixtures of Empirical Distributions: An IterativeApplication of the Age- Length Key. Biometrics, 43(1):23–35.Kimura, D. K. and Tagart, J. V. (1982). Stock Reduction Analysis, Another Solution to the CatchEquations. Canadian Journal of Fisheries and Aquatic Sciences, 39(11):1467–1472.King, J. R., McFarlane, G. A., Jones, S. R., Gilmore, S. R., and Abbott, C. L. (2012). Stockdelineation of migratory and resident Pacific hake in Canadian waters. Fisheries Research,114:19–30.Kvamsdal, S. F., Eide, A., Ekerhovd, N.-A., Enberg, K., Gudmundsdottir, A., Hoel, A. H., Mills,K. E., Mueter, F. J., Ravn-Jonsen, L., Sandal, L. K., Stiansen, J. E., and Vestergaard, N. (2016).Harvest control rules in modern fisheries management. Elementa: Science of the Anthropocene,4:000114.Lee, H.-H., Piner, K. R., Maunder, M. N., Taylor, I. G., and Methot, R. D. (2017). Evaluation ofalternative modelling approaches to account for spatial effects due to age-based movement.Canadian Journal of Fisheries and Aquatic Sciences, 74(11):1832–1844.Lett, C., Verley, P., Mullon, C., Parada, C., Brochier, T., Penven, P., and Blanke, B. (2008). ALagrangian tool for modelling ichthyoplankton dynamics. Environmental Modelling &Software, 23(9):1210–1214.Linton, B. C. and Bence, J. R. (2011). Catch-at-age assessment in the face of time-varyingselectivity. ICES Journal of Marine Science, 68(3):611–625.Liu, X. and Heino, M. (2013). Comparing Proactive and Reactive Management: Managing aTransboundary Fish Stock Under Changing Environment. Natural Resource Modeling,26(4):480–504.Liu, X., Lindroos, M., and Sandal, L. (2016). Sharing a Fish Stock When Distribution and HarvestCosts are Density Dependent. Environmental and Resource Economics; Dordrecht,63(3):665–686.Lloris, D., Matallanas, J., Oliver Reus, P. A., and Food and Agriculture Organization of the UnitedNations (2005). Hakes of the world (family Merlucciidae): an annotated and illustratedcatalogue of hake species known to date. Food and Agriculture Organization of the UnitedNations, Rome.Lo, N. C., Macewicz, B. J., and Griffith, D. A. (2011). Migration of Pacific Sardine (SardinopsSagax) Off the West Coast of United States in 2003–2005. Bulletin of Marine Science,87(3):395–412.99Magnuson-Stevens act (2007). Magnuson-Stevens Fishery Conservation and Management Act.Martell, S. and Stewart, I. (2014). Towards defining good practices for modeling time-varyingselectivity. Fisheries Research, 158:84–95.Martell, S., Stewart, I., and Wor, C. (2015). Exploring index-based PSC limits for Pacific halibut.In Report of Assessment and Research Activities 2015.Martell, S. J. D., Walters, C. J., and Wallace, S. S. (2000). The use of marine protected areas forconservation of lingcod (Ophiodon elongatus). Bulletin of Marine Science, 66(3):729–743.Maunder, M. N., Crone, P. R., Valero, J. L., and Semmens, B. X. (2014). Selectivity: Theory,estimation, and application in fishery stock assessment models. Fisheries Research, 158:1–4.Maunder, M. N. and Deriso, R. B. (2003). Estimation of recruitment in catch-at-age models.Canadian Journal of Fisheries and Aquatic Sciences, 60(10):1204–1216.Maury, O. and Gascuel, D. (2001). ‘Local overfishing’and fishing tactics: theoreticalconsiderations and applied consequences in stock assessment studied with a numericalsimulator of fisheries. Aquatic Living Resources, 14(4):203–210.McAllister, M. K. (1998). Modeling the effects of fish migration on bias and variance inarea-swept estimates of biomass: a vector-based approach. Canadian Journal of Fisheries andAquatic Sciences, 55(12):2622–2641.McClure, M., Workman, G., Prager, M. H., Holt, K., Sampson, D., and Branch, T. A. (2014). JointU.S.-Canada Scientific Review Group Report.McKeown, B. A. (1984). Fish migration. Timber press.Melnychuk, M. C., Peterson, E., Elliott, M., and Hilborn, R. (2017). Fisheries managementimpacts on target species status. Proceedings of the National Academy of Sciences,114(1):178–183.Merten, W., Appeldoorn, R., and Hammond, D. (2016). Movement dynamics of dolphinfish(Coryphaena hippurus) in the northeastern Caribbean Sea: Evidence of seasonal re-entry intodomestic and international fisheries throughout the western central Atlantic. FisheriesResearch, 175:24–34.Mesnil, B. and Shepherd, J. G. (1990). A hybrid age- and length-structured model for assessingregulatory measures in multiple-species, multiple-fleet fisheries. ICES Journal of MarineScience, 47(2):115–132.Methot, R. D. and Dorn, M. W. (1995). Biology and fisheries of North Pacific Hake (M.productus). In Hake: biology, fisheries and markets, pages 389–414. Springer, [S.l.].Methot, R. D. and Wetzel, C. R. (2013). Stock synthesis: A biological and statistical frameworkfor fish stock assessment and fishery management. Fisheries Research, 142:86–99.100Miller, K. A. and Munro, G. R. (2004). Climate and Cooperation: A New Perspective on theManagement of Shared Fish Stocks. Marine Resource Economics, 19(3):367–393.Minte-Vera, C. V., Maunder, M. N., Aires-da Silva, A. M., Satoh, K., and Uosaki, K. (2017). Getthe biology right, or use size-composition data at your own risk. Fisheries Research,192:114–125.Misund, O. A., Vilhjálmsson, H., Jákupsstovu, S. H. í., Røttingen, I., Belikov, S., Asthorsson, O.,Blindheim, J., Jónsson, J., Krysov, A., Malmberg, S. A., and Sveinbjórnsson, S. (1998).Distribution, migration and abundance of norwegian spring spawning herring in relation to thetemperature and zooplankton biomass in the Norwegian sea as recorded by coordinated surveysin spring and summer 1996. Sarsia, 83(2):117–127.Morais, P. and Daverat (2016). An Introduction to Fish Migration. CRC Press, 1st edition edition.Mourato, B. L., Hazin, F., Bigelow, K., Musyl, M., Carvalho, F., and Hazin, H. (2014).Spatio-temporal trends of sailfish, Istiophorus platypterus catch rates in relation to spawningground and environmental factors in the equatorial and southwestern Atlantic Ocean. FisheriesOceanography, 23(1):32–44.Moxnes, E. (2003). Uncertain measurements of renewable resources: approximations, harvestingpolicies and value of accuracy. Journal of Environmental Economics and Management,45(1):85–108.Mucientes, G. R., Queiroz, N., Sousa, L. L., Tarroso, P., and Sims, D. W. (2009). Sexualsegregation of pelagic sharks and the potential threat from fisheries. Biology Letters,5(2):156–159.Munro, G. (2005). The conservation and management of shared fish stocks: legal and economicaspects. Number 465 in FAO Fisheries and Aquaculture Technical Papers. FAO.Munro, G. R. (1979). The Optimal Management of Transboundary Renewable Resources. TheCanadian Journal of Economics / Revue canadienne d’Economique, 12(3):355–376.Nakamura, H. (1969). Tuna distribution and migration. Fishing news books.Nelson, S. (2014). Overview: Freezer Trawlers in the BC Groundfish Trawl Fishery. Prepared for:Province of British Columbia, Corporate Governance, Policy and Legislation Branch.Nichol, D. G. and Chilton, E. A. (2006). Recuperation and behaviour of Pacific cod afterbarotrauma. ICES Journal of Marine Science: Journal du Conseil, 63(1):83–94.Nielsen, A. and Berg, C. W. (2014). Estimation of time-varying selectivity in stock assessmentsusing state-space models. Fisheries Research, 158:96–101.Nikolic, N., Morandeau, G., Hoarau, L., West, W., Arrizabalaga, H., Hoyle, S., Nicol, S. J.,Bourjea, J., Puech, A., Farley, J. H., Williams, A. J., and Fonteneau, A. (2017). Review ofalbacore tuna, Thunnus alalunga, biology, fisheries and management. Reviews in Fish Biologyand Fisheries, 27(4):775–810.101Okamura, H., McAllister, M. K., Ichinokawa, M., Yamanaka, L., and Holt, K. (2014). Evaluationof the sensitivity of biological reference points to the spatio-temporal distribution of fishingeffort when seasonal migrations are sex-specific. Fisheries Research, 158:116–123.Pine, W. E., Pollock, K. H., Hightower, J. E., Kwak, T. J., and Rice, J. A. (2003). A Review ofTagging Methods for Estimating Fish Population Size and Components of Mortality. Fisheries,28(10):10–23.Pope, J. (1972). An Investigation of the Accuracy of Virtual Population Analysis Using CohortAnalysis. Research Bulletin International Commission for the Northwest Atlantic Fisheries,9:65–74.Punt, A. E. (2010). Harvest Control Rules and Fisheries Management. In Handbook of FisheriesConservation and Management. Oxford University Press.Punt, A. E., Butterworth, D. S., de Moor, C. L., De Oliveira, J. A. A., and Haddon, M. (2016).Management strategy evaluation: best practices. Fish and Fisheries, pages n/a–n/a.Punt, A. E. and Donovan, G. P. (2007). Developing management procedures that are robust touncertainty: lessons from the International Whaling Commission. ICES Journal of MarineScience, 64(4):603–612.Punt, A. E., Dorn, M. W., and Haltuch, M. A. (2008). Evaluation of threshold managementstrategies for groundfish off the U.S. West Coast. Fisheries Research, 94(3):251–266.Punt, A. E., Haddon, M., and Tuck, G. N. (2015). Which assessment configurations perform bestin the face of spatial heterogeneity in fishing mortality, growth and recruitment? A case studybased on pink ling in Australia. Fisheries Research, 168:85–99.Quinn II, T. J., Fagen, R., and Zheng, J. (1990). Threshold Management Policies for ExploitedPopulations. Canadian Journal of Fisheries and Aquatic Sciences, 47(10):2016–2029.Reed, W. J. (1979). Optimal escapement levels in stochastic and deterministic harvesting models.Journal of environmental economics and management, 6(4):350–363.Ressler, P. H., Holmes, J. A., Fleischer, G. W., Thomas, R. E., and Cooke, K. C. (2007). Pacifichake, Merluccius productus, autecology: a timely review. U S National Marine FisheriesService Marine Fisheries Review, 69(1-4). ZOOREC:ZOOR14601004724.Robichaud, D. and Rose, G. A. (2004). Migratory behaviour and range in Atlantic cod: inferencefrom a century of tagging. Fish and Fisheries, 5(3):185–214.Rodríguez-Sánchez, R., Lluch-Belda, D., Villalobos, H., and Ortega-García, S. (2002). Dynamicgeography of small pelagic fish populations in the California Current System on the regime timescale (1931?1997). Canadian Journal of Fisheries and Aquatic Sciences, 59(12):1980–1988.Ruttan, L. M. (2003). Finding fish: grouping and catch-per-unit-effort in the Pacific hake (Merluccius productus ) fishery. Canadian Journal of Fisheries and Aquatic Sciences,60(9):1068–1077.102Sainsbury, K. J., Punt, A. E., and Smith, A. D. M. (2000). Design of operational managementstrategies for achieving fishery ecosystem objectives. ICES Journal of Marine Science: Journaldu Conseil, 57(3):731–741.Sampson, D. B. (2014). Fishery selection and its relevance to stock assessment and fisherymanagement. Fisheries Research, 158:5–14.Sampson, D. B. and Scott, R. D. (2012). An exploration of the shapes and stability ofpopulation–selection curves. Fish and Fisheries, 13(1):89–104.Schnute, J. T. and Richards, L. J. (1995). The influence of error on population estimates fromcatch-age models. Canadian Journal of Fisheries and Aquatic Sciences, 52(10):2063–2077.Secor, D. H. (2015). Migration Ecology of Marine Fishes. The Johns Hopkins University Press.Sibert, J. R., Hampton, J., Fournier, D. A., and Bills, P. J. (1999). An advection-diffusion-reactionmodel for the estimation of fish movement parameters from tagging data, with application toskipjack tuna (Katsuwonus pelamis). Canadian Journal of Fisheries and Aquatic Sciences,56(6):925–938.Sippel, T., Paige Eveson, J., Galuardi, B., Lam, C., Hoyle, S., Maunder, M., Kleiber, P., Carvalho,F., Tsontos, V., Teo, S. L. H., Aires-da Silva, A., and Nicol, S. (2015). Using movement datafrom electronic tags in fisheries stock assessment: A review of models, technology andexperimental design. Fisheries Research, 163:152–160.Spijkers, J. and Boonstra, W. J. (2017). Environmental change and social conflict: the northeastAtlantic mackerel dispute. Regional Environmental Change, 17(6):1835–1851.Sullivan, P. J., Lai, H.-L., and Gallucci, V. F. (1990). A Catch-at-Length Analysis thatIncorporates a Stochastic Model of Growth. Canadian Journal of Fisheries and AquaticSciences, 47(1):184–198.Sumaila, U. R. (1999). A review of game-theoretic models of fishing. Marine Policy, 23(1):1–10.Sumaila, U. R. (2013). Game theory and fisheries: essays on the tragedy of free for all fishing.Number 41 in Routledge explorations in environmental economics. Routledge, Abingdon,Oxon ; New York.Taylor, I. G., Grandin, C., Hicks, A. C., Taylor, N., and Cox, S. (2015). Status of the Pacific Hake(whiting) stock in U.S. and Canadian waters in 2015.Taylor, N. G., Hicks, A. C., Taylor, I. G., Grandin, C., and Cox, S. (2014). Status of the PacificHake (whiting) stock in U.S. and Canadian waters in 2014 with a management strategyevaluation.Thompson, G. G. (1994). Confounding of gear selectivity and the natural mortality rate in caseswhere the former is a nonmonotone function of age. Canadian Journal of Fisheries and AquaticSciences, 51(12):2654–2664.103Thorson, J. T. and Cope, J. M. (2015). Catch curve stock-reduction analysis: An alternativesolution to the catch equations. Fisheries Research, 171(Supplement C):33–41.Tong, Y., Chen, X., and Kolody, D. (2014). Evaluation of three harvest control rules for BigeyeTuna (Thunnus obesus) fisheries in the Indian Ocean. Journal of Ocean University of China,13(5):811–819.United States State Department (2004). Agreement with Canada on Pacific hake/whiting: messagefrom the President of the United States transmitting agreement between the government of theUnited States of America and the government of Canada on Pacific hake/whiting, done atSeattle, November 21, 2003, volume 108. USGPO.Vasilakopoulos, P., O’Neill, F. G., and Marshall, C. T. (2016). The unfulfilled potential of fisheriesselectivity to promote sustainability. Fish and Fisheries, 17(2):399–416.Walters, C. (1998). Evaluation of quota management policies for developing fisheries. CanadianJournal of Fisheries and Aquatic Sciences, 55(12):2691–2705.Walters, C. (2004). Simple representation of the dynamics of biomass error propagation for stockassessment models. Canadian Journal of Fisheries and Aquatic Sciences, 61(7):1061–1065.Walters, C. and Ludwig, D. (1994). Calculation of Bayes posterior probability distributions forkey population parameters. Canadian Journal of Fisheries and Aquatic Sciences,51(3):713–722.Walters, C. and Parma, A. M. (1996). Fixed exploitation rate strategies for coping with effects ofclimate change. Canadian Journal of Fisheries and Aquatic Sciences, 53(1):148–158.Walters, C., Pauly, D., and Christensen, V. (1999). Ecospace: Prediction of mesoscale spatialpatterns in trophic relationships of exploited ecosystems, with emphasis on the impacts ofmarine protected areas. Ecosystems, 2(6):539–554. WOS:000084534100007.Walters, C. and Punt, A. (1994). Placing Odds on Sustainable Catch Using Virtual PopulationAnalysis and Survey Data. Canadian Journal of Fisheries and Aquatic Sciences,51(4):946–958.Walters, C. J., Martell, S. J., and Korman, J. (2006). A stochastic approach to stock reductionanalysis. Canadian Journal of Fisheries and Aquatic Sciences, 63(1):212–223.Walters, C. J. and Martell, S. J. D. (2004). Fisheries Ecology and Management. PrincetonUniversity Press.Waterhouse, L., Sampson, D. B., Maunder, M., and Semmens, B. X. (2014). Using areas-as-fleetsselectivity to model spatial fishing: Asymptotic curves are unlikely under equilibriumconditions. Fisheries Research, 158:15–25.Webster, R. A., Clark, W. G., Leaman, B. M., and Forsberg, J. E. (2013). Pacific halibut on themove: a renewed understanding of adult migration from a coastwide tagging study. CanadianJournal of Fisheries and Aquatic Sciences, 70(4):642–653.104White, J. W., Nickols, K. J., Malone, D., Carr, M. H., Starr, R. M., Cordoleani, F., Baskett, M. L.,Hastings, A., and Botsford, L. W. (2016). Fitting state-space integral projection models tosize-structured time series data to estimate unknown parameters. Ecological Applications,26(8):2677–2694.Winter, A., Foy, R. J., and Trussell, M. (2007). Factors influencing the mortality of tagged walleyepollock captured using a trawl net. North Pacific Research Board Final Report, 506(23):2.105


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items