The ecological module of BOATS-1 . 0 : a bioenergetically-constrained model of marine upper trophic levels suitable for studies of fisheries and ocean biogeochemistry

Introduction Conclusions References


Introduction
Humans have harvested fish and marine resources since prehistoric times, but due to the development of modern fish capture technologies since the end of the Second World War, and to a strong increase in demand arising from increasing population, global wild harvest increased at an unprecedented rate following 1945.This strong appetite for marine resources has had important impacts on marine ecosystems.A significant fraction of fisheries are overexploited, and estimates of the fraction of collapses range from 7-13 to 25 % of all fisheries (Mullon et al., 2005;Branch et al., 2011).Large finfish biomass is thought to be significantly depleted relative to its pre-harvest state (Myers and Worm, 2003), numerous species of finfish and invertebrates have witnessed range re-D. A. Carozza et al.: BOATS-1.0 ductions (local extinctions) (McCauley et al., 2015), and an index of marine finfish biomass indicates an aggregate loss of 38 % over many species (Hutchings et al., 2010).Despite increasing harvesting effort (Watson et al., 2013b), annual wild harvest appears to have peaked globally in the early 1990s (Watson et al., 2004;Pauly, 2007;FAO, 2014) at an annual rate that has been recently estimated at 130 million tonnes (Mt) per year (Pauly and Zeller, 2016), since which time it appears to have declined.As older coastal fisheries have become increasingly depleted (Jackson, 2001), harvest has extended to more taxa as well as further from the coast and deeper in the water column (Norse et al., 2012;Watson and Morato, 2013).
Anthropogenic climate change, on the other hand, is already altering nutrient dynamics and primary production through its effects on ocean temperature and circulation (Doney et al., 2012), with demonstrated consequences on the distributions of several fish populations (Pinsky et al., 2013;Walsh et al., 2015).Global climate models suggest that increased surface water stratification due to warming could decrease nutrient upwelling and so reduce net primary production (Sarmiento et al., 2004;Steinacher et al., 2010;Bopp et al., 2013).Warming can also directly influence fish biomass by affecting physiological rates that influence growth, mortality, reproduction, recruitment, and migration (Brander, 2010;Sumaila et al., 2011).Despite progress in identifying important mechanisms of biomass change, important uncertainties remain in constraining the overall impact and the spatial distribution of change in net primary production (Taucher and Oschlies, 2011) and fish biomass, with current analyses pointing toward heterogeneous spatial change in fish production and harvest potential (Cheung et al., 2010;Barange et al., 2014;Lefort et al., 2014).
Research in fisheries and fisheries economics often focusses on particular species, regions, and markets.In recent years, generalized, spatially resolved models of the marine ecosystem applicable to the global domain have been developed, but most are not directly coupled with predictive models of fishing activity (Jennings et al., 2008;Lefort et al., 2014;Watson et al., 2014).Our intention is to model fisheries and economic harvesting as parts of an integrated system that is bioenergetically constrained, and based on fundamental physical, ecological, and economic principles.The ecological module of the BiOeconomic mArine Trophic Sizespectrum model (BOATS) aims to represent commercial organisms as a set of super-organism populations (that we refer to as groups) that grow, reproduce, and die, taking into account their dependence on local environmental variables in the framework of a two-dimensional grid of the global ocean.The approach is structurally simpler than that of Christensen et al. (2015), and bears similarity with that of Jennings and Collingridge (2015), but unlike these models the BOATS ecological model explicitly treats life history and reproduction, similar to Maury et al. (2007).
The true ecological structure of marine communities is very complex, and includes many species-level ecological dynamics that are not understood at a mechanistic predictive level.A typical oceanic food web consists of dozens or more of interacting species, whose sizes span several orders of magnitude and whose lifetimes range from days to decades.Instead of attempting to model such species-level characteristics, we rely on the simple principle that the overall growth of organisms within a community depends on the availability of energy from net primary production, relative to the total consumption of energy by the metabolic activity of the community.Since one of our primary goals is to predict fishery harvest through coupling with an economic model, we define our community as including all commercially harvested organisms, including pelagic, demersal, and benthic species, both finfish and invertebrates (see discussion of size-based groups in the next section), referring to all as "fish" for simplicity.
In this paper, we describe the ecological module of the BOATS model.In a companion paper (Carozza et al., 2016), the ecological module is coupled to an economic harvesting module and extended to a two-dimensional global grid, in order to explore the spatial distribution of harvest as well as the parameter uncertainty.Here, we present in detail the equilibrium biomass at two ocean sites using a single set of parameter values, and conduct a sensitivity study to illustrate how the model biomass and the size structure of marine communities depend on net primary production and temperature.

Fish ecology model
The ecological module of BOATS uses the McKendrickvon Foerster (MVF) model (McKendrick, 1926;von Foerster, 1959), a widely used continuous-time model for an ageor size-structured population, to represent the evolution of biomass.Populations of fish biomass (all of the organisms in a group) are organized by size and are described by a continuous biomass distribution that we refer to as a biomass spectrum.Fish begin in the smallest size class and grow over time into adjacent (larger) size classes.In each size class, fish biomass evolves in time as the biomass growth less the natural mortality.
Biomass growth is determined by the net primary production that is transferred to fish from phytoplankton at the base of the food web, but cannot exceed the empirical maximum physiological fish growth rates that depend on the individual fish mass and temperature.As such, the local net primary production supports a maximum possible production rate of fish biomass.If actual production within the resolved fish spectra falls below this, due to a shortfall in the availability of biomass that can grow larger, the surplus net primary production is assumed to be taken up outside the resolved fish spectra, by non-commercial species (e.g.non-commercial invertebrates).The natural mortality in each mass class rep-resents biomass losses due to predation, by organisms both within and outside of the community, as well as other natural causes.The mortality formulation depends on an empirical relationship that considers the individual mass of the fish, the asymptotic mass of the fish (the maximum theoretical mass), and the temperature.The addition of new biomass into the smallest mass class, referred to as recruitment, is determined as a function of the net primary production and the production and survival of eggs.
BOATS is designed with the global ocean in mind, yet for ease of reading we present it for a single patch of the ocean, or in other words, for a single grid point on a twodimensional grid.By then applying BOATS to each oceanic grid cell independently, we represent fish biomass and harvest on a two-dimensional global grid.We force biomass using two-dimensional grids of vertically integrated net primary production (NPP) and vertically averaged temperature derived from satellite ocean colour and direct temperature measurements, respectively (Sect.2.8).At each grid point, we therefore simulate biomass spectra that are independent of the adjacent grid points.Hence, we do not take active or passive movement of fish, larvae, or eggs between adjacent grid points into account.These are complex processes that have been shown to play a role in determining fish biomass distributions (Watson et al., 2014).In BOATS, we assume that fish are present where there is NPP to provide food.Given that the model grid cells are 1 • × 1 • , we only effectively ignore nonlocal movements that occur over spatial scales that are larger than approximately 100 km × 100 km.This could bias our results in parts of the ocean where the advection of fish biomass is strong relative to the time step and spatial grid scale, such as in the Gulf Stream.This is especially true for larvae, but would likely pose less of a problem for larger fish since they swim faster than strong oceanic currents.Due to the movement of plankton by currents, a bias could also result from the difference in the locations at which plankton and fish production occur.We expect this to have a small impact on our results given our relatively coarse spatial resolution.Movement induced by ocean circulation and fish behaviour could be implemented in the future, with existing advection and diffusion algorithms (Faugeras and Maury, 2005;Watson et al., 2014).
In the current implementation of the model, we consider three independent populations of fish at every grid point, and so resolve three biomass spectra.These populations, which we refer to as groups, are defined by their asymptotic sizes as small, medium, and large fish, which allows for a very crude representation of biodiversity (Andersen and Beyer, 2006;Maury and Poggiale, 2013).There is no growth from one group to another; in other words, the small group consists of fish that remain small throughout their life history, such as anchovies and sardines, and so are distinct from the juveniles of the medium and large groups.The asymptotic mass, the mass at which all energy is allocated to reproduction and therefore the mass at which growth stops, characterizes each group.We employ groups since they allow us to make use of well-studied growth and mortality characteristics of fish of different asymptotic size (Andersen and Beyer, 2006;Maury and Poggiale, 2013).We work with a finite number of groups as opposed to a continuum (as in Andersen and Beyer, 2006;Maury and Poggiale, 2013), to directly compare our harvest results to the Sea Around Us Project (SAUP) harvest database (Watson et al., 2004;Pauly, 2007), using the three asymptotic masses (Sect.2.9) from the functional group definitions of the SAUP harvest database.Our group formulation combines functional groups (pelagic, demersal, and benthic, for example).Such an assumption may not be appropriate for particular aspects of benthic ecosystems, which have been shown to require more than a representation of size structure to adequately represent core ecosystem features (Duplisea et al., 2002;Blanchard et al., 2009).Nevertheless, for our global-scale model, we feel justified in using such a group formulation since Friedland et al. (2012) found little difference in how the biogeochemical attributes and harvest of pelagic and demersal fisheries reacted to primary production and trophic transfer efficiencies.Alternative group formulations remain a promising avenue of research in global fisheries modelling, one that could be pursued in future work (Blanchard et al., 2009;Maury, 2010).
Although we use the classical MVF model, we implement empirical relationships whenever possible to determine fundamental rates such as growth and mortality, since our goal is to represent fish biomass at the global scale, while limiting the model complexity and number of parameters.As opposed to determining both growth and mortality from explicit predation, as in Maury et al. (2007), Blanchard et al. (2009), Hartvig et al. (2011), and Maury and Poggiale (2013), NPP and the size distribution of phytoplankton set growth rates for all mass classes of fish through a trophic transfer of energy from phytoplankton to fish.To guarantee that growth rates do not exceed realistic values, a von Bertalanffy growth formulation that is based on field observations acts as an upper limit to the growth rate (von Bertalanffy, 1949;Hartvig et al., 2011;Andersen and Beyer, 2013).Mortality is based on an empirical parameterization that depends on mass and asymptotic mass, but also on the constant allometric growth rate in the empirical limit (Gislason et al., 2010;Charnov et al., 2012).
BOATS continues in a tradition of studies that model the global fishery by applying ecological principles to spatially resolved environmental properties.This line of research can be traced to the work of Ryther (1969), who estimated the potential global fish production and harvest based on NPP and simple trophic scaling relationships.More recently, Pauly and Christensen (1995), Chassot et al. (2010), Watson et al. (2013a), andRosenberg et al. (2014) Maury, 2010) was used to study tuna dynamics in the Indian Ocean (Dueri et al., 2012), as well as the impact of climate change on biomass and the spatial distribution of pelagic fish at the global scale (Lefort et al., 2014).Moreover, Blanchard et al. (2009Blanchard et al. ( , 2012) ) considered the impact of future environmental change in large marine ecosystems (LMEs) and Exclusive Economic Zones, while Woodworth-Jefcoats et al. (2012) examined the impact of climate change in three regions of the Pacific Ocean.

Biomass evolution: the McKendrick-von Foerster (MVF) model
The MVF model, a first-order advection-reaction partial differential equation, was first presented by McKendrick (1926) for use in epidemiology, but was later more formally derived for use in the study of cellular systems by von Foerster (1959).Since it provides a natural framework for representing aspects of size dependency and fish life history, and generates biomass spectra that resemble those found in the field (Sheldon et al., 1972;Blueweiss et al., 1978;Brown et al., 2004;Marquet et al., 2005;White et al., 2007), the MVF model has seen a wide variety of applications in marine ecosystems and fisheries.Ecosystem models that have applied the MVF approach to large-scale fisheries studies generally make use of the classical size-structured equation, but differ in the formulations used to calculate growth, mortality, and reproduction, and differ in the structural organisation of fish groups.
Although the MVF model can be expressed by a variety of variables, it is usually presented in terms of the number of fish (the abundance) that evolve in time as a function of the fish age.As an alternative to age, size (measured as length or mass) is also used as an organizing variable, since it can be more descriptive than age for certain applications.Since fish growth (von Bertalanffy, 1949;Andersen and Beyer, 2013), natural mortality (Pauly, 1980;Gislason et al., 2010;Charnov et al., 2012), and harvest (Rochet et al., 2011) are generally size-dependent, we employ size in lieu of age.Moreover, we describe size in terms of mass as opposed to length, although there is a strong relationship between fish mass and length (Froese et al., 2013).
The MVF model uses a spectral framework to describe fish populations; that is, it describes the biomass of fish of mass m at time t by a continuous spectrum f (m, t) such that the fish biomass in the mass interval [m, m + dm] is f (m, t) dm.Although abundance is typically used in applications of the MVF model, and has been used in marine ecosystem applications, see for example Andersen and Beyer (2006); Blanchard et al. (2009), or Datta et al. (2010), we use biomass to compare our results more directly with the harvest data that we use to evaluate BOATS.Regardless, since the abundance n and biomass f spectra are related by f (m, t) = n(m, t)m, in the continuous case, using one form over the other does not influence the model dynamics.We note that, in the numerical implementation of the model, there will be a small difference between the two since we use the geometric mean to represent a discretized range of masses (Sect.2.9).Hence, as fish grow they jump from one geometric mean to next, which may result in an accumulation of biomass.
Fish biomass evolves in time as where f k (m, t) is the biomass spectrum in grams of wet fish biomass (gwB) per square metre of ocean surface per unit of the mass class (gwB m −2 g −1 ), for an individual fish of mass m, at time t, belonging to group k.In Appendix A, we derive the biomass form of the MVF model used in Eq. ( 1).From the definition of the biomass spectrum above, we have that the cumulative biomass at time t of individuals of mass ranging from 0 to m is the integral F k (m, t) = m 0 f k (m , t)dm .In this paper, spectral variables such as the biomass spectra f k (m, t) are written in lower case, whereas cumulative variables that are integrated over size are written in upper case.
Fish biomass is controlled by growth, mortality, reproduction, and recruitment (note that we present harvest in the companion paper, Carozza et al., 2016).Term 1 on the right hand side of Eq. (1) represents the somatic growth in fish biomass that occurs at a rate γ S,k (m, t) (g s −1 ).This term results from fish growing from one interval of mass, which in the discrete case is called a mass class, into the adjacent mass class (for example from a class of 1 to 2 kg to a class of 2 to 3 kg).Since the MVF model is founded on the conservation of numbers of fish (Appendix A), term 2 represents the biomass accumulation that occurs from fish growing in size.Term 3 of Eq. ( 1) represents the natural mortality k (m)f k (m, t) (gwB m −2 g −1 s −1 ), or all non-harvesting sources of fish mortality, which includes losses to predation as well as non-predation losses such as parasitism and disease, senescence, and starvation (Pauly, 1980;Brown et al., 2004).Although we do not consider harvest mortality in this paper, in the full BOATS model (described by Carozza et al., 2016, in review) it is represented by another loss term on the right hand side of equation Eq. (1).The growth rate γ S,k (m, t) (Eq.22) and the mortality rate k (m) (Eq.26) depend on both mass and temperature.
Since the time evolution equation of the MVF model is a first-order partial differential equation, we specify an initial condition (Eq.2) and a boundary condition (Eq.3).The initial condition, or the fish biomass spectrum at the starting time f k,m,0 , is discussed in Sect.3.1.The boundary condition, which is defined at the lower mass boundary m 0 , de- and processes of the ecological module of BOATS.Net primary production (NPP) and temperature (T ) force the model and are used to calculate the fish production spectrum, by assuming a transfer of energy from phytoplankton to successive sizes of fish that depends on the trophic efficiency and the predator to prey mass ratio.From fish production, we calculate the size-dependent growth rate of biomass in three independent groups that represent small, medium, and large commercial fish.Mortality rates are calculated as a function of size and asymptotic size, and also depend on temperature.Adult fish, the largest sizes in each spectrum, allocate energy to reproduction, of which a fraction is returned to the smallest mass class of the corresponding spectrum, representing recruitment of juveniles.
termines the flux of biomass that is added to the biomass spectrum at the initial size class, and depends on the energy allocated to reproductive biomass, the recruitment, and the NPP.This term is detailed in Sect.2.4 and summarized in Eq. ( 29).A schematic of the ecological module of BOATS, with the main model components and processes, is presented in Fig. 1.

Temperature dependence
Organismal body temperature is a fundamental driver of physiological processes since it strongly controls rates of metabolic activity and therefore strongly influences growth, mortality, and reproduction rates (Brown et al., 2004).To model temperature dependence, which we represent by the function a(T ), we apply the van't Hoff-Arrhenius equation where T r (K) is the reference temperature of the process in question (growth or mortality, for example), k B (eV K −1 ) the Boltzmann constant, and ω a (eV) the activation energy of metabolism.Although there is at present no mechanistic derivation of the relationship between metabolic rate and temperature at the level of an entire organism, we interpret the exponential temperature dependence of Eq. ( 4) as an empirical parameterization of this complex relationship with strong observational constraints (Clarke, 2003(Clarke, , 2004;;Marquet et al., 2005;Vandermeer, 2006).For all temperature-dependent rates, we use the average water temperature from the upper 75 m of the water column (Jennings et al., 2008), since it is representative of an average mixed layer depth and so identifies the average temperature at which photosynthesis takes place (Dunne et al., 2005), and since it is representative of the depths at which many fish live and are harvested (Morato et al., 2006;Watson and Morato, 2013).We further assume that fish adopt exactly the water temperature.Given that the greater majority of marine organisms are ectotherms, we feel that this is a reasonable assumption.Taking the average of the upper 75 m of the water column could create biases in regions with strong vertical temperature gradients, since different components of the ecosystem could live at substantially different temperatures, or in regions that are dominated by bottom dwellers in regions deeper than 75 m.However, given that many commercial fish spend significant time near the surface, but actively travel throughout the water column, we feel that this depth is an appropriate first approximation of the average temperature felt by the community.Note that the temperature we apply is generally not accurate for mesopelagic ecosystems, which could make up a large part of marine biomass (Irigoien et al., 2014), but since the majority of these ecosystems have not been commercially exploited, they are not included in our modelled community.

Energy allocation to growth
Fish growth rates are key mass-dependent quantities that characterize each fish group and are limited by the energy available to consumers, and, ultimately, by the photosynthetic primary production.We assume that there is a constant energetic content of biomass (Krohn et al., 1997;Maury et al., 2007), and so treat biomass and energy as equivalents.We envision that energy is supplied to a fish of mass m by the transfer of biomass through the food web by means of predation.Following macroecological theory, this complex process is parameterized by assuming that a fraction of the energy from NPP is transferred up through the food web to become fish biomass production, depending on the average trophic efficiency, the average predator to prey mass ratio, and the phytoplankton size (Ernest et al., 2003;Brown et al., 2004) (Eq. 8).Individual fish then allocate this energy input to either somatic growth (that is, the formation of additional biomass, which we from here on refer to simply as growth γ S,k (m, t), g s −1 ) or to formation of reproductive biomass γ R,k (m, t) (g s −1 ), and so we have that where ξ I,k (m, t) is the input of energy to a fish at time t in group k.We rearrange to write the growth rate as It is important to recognize that the individual fish growth rate cannot exceed a biologically determined maximum rate, no matter how much food is available.This aspect of fish growth is based on empirical observations and allometric arguments, and founded on the work of von Bertalanffy (1938Bertalanffy ( , 1949Bertalanffy ( , 1957) ) and expanded upon by many others including Paloheimo and Dickie (1965), West et al. (2001), andLester et al. (2004).To take this growth rate limitation into account, we assume that the realized input energy ξ I,k (m, t) cannot exceed that supplied by NPP through the trophic scaling, or that determined by empirical growth limits, and so have that the energy input is where ξ P,k (m, t) is the energy input to fish from NPP as transferred through the food web, and ξ VB,k (m, t) is that input from a purely empirical allometric framework following von Bertalanffy (1949).Essentially, ξ VB,k (m, t) describes the maximum growth rate of fish in the case that food is extremely abundant.
We define φ ,C as the fraction of NPP that is potentially available to the sum of all commercial fish groups.In the present work, we assume that φ ,C is equal to 1, and therefore omit it from the equations.This simplifying assumption implies that the entire global ecosystem of animals larger than 10 g would have consisted of potentially commercial species prior to fish harvesting (including bycatch).Obviously this is incorrect, in that the existence of noncommercial animals larger than 10 g requires that φ ,C < 1 in the real world.However, given the weak observational constraints on biomasses of non-commercial animal species at the large scale, and the fact that the species composition of all marine ecosystems has been heavily altered by human activity, it is very difficult to estimate the true value of φ ,C .Despite this difficulty, sensitivity tests revealed that biomass and harvest in the model are approximately linear with φ ,C (not shown).Since we constrain the parameters in BOATS by comparing linear correlations of modelled and observed harvest (see Sect. 3, Table 1, and Carozza et al., 2016), and given the linearity of modelled harvest vs. φ ,C , the value of φ ,C would have a negligible effect on the spatial correlation criterion used for the optimized parameter choices.Nevertheless, it should be pointed out that using alternate values of φ ,C would change the predicted biomass and harvest, all else being equal.
We further assume that each of the three fish groups has access to an equal fraction of the available primary production, φ ,C /3.By assuming that constant portions of the available photosynthetic energy are available to each of the commercial fish groups, all groups are assured to coexist stably.Ecologically, our assumption implies equal resource partitioning to each group, both when they are at the larval stage (through recruitment) and as juveniles and adults (through growth) (Chesson, 2000).This can be thought of as reflecting a separate ecological niche for each group that remains stable over time, and implies that excess NPP, which would result from growth-rate limitation of one group, is not available to other potentially commercial groups, but rather supplied to noncommercial species.Non-commercial species could include, among others, unharvested mesopelagic fish, planktonic invertebrates such as cnidarians, and benthic invertebrates such as amphipods and nematodes.Although this and the previous assumption are not strictly accurate representations of the marine ecosystem, we feel that they are commensurate with the simple three-group representation of the ecosystem and the scarcity of appropriate data constraints, and could be improved in future work.
Each individual fish receives an equal part of the fish production that is input to its mass class, which we here identify as an infinitesimal mass interval of width dm.Where φ C,k is the fraction of φ ,C that is available to group k, and π(m, t) the fish production distribution, the individual fish production is therefore the fish production in the mass interval φ C,k π(m, t)dm divided by the number of individuals in the mass class n k (m, t)dm (Eq.8).Since the abundance spectrum n k (m, t) is equal by definition to f k (m, t)/m, the primary-production-based input of energy to each individual fish is Since we assume that the NPP that is transferred up through the trophic web is uniformly input to all individuals in a given mass class, if the biomass in a mass class falls due to a removal (such as harvesting, for example) then this is equivalent to a decrease in the number of individuals in that mass class.This implies that more fish production would input to each individual, and so in such a scenario ξ P,k would increase.This input of energy depends on the biomass (also referred to as density dependence) and the fish production.The fish production term depends on temperature through the representative mass of phytoplankton m ψ (t) (Eq.25), which is a function of the temperature-dependent large fraction of phytoplankton L (t) (Dunne et al., 2005).
In conditions that are not limited by food availability, the standard von Bertalanffy (somatic) growth rate equation is where the Am b term represents the energy input from food intake after assimilation and standard metabolism, and k a m and k r m represent the energy used in activity and reproduction, respectively (von Bertalanffy, 1949;Paloheimo and Dickie, 1965;Chen et al., 1992;Andersen and Beyer, 2013).The allometric growth rate (not to be confused with the somatic growth rate γ S,k ), which we write as A = A 0 a A (T ), Table 1.Ecological model parameters.Assumption (I) (Brown et al., 2004;Savage et al., 2004;Andersen and Beyer, 2013); assumption (II) value of slope sufficiently large to have abrupt increase in allocation of reproduction from 0 to 1; assumption (III) (Beverton, 1992;Charnov et al., 2012); assumption (IV) (Jennings et al., 2008;Barnes et al., 2010;Irigoien et al., 2014).β truncated since we only consider fish up to 100 kg; assumption (V) Equal partitioning of net primary production to each group; assumption (VI) (Dahlberg, 1979;Andersen and Pedersen, 2010;Pulkkinen et al., 2013).Assumption (VII) (Duarte and Alcaraz, 1989;Cury and Pauly, 2000;Freedman and Noakes, 2002;Maury et al., 2007).The † symbol in the first column identifies parameters that were considered in the tuning procedure of the companion paper (Carozza et al., 2016).∂F /∂p is the rate of change of equilibrium biomass (calculated over the three groups) with respect to change in a parameter p. is the growth constant A 0 (g 1−b s −1 ) modulated by the van't Hoff-Arrhenius temperature dependence for growth a A (T ) (Eq. 4).
The energy input we wish to resolve is that for both growth and reproduction, and so we add the reproduction term k r m to both sides of Eq. ( 9) to find the energy input to be Although the interpretation of the terms in Eq. ( 10) do not exactly correspond to von Bertalanffy's original interpretation of a balance between anabolic growth and catabolic decay, we refer to this equation as the von Bertalanffy energy input ξ VB,k .We consider different values of the activation energy of metabolism for growth ω a,A and mortality ω a,λ (Eq.4), which result in different temperature dependence curves a A (T ) and a λ (T ).The parameter b (unitless) is the allometric scaling constant, and k a (s −1 ) is the mass specific investment in activity.We follow Andersen and Beyer (2013) and define a new constant a = k a /(k a + k r ), which when combined with the idea that there is zero growth at the asymptotic mass m ∞,k (Munro and Pauly, 1983;Chen et al., 1992;Andersen and Beyer, 2013), allows us to express the mass specific investment in activity as k a = A a m b−1 ∞,k .At each group's asymptotic mass, we therefore have that 7) for the input of energy to growth and reproduction is therefore the minimum of a term that depends on biomass and one that does not.Applying the definition of the fish production spectra that we introduce in the next section (Eq.24), we have a change in growth regime when f k is such that When biomass is low enough that this equation holds, NPP no longer influences the input energy, fish will grow at their maximum physiological rate, and any unused energy available to fish production is assumed to be transferred to unresolved parts of the ecosystem.For low productivity systems, the model could overestimate biomass since a larger fraction of primary production will be transferred to commercial species.However, in high productivity systems, the allometric limit is more likely to set growth rates; therefore a larger fraction will be transferred to the non-commercial groups.
That said, the potential for, and the magnitude of, such a feature will depend on the particular values of the growth rates at the site in question (Eq.11).

Energy allocation to reproduction
We assume that the energy allocated to reproduction γ R,k (m, t) is proportional to the total input energy ξ I,k (m, t) where k (m) is the mass-dependent fraction of input energy that is allocated to reproduction.From Eq. ( 6), we write the growth rate as We now derive an expression for k (m).Following Hartvig et al. (2011), we assume that the allocation to reproduction is proportional to mass (Blueweiss et al., 1978;West et al., 2001;Lester et al., 2004;Andersen and Beyer, 2013), and that it also scales with a size-dependent rate s k (m) that defines the size-structure of the transition to maturity (Eq.23).This gives us where k max r is a normalizing constant.Combined with Eq. ( 13), we have that where ξ I,k is a representative input energy that we employ to guarantee that the allocation to reproduction does not change with biomass.For the representative input energy, we take the maximum possible value; that is, the von Bertalanffy input energy described in Eq. ( 10), and so have that ξ I,k = ξ VB,k .
We therefore determine k (m) for the energy input regime that is not limited by fish production, and find that We determine k max r by applying the definition of the asymptotic mass, namely that it is the mass at which energy is only allocated to reproduction and so k (m ∞,k ) = 1.This gives and so we have that We replace this value of k max r into Eq.( 17) to find that Applying Eq. ( 10) for ξ VB,k , and noting that s k (m ∞,k ) is essentially equal to 1, we find that Bringing this development together with Eq. ( 14), the individual fish growth rate is As in Hartvig et al. (2011), we assume that the mass structure of the allocation of energy to reproduction s k (m) is a sharply transitioning function that shifts from near zero to near one around the mass of maturity m α,k .Based on Beverton (1992) and Charnov et al. (2012), we further assume that the mass of maturity is proportional to the asymptotic mass m ∞,k such that m α,k = ηm ∞,k (Table 1).Although other functional forms are plausible, s k (m) must have a transition in mass that is proportional to m ∞,k (or to the maturity mass) (Hartvig et al., 2011), and so we use the functional form used by Hartvig et al. (2011), where the parameter c s determines how quickly the transition from zero to one takes place (Fig. 2).For reference, we calculate the reproduction allocation mass scale, the range over which the majority of the change in reproduction allocation takes place, as the inverse of the derivative evaluated at the maturity mass, ( ds k dm | m=m α,k ) −1 , which we find to be 4m α,k c s .

Fish production spectrum
We model the biomass production of fish by assuming that both phytoplankton and fish production are part of the same energetic production spectrum (Sheldon et al., 1972;Ernest et al., 2003;Brown et al., 2004).Unlike in the approaches of Maury et al. (2007) and Hartvig et al. (2011), among others, we do not model the growth and decay dynamics of phytoplankton biomass.Instead, we represent fish production over a spectrum of individual fish masses, π(m, t) Representative mass of phytoplankton g ( 24), ( 25) Fraction of large phytoplankton production Unitless (25) R P (m 0 , t) Primary-production-determined recruitment gwB m −2 s −1 (27) R e,k (m 0 , t) Egg production and survival determined recruitment gwB m −2 s −1 (28) R k (m 0 , t) Overall recruitment gwB m −2 s −1 (29) (mmol C m −2 g −1 s −1 ).Following Brown et al. (2004) and Jennings et al. (2008), we base this formulation on (1) the NPP ψ (t) (mmol C m −2 s −1 ) (Sect.2.8), (2) the representative size at which NPP takes place m ψ (t) (g) (Jennings et al., 2008), and (3) the trophic scaling exponent τ that indicates how efficiently energy is transferred through the trophic web, where τ depends on the trophic efficiency α and the predator to prey mass ratio β, and is equal to log(α)/ log(β) (Brown et al., 2004).The fish production spectrum follows As in Brown et al. (2004), we assume that α and β, and hence τ , are constant.From the expression for fish production detailed in Eq. ( 24), we determine the individual fish growth rate using Eq. ( 22).
Although variability in the trophic scaling τ , that could depend on environmental or ecosystem characteristics, is potentially of significant importance, we take here the simple assumption that the trophic scaling is globally constant, as other authors have (Brown et al., 2004;Jennings et al., 2008).We note that, using a large database of individual prey eaten by individual predators, Barnes et al. (2010) found that the predator to prey mass ratio increases with predator mass.Given that we apply an average value of β, and assuming that all else remains equal, the work of Barnes et al. (2010) implies that we would underestimate β for large m and overestimate β for small m, and so (by Eq. 24) we underestimate D. A. Carozza et al.: BOATS-1.0 π k for large m and overestimate π k for small m.Essentially, a mass-dependent β would tend to decrease the steepness of biomass spectra relative to what is shown here.It is also commonly assumed that the trophic efficiency α is constant (Brown et al., 2004;Jennings et al., 2008;Tremblay-Boyer et al., 2011).Based on acoustic biomass estimates and modelling work, Irigoien et al. (2014) suggests that trophic efficiency can instead be significantly different in low and high productivity regions, at different levels in the food web (from phytoplankton to mesozooplankton and from mesozooplankton to fish) and that it can also depend on environmental parameters such as temperature (through its influence on organismal metabolic rates) and water clarity (which affects visual predation).Quantifying variability in τ is an important target for future work.
The production spectrum is the product of two terms.The first is the initial value determined at the representative phytoplankton mass m ψ (t), which corresponds to the NPP normalized by the representative phytoplankton size.The fish production spectrum then follows a power law dependence in m with a scaling exponent of τ − 1.This mass scaling represents larger phytoplankton (larger m ψ (t)) being trophically closer to fish than smaller phytoplankton, thereby permitting more energy to be transferred from phytoplankton to fish (Ryther, 1969).The power law dependence that we use is based on Kooijmann (2000) and Brown et al. (2004).The model is forced with observations of NPP, and so we run the model in units of mmol C. For analysis and presentation, we convert to grams of wet biomass (gwB) by assuming that there are 12 g C per mol C, and that there are 10 gwB for every g of dry carbon (Jennings et al., 2008).
Phytoplankton mass ranges over several orders of magnitude (Jennings et al., 2008).We take a simple approach and express the spectrum of phytoplankton as a single representative mass at which NPP takes place.Due to the wide range of phytoplankton mass, we calculate the representative mass as and so take the geometric mean of the mass of a typical large, m L , and a typical small, m S , phytoplankton, weighted by the fraction of production due to large or small phytoplankton, L (t) and 1 − L (t), respectively.We calculate this fraction using the phytoplankton size structure model of Dunne et al. (2005), which resolves small and large phytoplankton and assumes that small zooplankton are able to successfully prey upon increasing production of small phytoplankton, but that large zooplankton are unable to do so as effectively for large phytoplankton production.Dunne et al. (2005) propose an empirical relationship for the large fraction of NPP L (t) in terms of temperature T C (t) ( • C) and the NPP, the Eppley factor e k E T C (t) where k E ( • C −1 ) is the Eppley temperature constant for phytoplankton growth, and * (mmol C m −3 d −1 ) the productivity normalized to a temperature of 0 • C. The Dunne et al. (2005) model resolves a high fraction of the variability in phytoplankton community structure (Agawin et al., 2000), and provides a mechanism to explain how the fraction of large phytoplankton biomass increases with increasing phytoplankton biomass.Although we use this particular formulation for the large fraction in Eq. ( 25), future work could examine alternatives (Denman and Pena, 2002).

Natural mortality
The natural mortality term represents all forms of natural (non-fishing) mortality.It mainly consists of predation, but also includes non-predatory sources of mortality such as parasitism, disease, and senescence (Pauly, 1980).This term is of first-order importance in determining energy flows in marine food webs, and so also in determining biomass.In pursuing our principle of using empirical parameterizations to represent complex processes that are incompletely understood, we follow the work of Gislason et al. (2010) andCharnov et al. (2012) and take the mortality rate to be where λ = e ζ 1 (A 0 /3)a λ (T ) (see Appendix B for a full derivation of this form).ζ 1 is a parameter estimated from mortality data (Gislason et al., 2010), A 0 (g 1−b s −1 ) is the growth constant from Eq. ( 10), and a λ (T ) is the van't Hoff-Arrhenius exponential for mortality as described in Eq. ( 4).Charnov et al. ( 2012) provided a mechanistic underpinning for Eq. ( 26) by calculating the optimal number of daughters per reproducing female over that female's lifetime.Unlike other empirical mortality rate frameworks, such as that of Savage et al. (2004), the mass dependence m −h does not depend on the allometric growth scaling b, and so the mass dependence of the mortality rate is not determined by internal biological parameters, but by predation and competition (Charnov et al., 2012).The losses due to natural mortality, term 3 in Eq. ( 1), are linearly proportional to biomass as in Gislason et al. (2010), and in keeping with the classical MVF model.
It is important to highlight the fact that unlike some other models, we do not adopt an explicit representation of predation-dependent mortality (Maury et al., 2007;Blanchard et al., 2009;Hartvig et al., 2011).The mortality rate only depends on the organism mass, asymptotic mass, and temperature, and is linear in biomass.This choice is motivated by the wide range of predator-prey mass ratios in marine ecosystems (Barnes et al., 2010), and the complexity and non-stationarity of food web relationships.In applying this parameterization, we avoid the complication of choosing a difficult-to-constrain prey selectivity function, and benefit from applying mortality rates that are directly founded in observed rates.Without necessarily losing realism, this parameterization simplifies the complicated dynamics that result from more sophisticated prey selectivity formulations (Andersen and Pedersen, 2010).
Since the prey mortality rate does not depend on the predator biomass, we do not resolve top-down trophic cascades (Andersen and Pedersen, 2010;Hessen and Kaartvedt, 2014).At present, a scarcity of data hinders a formal verification of generalized trophic cascades in the open ocean, which would be desirable for the formulation of their impact within the BOATS framework.However, we do represent bottom-up effects through the growth formulation described in Eq. ( 1), since a change in biomass in one size class is carried upward through the trophic web as fish grow to larger mass classes.

From reproduction to recruitment
Fish reproduction and recruitment comprise a set of complex ecological processes that result in new fish biomass entering a fishery (Myers, 2002).This first involves fish allocating energy to reproduction and releasing eggs and sperm during spawning.Fertilized eggs must then survive predation until they hatch to become larvae, when they must again survive predation until they grow into juveniles (Dahlberg, 1979;McGurk, 1986;Myers, 2001).The end of the juvenile stage is generally defined as when fish reach sexual maturity or when they begin interacting with other adult members of the fishery (Kendall et al., 1984).The definition of a recruit is more nuanced since it generally depends on the fishery in question and can be based on a particular size or age, the size or age of sexual maturity, or the size or age at which fish can be caught (Myers, 2002).For the model, we refer to recruitment as the flux of new biomass into the lower boundary mass (m 0 ) of each group.
Recruitment is driven by biomass-dependent (densitydependent) processes, such as predation and disease, as well as by biomass-independent (density-independent) processes such as environmental change.These processes strongly and nonlinearly affect mortality throughout the egg, larval, and juvenile stages (Dahlberg, 1979;McGurk, 1986;Myers, 2002).To model the number of recruits that result from a given spawning stock of biomass, one must make assumptions on the nature of these processes.The widely used stockrecruitment models of Ricker (1954), Schaefer (1954), and Beverton and Holt (1957), and the generalization of these models by Deriso (1980) and Schnute (1985), make such assumptions and operate in terms of the spawning stock biomass; that is, the biomass that is of reproductive age.
We model recruitment by considering both the NPP and the production and survival of eggs by adult fish.Our formulation is based on the Beverton-Holt stock recruitment relationship (which employs a Holling Type 2 functional form, Holling, 1959), as used by Beverton and Holt (1957) and Andersen and Beyer (2013), with NPP setting the upper limit and the half-saturation constant (Eq.29).This form allows for an approximately linear decrease to zero recruitment as the spawning stock biomass goes to zero, but sets an upper limit that depends on the NPP when the spawning stock biomass is large, in order to represent the role of food availability in determining larval survival.
The flux of biomass out of a mass class is the growth rate multiplied by the biomass in that mass class (Eq.1).Since the recruitment is also a flux of biomass (one that occurs at the lower mass boundary), to define it in terms of NPP R P,k (m 0 , t) (gwB m −2 s −1 ), we apply Eq. ( 8) and find that where m 0 is the lower bound of the smallest mass class, and π is the fish production spectrum from Eq. ( 24).Alternatively, the recruitment from the production and survival of eggs to recruits, R e,k (m 0 , t) (gwB m −2 s −1 ), depends on the energy allocated to reproduction, γ R,k (t) (Eq.13), by all n k individuals over all mass classes, which we write as The model biomass includes both males and females, which are assumed to mature at the same mass (Beverton, 1992).As in other model studies (Maury et al., 2007;Andersen and Pedersen, 2010;Andersen and Beyer, 2013), males and females of reproductive age continually reproduce, yet only the female contribution is counted in the flux into the smallest mass class, since the male contribution to a fertilized egg is negligible compared to that of the female.Hence, when the integral part of Eq. ( 28) is multiplied by the fraction of females, φ f , we have the biomass of eggs produced.Dividing by the mass of an egg m e therefore gives the number of eggs produced, which when multiplied by the survival fraction s e , expressing the probability that an egg becomes a recruit, gives the number of recruits.From the number of recruits produced per unit time, we multiply by the mass of a recruit, m 0 , to determine the biomass flux of recruits.
Applying the same form as the stock-recruitment model developed by Beverton and Holt (1957) (see Andersen and Beyer, 2013) we take the overall recruitment R k (m 0 , t) Following Andersen and Beyer (2013), we take the halfsaturation constant (the value of R e,k (m 0 , t) at which the overall recruitment is one half of the maximum recruitment allowed by productivity) to be R P,k (m 0 , t). Figure 3 shows how the overall recruitment R k (m 0 , t) changes as a function of R P,k (m 0 , t) and R e,k (m 0 , t).As is the case for a Holling Type 2 functional form, as biomass and therefore also the egg-and survival-based recruitment R e,k (m 0 , t) increases, the overall recruitment saturates toward the primary production-based limit R P,k (m 0 , t).This indicates that for sites with high biomass, NPP limits recruitment.At the other extreme, when R e,k (m 0 , t) is small relative to R P,k (m 0 , t), the recruitment is approximately linear in R e,k (m 0 , t) and so has a weak dependence on R P,k (m 0 , t) such that at low biomass the egg production and survival limits recruitment.
Tables 1 and 2 detail the fish model parameters and variables, respectively.The group and mass class structure, and the numerical discretization of the continuous biomass spectra, are presented in Sect.2.9 and Sect.2.10, respectively.

Environmental forcing: temperature and net primary production
The ecological model requires temperature and NPP information as forcing input to calculate the time evolution of biomass (Eq.1).These variables can be provided by an ocean general circulation model that includes a lower trophic level model.Here, we instead use observational estimates, which would be expected to provide a more realistic simulation.For temperature, we use the World Ocean Atlas 2005 (Locarnini et al., 2006), which brings together multiple sources of in situ quality-controlled temperature interpolated to monthly climatologies on a 1 • ×1 • grid.We discuss our usage of temperature in Sect.2.2, and as discussed above, use the average water temperature from the upper 75 m of the water column to force temperature-dependent rates.For NPP, we take the average of three satellite-based estimates (Behrenfeld and Falkowski, 1997;Carr et al., 2006;Marra et al., 2007) to capture some of the variability that exists in different NPP models (Saba et al., 2011).We note that satellite-based estimates suffer from a range of shortcomings, including lack of pro-ductivity sources other than phytoplankton (e.g.seagrass and corals), and biases in coastal regions and estuaries (Smyth, 2005;Saba et al., 2011).Although overall minor (see Crossland et al., 1991;Duarte and Chiscano, 1999), these uncertainties will carry through to the modelled biomass and harvest.

Group and mass class structure
Fish span several orders of magnitude in mass, and we therefore discretize the mass spectra into logarithmic mass classes.In order to directly compare our results with the Sea Around Us Project (SAUP) harvest database (Watson et al., 2004;Pauly, 2007), we consider three fish groups each with a different asymptotic mass.We first convert the maximum lengths used in the SAUP (30 cm for the small group, 90 cm for medium group, and up to our maximum resolved length for the large group) to asymptotic length assuming that the maximum length is 95 % of the asymptotic length (Taylor, 1958;Froese and Pauly, 2014), and then apply a length-weight relationship of the form m = δ 1 l δ 2 (Froese et al., 2013) to calculate the asymptotic mass.This results in asymptotic masses of 0.3, 8.5, and 100 kg for the small, medium, and large groups, respectively.
Although the asymptotic masses differ, all three groups have the same mass class structure, with lower and upper bounds of m 0 = 10 g and m u = 100 kg, respectively.Since the groups have different asymptotic masses m k,∞ , there are therefore fewer resolved mass classes for groups with smaller asymptotic mass.We define the mass classes by dividing the mass spectrum into N M classes with lower bounds m i,L such that where i is the index of the mass class that ranges from 1 to N M .Based on this definition, we describe a mass class as an interval We divide the spectrum into 50 mass classes (N M = 50).Although we use fewer mass classes than some other studies (Maury et al., 2007;Hartvig et al., 2011), we have tested higher temporal and spatial resolutions and find that our interpretations would not be influenced by our choice of temporal or spatial resolution.
When we calculate and present mass-dependent quantities, we consider a mass m i that represents the average or central value of its class.For this, we apply the geometric mean of the lower and upper bounds of a mass class, which we calculate as since the upper bound of a mass class is the same as the lower bound of the adjacent class.

Numerical methods
The biological part of our model is a system of three nonlinear first-order (in mass) partial differential equations that describe the evolution of the biomass spectra of three fish groups.Each equation is forced with the same net primary production and temperature information, and the equations do not interact with one another.Here, we use the standard notation of a subscript i to describe a mass cell, and a superscript n to describe a temporal cell.The notation k, as in the main text, refers to a fish group.For example, f n+1 k,i represents the biomass spectral value f of group k, at mass class i at time n + 1.
Since the McKendrick von-Foerster model is an advective equation in biomass, as is true of advective equations, transport errors are a concern (Press et al., 1992).To limit such errors, and because growth is always defined to be positive (or zero), we apply an upwind scheme (Maury et al., 2007;Hartvig et al., 2011).This numerical scheme uses only biomass information that is upwind of the cell of interest; that is, it only uses biomass information at cells i and i − 1 to integrate and determine the biomass at cell i at the next timestep.We use a forward difference scheme for the temporal rate of change, and explicitly calculate the growth (γ ) and mortality ( ) rates; that is, we use the current temporal state of biomass f n k,i to update the biomass, as opposed to using the future biomass state f n+1 k,i as in an implicit scheme, and integrate biomass as The model is stable and converges as we decrease t.

Results and discussion
Here we describe the behaviour of the fish ecology model, and make use of a simplified version of the model as a reference point and initial biomass condition.We consider two model grid points that correspond to individual patches of ocean at a cold-water site in the East Bering Sea (EBS) LME (64 • N, 165 • W) and a warm-water site in the Benguela Current (BC) LME (20 • S, 12 • E), and describe the resulting biomass spectra and other model variables.We discuss the results from a sensitivity test that considers the role of NPP (ranging from 50 to 2000 mg C m −2 d −1 ) and temperature (ranging from −2 to 30 • C) on biomass.For these simulations, we use a 15 day timestep and constant forcing of annually averaged NPP and temperature.We do not use these sites for a thorough data-based model validation, which is difficult at this time due to a lack of suitable fish biomass data.The parameter values used here  are taken from an extensive data-model comparison that employs the global implementation of the model, and is fully described in the companion paper (Carozza et al., 2016).In that study, we take a Monte Carlo approach with over 10 000 parameter sets to find parameter combinations that best fit observed harvest at the LME-scale, considering the full range of the uncertain parameter space for the 13 most important parameters.Of these 13 parameters, 2 are economic, with the remaining 11 ecological parameters being identified with a dagger symbol in Table 1.Beyond the validation to harvest at the LME-scale in the companion paper, more specific validation could be done in the future with suitable data sets when they become available (that is, size aggregated, regional-scale, species-comprehensive biomass assessments).

Initial biomass state
To begin our results and analysis section, we make a series of simplifying assumptions in order to derive an analytical biomass spectrum f k,m,0 , which we use as a reference point for evaluating aspects of the full model.Since this analytical biomass state is a reasonable approximation of the full model, we also use it as an initial biomass condition for our simulations.
Beginning with the evolution of biomass in Eq. ( 1), we assume that the input energy expressed in Eq. ( 7) is solely controlled by NPP, so that ξ I,k (m, t) = ξ P,k (m, t) = φ C,k π(m, t)m/f k (m, t), and that there is no allocation of energy to reproduction, so that k (m) = 0.These two assumptions result in a growth rate of γ S,k (m, t) = φ C,k π(m, t)m/f k (m, t), which allows us to calculate the equilibrium biomass spectrum ( ∂ ∂t f k (m, t) = 0) in terms of the fish production spectrum (Eq.24) and the mortality rate (Eq.26).We consider constant forcing and so apply the annual average NPP ψ and temperature (which are contained in the mortality rate λ and representative phytoplankton mass m ψ terms), and find that the equilibrium biomass spectrum of each each group is As expected from the MVF model, biomass follows a power law spectrum with respect to mass.Given that the power law scaling exponent is τ + h − 1, biomass scales as a function of the trophic and mortality scalings, which we assume are constant.On the other hand, the intercept of the spectrum (in logarithmic space, when m = m 0 = 10 g) depends on a variety of parameters such as the NPP and trophic efficiency, as well as the natural mortality rate and the representative phytoplankton mass.Unlike the mass scaling, the intercept is also group dependent through the fraction of primary production allocated to each group and the asymptotic mass.

Biomass equilibrium
As in other studies, we use features of the modelled biomass spectra, shown in Fig. 4, to interpret the model results.Work on marine ecosystems indicates that biomass spectra, when plotted in log-log space, are approximately linear over most of the size range and have slopes that range from −1.0 to −1.2 (Blueweiss et al., 1978;Brown et al., 2004;Marquet et al., 2005;White et al., 2007).Ignoring harvest, group biomass spectra generally decrease with size, except at the maturity mass at which energy begins to be allocated to reproduction (Fig. 2), where there is a decrease in the growth rate and thereby an accumulation of biomass that may result in a local maximum or a local decrease of the spectrum slope (Andersen and Beyer, 2013).As expected from Eq. ( 33), the group intercepts differ, but only by little since in our formulation the only difference arises from the weak asymptotic mass dependence m h+b−1 ∞,k in the mortality term.Biomass is larger at the cold-water site, despite it having a lower NPP (Fig. 5).In particular, large group biomass is larger at the cold-water site, which is consistent with the findings of Watson et al. (2014).
There is a nonlinear decrease in biomass at larger mass classes (Fig. 4).The shape of the biomass spectra are determined from the growth and mortality rates.Since the growth rate consists of NPP and allometric regimes (Eq.22), and the mortality rate of a single regime (Eq.26), any changes in the shape of the biomass spectra are determined by the growth rate.We generally find that the NPP regime (Eq.8) limits energy input in smaller mass classes, whereas the allometric regime (Eq.10) plays the limiting role in the largest mass classes.

Sensitivity tests
Total biomass (Fig. 5a) increases monotonically for increasing NPP, yet decreases monotonically for increasing temperature.Increasing temperature not only reduces the primaryproduction-based growth rate γ P by reducing the representative phytoplankton size (Eq.24), it also significantly drives up the mortality rate, generating a clear pattern of reduced biomass.Under the allometric regime of growth (Eq.10), higher temperature implies a greater growth rate, which on its own results in an increase in biomass (not shown).However, this feature is more than counterbalanced by the mortality rate increase, which results in an overall lower biomass for higher temperature.
We calculate the total biomass spectrum as the sum of the biomass of each mass class over all groups.We use the biomass value at the first mass class to define the intercept, and calculate the slopes based on the non-reproducing parts of the spectra (the mass classes that are smaller than the maturity mass m α,k ) since this is generally the linear part of the spectra (Maury and Poggiale, 2013), using linear regression on the log-transformed data (Xiao et al., 2011).The spectral intercept (Fig. 5b) depends on both NPP and temperature, monotonically increasing with increasing NPP, but nonlin-early changing in temperature due to the multiple sources of temperature dependence in the intercept (Eq.33).The biomass slope does not depend on NPP (Fig. 5c), as indicated in Eq. ( 33), and the resulting total slope values (grey curve in Fig. 5c), given the parameters used in this single realization of the model, are consistent with published values from marine ecosystems that range from −1.0 to −1.2 (Blueweiss et al., 1978;Brown et al., 2004;Marquet et al., 2005;White et al., 2007).However, we find flatter slopes for lower temperatures, to values as low as −0.9.This implies that our model would result in generally higher biomass than if the slope of the spectra fell between −1 and −1.2.Equation ( 33) also indicates that the slope is not a function of temperature.That equation applies for the small group (blue curve in Fig. 5c) over all temperatures, and for the medium group at low temperatures.However, when the input energy is determined by the von Bertalanffy limit, as is the case for high temperatures in the medium group and all temperatures for the large group, a rise in temperature steepens the biomass slope.Overall, NPP only influences spectra by shifting the intercept, whereas temperature both shifts the intercept and changes the slopes of biomass spectra when the input energy is set by the von Bertalanffy limit.
The model illustrates hypothetical inferences, based on the macroecological theory it uses, that need to be compared to suitable observations.Further validation of the model at specific locations and at the size-class level of detail remains a challenge because of the scarcity of suitable data sets.To further validate BOATS and comparable models, we require size-class-resolved observations at the ecosystem level, at a high enough resolution to detect variations in spectral properties, and at a sufficient number of sites so as to detect bulk variations due to different temperature and NPP.This type of detail at the ecosystem level is not available even in current stock assessment databases, and it should be considered an important target for future data syntheses.

Conclusions
We have described a new marine upper trophic level model for use in gridded, global ocean models.The model as described here is used as the ecological module of the BOATS model, designed to study the global fishery.In a companion paper, we discuss the economic module of the BOATS model and complete the model evaluation by comparing harvest simulations to global harvest observations.The approach could be readily adapted to other purposes, such as for use in studies of ocean biogeochemistry or ecology.
The model uses NPP and temperature to represent the firstorder features of fish biomass using fundamental marine biogeochemical and ecological concepts.When possible, we apply empirical relationships with mechanistic underpinnings to simplify complex ecological processes that are difficult to constrain.Phytoplankton community structure is represented by the proportion of large phytoplankton.Fish growth rates are determined by a parameterized trophic transfer of energy from primary production, but limited by empirical allometric estimates.The natural mortality rate is based on an empirical relationship that depends on the individual and asymptotic mass, and reproduction depends on the NPP and the fish biomass of reproductive age.The resulting biomass spectra, as defined here, include all commercially harvested organisms longer than 10 cm (greater than 10 g).
We presented simulated biomass spectra at a warm-and a cold-water site, and performed a sensitivity test of the model forcing variables to examine key model variables.We find that the structure of modelled biomass spectra is broadly consistent with observations, and biomass slopes match observations over a wide range of NPP and temperature.Although the model employs a limited number of parameters compared to similar modelling efforts, it retains reasonably realistic representations of biological and ecological processes, and is computationally efficient, which allows for extensive sensitivity studies and parameter-space analyses even when implemented globally.Due to its dynamical generality and conceptual simplicity, the ecological module of BOATS is well-suited for global-scale studies where the resolution of species or functional-groups is not necessary.

D. A. Carozza et al.: BOATS-1.0 Appendix A: Biomass version of the McKendrick-von Foerster (MVF) model
The MVF model equation is an expression of the conservation of the number of fish (Kot, 2001), and in terms of abundance is written as where γ (m, t) is a characteristic velocity of growth (Kot, 2001), which we assume is equivalent to the individual growth rate dm dt , and (m, t) is the instantaneous natural mortality rate.For ease of reading, we ignore the mass and time dependencies and write f = f (m, t), γ = γ (m, t), = (m, t), and n = n(m, t).The biomass spectrum f (m, t) is defined as n(m, t)m, and so n(m, t) = f (m, t)/m.Substituting this expression into Eq.(A1), we have that This result is similar in structure to its abundance-based counterpart in Eq. (A1), aside from the extra term γf m , which is equivalent to γ n.This new term is a direct consequence of the conservation of the number of fish written in terms of biomass, and represents the increase in biomass that occurs as a given number of fish grow into a larger mass interval at the rate γ .

Appendix B: Derivation of natural mortality formulation
We apply the empirical model of natural mortality from Gislason et al. (2010) to derive Eq. ( 26).The natural mortality rate is model 2 of Table 1 from Gislason et al. (2010), where is the natural mortality rate, l is the organism length, l ∞ is the asymptotic organism length, K is the von Bertalanffy growth parameter that is equivalent to A 3 m b−1 ∞ , and T is temperature.The variable A = A 0 a λ (T ) is the growth constant A 0 scaled by the van't Hoff-Arrhenius exponential function for mortality, and b is the allometric scaling constant (Eq.10).Gislason et al. (2010) found that the ζ 4 parameter was not statistically significant, and so we rewrite the natural mortality rate ignoring the temperature term as We apply the length-weight relationship l = (m/δ 1 ) 1/δ 2 taking δ 2 = 3 (Froese et al., 2013) to write the equation in terms of mass, and find that

Figure 1 .
Figure1.Schematic diagram of the main modules, components, and processes of the ecological module of BOATS.Net primary production (NPP) and temperature (T ) force the model and are used to calculate the fish production spectrum, by assuming a transfer of energy from phytoplankton to successive sizes of fish that depends on the trophic efficiency and the predator to prey mass ratio.From fish production, we calculate the size-dependent growth rate of biomass in three independent groups that represent small, medium, and large commercial fish.Mortality rates are calculated as a function of size and asymptotic size, and also depend on temperature.Adult fish, the largest sizes in each spectrum, allocate energy to reproduction, of which a fraction is returned to the smallest mass class of the corresponding spectrum, representing recruitment of juveniles.

Figure 2 .
Figure 2. Mass dependence of reproduction by group.The mass scaling function s k (m) (thin lines, Eq. 23) determines the mass dependence of the allocation of energy to reproduction.k (m) (thick lines, Eq. 21) is the fraction of energy allocated to reproduction.

Figure 4 .
Figure 4. Steady state biomass spectra at two sites.Black solid, dashed, and dash-dot curves represent the small, medium, and large group biomass, respectively, whereas the grey curves represent the total of the three groups.The model is forced at two sites with annual average net primary production (NPP) and temperature (T ) with a timestep of 15 days.Simulations are for a (a) cold-water site in the East Bering Sea LME (64 • N, 165 • W) and a (b) warm-water site in the Benguela Current LME site (20 • S, 12 • E).

Figure 5 .
Figure5.Model sensitivity to net primary production (NPP) and temperature (T ).(a) Total biomass in terms of NPP and T , (b) intercept of total fish spectrum in terms of NPP and T , and (c) group and total slope of the non-reproducing part of the fish biomass spectra.In (c), since the slopes of the biomass spectra do not depend on NPP, the slopes are lines that depend only on temperature.Red and blue circles in (a) and (b) represent the NPP and T of the warm-and cold-water sites, respectively, used in Fig.4.All total spectral intercepts and slopes are calculated by adding the biomass in each mass class over all three groups.The intercept is the spectral biomass of the first mass class, and the slope is calculated from the mass classes that are smaller than the maturity mass m α,k (the non-reproducing mass classes).