Global Change Biology Global Change Biology (2014) 20, 1559–1584, doi: 10.1111/gcb.12489
Mechanistic insights into the effects of climate change on larval cod T R O N D K R I S T I A N S E N 1 * , C H A R L E S S T O C K 2 , K E N N E T H F . D R I N K W A T E R 1 and ENRIQUE N. CURCHITSER3 1 Institute of Marine Research, Bergen, 5817 Norway, 2NOAA Geophysical Fluid Dynamics Laboratory, Princeton University, Forrestal Campus, 201 Forrestal Road, Princeton, NJ 08540-6649, USA, 3Department of Environmental Sciences/Institute of Marine and Coastal Sciences, Rutgers University, New Brunswick, NJ 08901 USA
Abstract Understanding the biophysical mechanisms that shape variability in fisheries recruitment is critical for estimating the effects of climate change on fisheries. In this study, we used an Earth System Model (ESM) and a mechanistic individual-based model (IBM) for larval fish to analyze how climate change may impact the growth and survival of larval cod in the North Atlantic. We focused our analysis on five regions that span the current geographical range of cod and are known to contain important spawning populations. Under the SRES A2 (high emissions) scenario, the ESMprojected surface ocean temperatures are expected to increase by >1 °C for 3 of the 5 regions, and stratification is expected to increase at all sites between 1950–1999 and 2050–2099. This enhanced stratification is projected to decrease large (>5 lm ESD) phytoplankton productivity and mesozooplankton biomass at all 5 sites. Higher temperatures are projected to increase larval metabolic costs, which combined with decreased food resources will reduce larval weight, increase the probability of larvae dying from starvation and increase larval exposure to visual and invertebrate predators at most sites. If current concentrations of piscivore and invertebrate predators are maintained, larval survival is projected to decrease at all five sites by 2050–2099. In contrast to past observed responses to climate variability in which warm anomalies led to better recruitment in cold-water stocks, our simulations indicated that reduced prey availability under climate change may cause a reduction in larval survival despite higher temperatures in these regions. In the lower prey environment projected under climate change, higher metabolic costs due to higher temperatures outweigh the advantages of higher growth potential, leading to negative effects on northern cod stocks. Our results provide an important first large-scale assessment of the impacts of climate change on larval cod in the North Atlantic. Keywords: climate change, cod, earth system modeling, fish, Gadus morhua, individual-based modeling, mechanistic modeling, phytoplankton production, recruitment, zooplankton Received 30 August 2012 and accepted 19 November 2013
Introduction Global climate change is expected to alter the physical environment as well as the biological structure and functioning of the world’s ocean ecosystems (Rosenzweig et al., 2007; Doney et al., 2012). Changes in physical features, such as ocean temperature, stratification, and currents, will have considerable impacts on marine ecosystems (Hays et al., 2005; Doney et al., 2012). For example, projected reductions in primary productivity (Steinacher et al., 2010) could have effects throughout the foodweb including the recruitment of fish stocks (Brander, 2007a). Climate change is also expected to *Present address: National Oceanic and Atmospheric Administration, Silver Spring, Maryland, USA Correspondence: Trond Kristiansen, tel. +4797701109, fax +4755238531, e-mail: [email protected]
© 2013 John Wiley & Sons Ltd
affect fish stocks by causing major geographic shifts in species distribution and abundance over the next 50–100 years (Barker & Knorr, 2007; Brander, 2007b; Cheung et al., 2009b). Together these changes could have major direct and indirect impacts on fisheries recruitment, although considerable uncertainty exists with respect to the magnitude of these changes and the mechanisms that underlie them (Brander, 2007a). Observations and climate model simulations suggest that climate change impacts in the North Atlantic will include several physical changes that could affect biological productivity. Observed sea surface temperatures in the North Atlantic over the twentieth century indicate a warming trend everywhere with exception of the far northwestern Atlantic (Deser et al., 2010b). Climate model simulations from the fourth IPCC assessment report project these trends to continue under high CO2 emissions scenarios (Steinacher et al., 2010; Capotondi et al., 2012). Warming and decreased surface salinity 1559
1560 T . K R I S T I A N S E N et al. are projected to increase surface ocean stratification over most of the North Atlantic during the second half of the twenty-first century (Capotondi et al., 2012). Enhanced stratification limits the vertical exchange of cold nutrient-rich water to the surface layer and is projected to decrease primary productivity in most areas of the North Atlantic (Steinacher et al., 2010). The prominent exceptions to this trend occur in extremely highlatitude ecosystems, where increased stratification and warming is projected to increase primary productivity by alleviating strong light limitation and reducing ice cover (Bopp et al., 2001; Steinacher et al., 2010; Doney et al., 2012). These physical changes can have diverse impacts on cod recruitment across the North Atlantic basin. Atlantic cod (Gadus morhua) is a benthic fish species with a pelagic larval stage that is distributed across the North Atlantic in more than 25 different cod populations (Drinkwater, 2005). Past natural climate variability in the North Atlantic basin has caused oscillating cold and warm conditions in cod spawning grounds (Drinkwater, 2006; Sundby & Nakken, 2008). Warm conditions have been found to have a positive effect on survival and recruitment of Atlantic cod through enhanced zooplankton productivity, faster growth rates, and elevated survival (Ottersen et al., 2001). However, the response to ocean warming may depend on the geographical location of the population, the availability of prey resources, and the species involved (Drinkwater, 2005; Cheung et al., 2009a, 2013; Kristiansen et al., 2011). For example, increased ocean temperatures have resulted in stronger recruitment of Barents Sea stocks of cod, haddock, and herring (Ottersen & Loeng, 2000), while increased temperatures in the North Sea have had negative effects on cod recruitment through changes in the prey composition, distribution, and abundance (Beaugrand et al., 2003). In the Barents Sea, warmer temperatures were associated with an increased inflow of water rich in zooplankton, which increased growth and survival through the larval phase (Sundby, 2000). However, in the North Sea, warmer temperatures caused a shift in zooplankton species composition that was lower quality for larval fish, leading to starvation (Beaugrand et al., 2003). Although past conditions in the North Atlantic provide key insights into interactions between cod and climate, we still have limited knowledge about how the physical and biological changes that result from climate change may affect cod recruitment. Cod production, or growth and recruitment, is strongly linked to the mortality rate during larval and juvenile stages (Sundby et al., 1989). Successful feeding, growth, and survival during this period depends on the combination of physical and biological environmental conditions
such as turbulence, prey abundance, light level, and predator abundance (Leggett & Deblois, 1994) among others. For example, although warmer temperatures increase fish metabolic rates (e.g., K€ oster et al., 2003; O’Connor et al., 2007) enabling higher and faster growth potential, they also require higher prey abundance to sustain higher metabolic rates. In addition, the overlap in time and space between larval fish and phytoplankton and zooplankton (Cushing, 1990; Edwards & Richardson, 2004; Durant et al., 2007) can have a significant impact on fisheries recruitment patterns. Because no single variable determines the survival of larval fish, understanding future dynamics will require considering several mechanisms that all have an effect on fish growth, feeding, and survival (Clark et al., 2003; Beaugrand & Kirby, 2010; Kristiansen et al., 2011). While empirical relationships between fish populations and environmental variables (e.g., bioclimate envelope models or niche models) provide the basis for initial assessments of future responses (e.g., Cheung et al., 2009a), such approaches allow for only limited exploration of the underlying interactions between processes. For example, the behavioral response of larval fish to vertical variations in food abundance can only be considered using models that account for the dynamic mechanistic interactions between larvae and the biophysical conditions of the ocean (DeAngelis & Gross, 1992). In addition, the static parameterization of niche models may limit their robustness for climate change applications as ecosystems evolve into novel states. The continued development and application of ecosystem models that can account for ecology and physiology is essential for improving confidence in projections of climate change impacts on living marine resources (Stock et al., 2011). In this study, we forced a highly mechanistic individual-based cod larvae model (Fiksen & MacKenzie, 2002; Kristiansen et al., 2007, 2009b) with climate change projections from a global Earth System Model (ESM). This combination of models provided a powerful tool for elucidating the mechanisms underlying the response of cod larvae to climate changes. Our modeling framework quantified physical changes, the planktonic ecosystem response, and the larval behavioral response within five habitats, or areas, in the North Atlantic that are historically important spawning and nursery grounds for larval and juvenile stages of Atlantic cod and span the current geographical range of the species. Using our modeling framework, we analyzed how survival through the larval stages may be affected by large-scale climate change from 1950 to 2099 and how changes in survival may differ across habitats in the North Atlantic. We also assessed the relative importance of temperature changes and productivity changes in driving the predicted responses. © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
C L I M A T E C H A N G E E F F E C T S O N F I S H E R I E S 1561 depth. An overview of each of the model components is provided below and the relationships between the components are summarized in Fig. 1. This is followed by a description of the implementation of the modeling framework for each study site and the experiments conducted.
Materials and methods This study combined projected changes in physical and planktonic ecosystem dynamics from an ESM with an IBM for cod larvae to project how changes in climate will affect survival and growth rates of larval cod. A zooplankton model was also used to translate changes in primary production projected by the ESM into changes in the larval cod prey field. The IBM is a two-dimensional (time and depth) vertical model and uses physical and plankton forcing specific to five different regional domains. Each of these domains encompasses large internal spatial variability that we did not capture with this analysis where we focused only on variation with time and
Wind
Chl-a
This study used the Geophysical Fluid Dynamics Laboratory (GFDL) prototype Earth System Model (ESM2.1). This model is built on the physical components of GFDL’s CM2.1 physical climate model (Delworth et al., 2006), with added marine biogeochemistry and land vegetation dynamics. The ocean
Phytoplankton production rates
Earth System Model GFDL ESM v.2.1
Temp.
The earth system model (NOAA GFDL ESM2.1)
Zooplankton model
Light
Prey
Forcing: ESM + zooplankton models
Forcing data for IBM
Ingest
Encounter IBM: Mechanistic feeding module
Capture
Approach
- Temperaturedependent growth Stomach - Temperature and food dependent growth
Update: - weight - stomach
IBM: Bio-energetics module
- Metabolism
Calculate mortality from starvation, invertebrates, and visual predators
IBM: Predation module
Calculate vertical behavior based on ingestion and mortality rates
IBM: Vertical behavior module
Fig. 1 A schematic of the modeling setup consisting of the Earth System Model (ESM), a prey model, and the individual-based model. Temperature, wind (used to calculate turbulence), and current velocity data were used to force the individual-based model directly. The phytoplankton values estimated by the ESM model were used as input to a prey model to estimate abundance of zooplankton at 12 different prey stages. The prey items were used as input to a mechanistic individual-based model that simulated the feeding ecology and bioenergetics of larval cod. Growth was either food-limited or temperature-limited (growth saturated) depending on how much prey the larva consumed in one time-step and the amount of food in the stomach. © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
1562 T . K R I S T I A N S E N et al. component of ESM2.1 is the Modular Ocean Model (Griffies et al., 2005) with 1° 9 1° longitude/latitude horizontal resolution except along the equator where resolution is refined to 1/3°. Fifty vertical levels are resolved, with 10 m resolution over the top 200 m. The TOPAZ (Tracers of Ocean Phytoplankton with Allometric Zooplankton) model provided the biogeochemical dynamics. TOPAZ includes major macro- and micronutrients (N, P, Si, and Fe) and three types of phytoplankton (5 lm ESD large phytoplankton, and diazotrophs). Phytoplankton growth and chlorophyll to carbon ratios are based on the algorithm of Geider et al. (1997). Particle remineralization and sinking is modulated by biogenic and lithogenic minerals (Dunne et al., 2005, 2007) and zooplankton grazing is modeled implicitly (Dunne et al., 2005). In the North Atlantic, historical ocean-ice hindcasts using the MOM4p1-TOPAZ coupling have successfully captured North Atlantic bloom dynamics as observed by SeaWiFS (Henson et al., 2009a) and interannual variations in bloom dynamics observed over the past 50 years in the Continuous Plankton Recorder (CPR) time series (Henson et al., 2009b). However, these simulations did not fully capture plankton regime shifts observed in the CPR record (Henson et al., 2009b). The ESM simulations used herein are identical to earlier work (Rykaczewski & Dunne, 2010; Polovina et al., 2011). The simulations were initialized using data from the World Ocean Atlas 2001 (Conkright et al., 2002) and Carbon Dioxide Information Analysis Center (Key et al., 2004). The model was spun up for 1000 years with a pre-industrial CO2 concentration of 286 ppm (Delworth et al., 2006). Simulations were then made based on the A2 scenario of the IPCC Special Report on Emissions Scenarios (Nakicenovic et al., 2000). This scenario calls for CO2 to increase to 850 ppm in 2100 and is considered a high carbon emission scenario. Focus on a single scenario from a single model prevented a full exploration of response uncertainty across models and scenarios. The intent of this study, however, was to elucidate the potential mechanistic responses of larval cod to characteristic climate-change-driven changes in stratification, temperature, and productivity across the North Atlantic shared by most climate projections (Steinacher et al., 2010; Capotondi et al., 2012). The projected trends will be compared against an ensemble of earth system projections from the fifth IPCC assessment. Six output variables were taken from ESM2.1 to force the zooplankton and larval cod modules: temperature, primary production from large and small phytoplankton, chlorophyll, and wind, and surface light. Profiles of temperature and phytoplankton conditions in the upper 60 m of the water column were used to simulate a dynamically changing environment where the larval fish feeds and grows. The wind conditions at the surface were used to estimate the turbulence level of the water column (MacKenzie & Leggett, 1993), while chlorophyll values were used to determine the attenuation of light with depth (Riley, 1956). For temperature, the difference between the modeled (1950–1999) climatology and the World Ocean Atlas (Locarnini et al., 2010) was used to bias correct ESM projections (Table 1). The 1950–1999 average ocean temperature
estimated by the ESM model for the upper 30 m of the water column was within 1 °C of WOA values for all sites except the Georges Bank/Gulf of Maine (GB/GOM) region, where it is 3.2 °C too warm (Table 1). This warm bias is shared by many coarse resolution climate models (Randall et al., 2007) and is generally attributed to the separation of the Gulf Stream in the models being too far to the north. The magnitude of seasonal changes in temperature was reasonably captured by the model (Table 1), so no bias correction for seasonal temperature changes was applied. The modeled average (2002–2011) primary productivity (g C m2 d1) was within the range of observed 14C-based methods or satellite-based estimates (see ESM vs. Eppley-VGPM; Table 2, Fig. S4), so no bias correction was used.
Zooplankton model A zooplankton model was created to simulate a generic species with six nauplii (I–VI) and six copepodite (I–VI) stages, where the abundance and distribution was calculated based on the primary production estimates obtained from the ESM. The zooplankton model estimated the mesozooplankton biomass (Z, mg C m3) based on estimates of mesozooplankton production and mortality. The mesozooplankton production estimate was derived from small and large phytoplankton primary productivity (PPSP and PPLP, converted to mg C m3 day1) taken from the ESM. Consumption by microzooplankton was assumed to be the dominant loss mechanism for small phytoplankton (Calbet & Landry, 2004) and consumption by mesozooplankton was assumed to be the dominant loss mechanism for large phytoplankton and microzooplankton. Therefore, microzooplankton was not considered a prey for larval fish, but only regarded as a prey resource for mesozooplankton. Assuming a zooplankton gross growth efficiency (gge) of 0.25 (Straile, 1997) then yielded the following estimate of mesozooplankton production (PZ, mg C m3 day1): PZ ¼ gge2 PPSP þ gge PPLP
ð1Þ
PZ was added to the total mesozooplankton biomass in the water column (Z). Mesozooplankton mortality (mZ, mg C m3 day1) was assumed to be temperature and density dependent (Ohman & Hirche, 2001): ðTT0 Þ=10
mZ ¼ Q10
mZ0 Z2
ð2Þ
Here, T was the temperature, Q10 = 2.0 was the factor by which the mortality rate change for a 10 °C rise in water temperature, T0 = 6.4 (Ohman & Hirche, 2001), and mZ0 was set to be consistent with observed mesozooplankton mortality rates (0.2 day1, Ohman et al. (2004)) at typical mesozooplankton concentrations in the North Atlantic. The density-dependent (i.e., quadratic) formulation was assumed to arise from a combination of mesozooplankton predators aggregating over areas of increased mesozooplankton biomass and adult mesozooplanton consuming earlier life stages (Ohman & Hirche, 2001). This formulation produced average zooplankton biomass that was comparable with observed values (see Table 2 for details; Moriarty & O’Brien, 2012). The zooplankton
© 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
C L I M A T E C H A N G E E F F E C T S O N F I S H E R I E S 1563 Table 1 Column two shows the climatological temperature from the World Ocean Atlas 2009 (Locarnini et al., 2010), while column three shows the calculated climatology (1950–1999) from Earth System Model (ESM) for the five areas Georges Bank, West Greenland, Iceland, North Sea, and Lofoten. Column four shows the temperature difference (DT, °C) between the model (ESM) and the observed (WOA). DT was used to bias correct the ESM temperature data that forced the individual-based model (IBM) used in this study. Data for both WOA and ESM were depth averaged for the upper 30 m of the water column over the domain (Fig. 1). Columns 6 to 9 show the seasonal and annual temperature (°C) and salinity (psu) from the World Ocean Atlas (climatology) and the ESM averaged between 1950 and 1999. The World Ocean Atlas shows standard deviations between stations found within the regions and therefore indicates a different standard deviation compared to the standard deviation show for the ESM values, which depicts variability between years.
Location Georges Bank
West Greenland
Iceland
North Sea
Lofoten
Clim (WOA, °C)
Clim (ESM, °C)
10.2
13.39
3.9
8.72
9.88
7.53
4.4
7.87
9.28
6.70
DT (ESMWOA, °C) 3.2
0.5
0.86
0.6
0.83
Period Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND
biomass was divided into zooplankton length intervals of 133 lm ranging from 100 to 1600 lm according to the algorithm described in Daewel et al. (2008). The result is a generic prey item consisting of 12 different size categories from nauplii to copepodite.
The individual-based model The IBM was configured for Atlantic cod and has been thoroughly tested against several observational datasets in prior studies. For example, the IBM was able to reproduce growth and feeding patterns as observed in a macrocosm over the course of the first 47 days posthatching. The IBM adequately © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
WOA temp (°C)
WOA salt (psu)
ESM temp (°C; bias corrected)
ESM salt (psu)
10.24 0.58
32.84 0.09
10.62 3.52
35.67 0.25
5.85 8.25 14.98 11.87 3.90
0.05 0.56 1.72 0.10 0.02
32.88 32.73 32.80 32.97 34.31
0.04 0.11 0.14 0.07 0.06
6.42 13.47 13.34 9.38 1.42
1.16 2.59 1.95 1.75 2.70
35.68 35.52 35.70 35.78 34.15
0.22 0.24 0.25 0.21 0.31
2.89 2.99 5.92 3.80 8.73
0.03 0.03 0.16 0.05 0.1
34.60 34.49 34.09 34.06 35.18
0.02 0.04 0.13 0.05 0.0
0.99 4.34 2.79 0.46 7.10
1.08 1.81 1.87 1.21 1.96
34.39 34.25 33.95 34.02 35.36
0.23 0.26 0.27 0.24 0.08
7.29 8.14 10.89 8.59 9.88
0.02 0.12 0.22 0.02 0.33
35.17 35.18 35.21 35.16 34.60
0.0 0.0 0.01 0.0 0.10
5.61 9.58 7.52 5.69 9.24
0.50 1.40 1.58 0.35 3.29
35.41 35.35 35.32 35.38 27.50
0.03 0.09 0.10 0.06 0.29
5.87 8.72 14.23 10.71 7.53
0.02 0.46 0.86 0.01 0.20
34.64 34.58 34.55 34.61 34.42
0.09 0.12 0.10 0.10 0.09
5.22 11.66 12.14 7.94 8.35
1.21 1.86 1.59 1.90 2.33
27.48 27.27 27.53 27.74 34.59
0.24 0.19 0.25 0.26 0.20
5.7 6.43 10.14 7.87
0.07 0.16 0.79 0.08
34.50 34.54 34.37 34.27
0.04 0.06 0.16 0.09
6.01 10.76 9.61 7.01
0.56 1.90 1.52 0.88
34.73 34.66 34.49 34.49
0.15 0.20 0.19 0.16
reproduced observed growth patterns when the observed environment (prey resources, temperature, and turbulence) was used as input to the IBM (Kristiansen et al., 2007). The IBM has also been tested against a very detailed observational dataset on GB and proved to be able to simulate growth, feeding, behavior, and prey selection comparable with observations for two separate years 1993 and 1994 (Kristiansen et al., 2009b). The IBM consisted of 4 modules: (1) a feeding module, (2) a growth module, (3) a predator module, and (4) a behavioral module. These modules consisted of a number of functions that were estimated sequentially (Fig. 1). All processes were estimated for each time-step as responses may vary with time of day, the depth position in the water column, and with
1564 T . K R I S T I A N S E N et al. Table 2 Table showing average (variable period see below) observed vs. average (1958–1999) modeled zooplankton biomass log10 (mg C m3), and annual and seasonal average modeled phytoplankton productivity (g C m2 d1)
Location
Period
Georges Bank
Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND
West Greenland
Iceland
North Sea
Lofoten
Zooplankton biomass log10 (mg C m3)
Phytoplankton production (g C m2 d1)
ESM
COPEPOD
ESM
Eppley-VGPM
0.85 0.22 (n = 5628)
0.79 0.034 (n = 144)
0.73 0.046 (n = 114)
0.65 0.29 (n = 72)
0.24 1.19 1.25 0.47 0.40
0.06 0.16 0.09 0.18 0.014 (n = 144)
0.27 0.73 1.04 0.61 0.35
0.12 0.26 0.36 0.26 0.029 (n = 85)
0.45 0.078 (n = 43)
0.04 0.75 0.54 0.16 0.40
0.02 0.31 0.08 0.08 0.007 (n = 144)
0.04 0.47 0.42 0.10 0.54
0.04 0.26 0.16 0.11 0.077 (n = 85)
n/a (n = 0)
0.01 0.86 0.69 0.12 0.54
0.00 0.5 0.06 0.08 0.010 (n = 144)
0.05 0.43 0.91 0.10 0.58
0.05 0.27 0.41 0.11 0.028 (n = 95)
0.37 0.26 (n = 165)
0.14 0.91 0.89 0.26 0.43
0.06 0.08 0.05 0.12 0.011 (n = 144)
0.11 0.68 0.63 0.36 0.64
0.13 0.24 0.22 0.15 0.055 (n = 76)
0.01 0.89 0.75 0.08
0.01 0.42 0.14 0.07
0.00 0.67 0.72 0.18
0.00 0.30 0.27 0.23
1.1 0.13 (n = 470)
0.68 0.14 (n = 470)
1.1 0.15 (n = 470)
0.74 0.06 (n = 470)
0.80 0.11 (n = 470)
For both the modeled and observed phytoplankton, the values were averaged for the period 2002–2011. Both observed and modeled average and standard deviation values (n equals number of samples) were calculated based on annual average values for the described time periods over the geographical regions show in Fig. 2. The seasonal averages of observed and modeled phytoplankton productivity were estimated by averaging over the seasons January, February, March (JFM), April, May, June (AMJ), July, August, September (JAS), and October, November, December (OND). All modeled values were taken from the ESM model. Observed zooplankton within each of the 5 regions were extracted from the COPEPOD database (Moriarty & O’Brien, 2012; http:// www.st.nmfs.noaa.gov/plankton/) and averaged over the entire period when observations were available: 1958–2001 for Georges Bank, 1963 for West Greenland and Iceland, North Sea, and 1958–2005 for Lofoten. Eppley-VGPM (Eppley, 1972; Behrenfeld & Falkowski, 1997) primary production estimates were extracted for each region, annually averaged, before average and standard deviation values for 2002–2011 were calculated. The Eppley-VGPM estimates were calculated using observed chlorophyll and temperature data (MODIS), SeaWiFS PAR, and estimates of euphotic zone depth were available from the Oregon State University (www.oregonstate.edu/ocean.productivity/index.php). Also see Fig. S4 for a seasonal comparison between phytoplankton production estimates from ESM and Eppley-VGPM. larval ontogeny. The mechanistic approach relied on a realistic representation of the physical and biological environment, a fundamental understanding of the biology of the prey and the larvae, and the interaction between these components.
Feeding module. The feeding module calculated the number of prey items encountered, approached, captured, and ingested by each larva within one time-step. The encounter
rate was estimated based on the larva’s ability to visually perceive the prey and the prey density (Aksnes & Giske, 1993). The prey differed in length and weight from smaller nauplii to larger copepodites where the smallest prey items were visually hard to detect but easy to catch, and the larger prey were more visible but harder to catch (Kristiansen et al., 2009b). These prey items were visually detected by the larva using a functional response that depended on (1) the eye sensitivity of © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
C L I M A T E C H A N G E E F F E C T S O N F I S H E R I E S 1565 the larva that improved through ontogeny, (2) the image size of the prey, and (3) the light level (Aksnes & Utne, 1997). The light conditions at the depth of the larva were estimated based on how the surface light attenuated with depth, as well as beam attenuation (loss of light) between the predator and the prey. The attenuation coefficient was calculated assuming a background attenuation of 0.1 m1, which increased depending on the chlorophyll level in the surface layer (Riley, 1956). The chlorophyll and light values were taken from the ESM model. The biological characteristics of the prey, such as its contrast to the background and its image area, affected the visibility to the larva (Fiksen & MacKenzie, 2002). The movement of the larva during feeding was modeled as a pause-travel pattern, where the larva searched for prey in the visible half-sphere in front of its snout during pause (O’Brien et al., 1989). If a prey was located within the field of perception, or a prey item was moved into the field of perception due to turbulence, the larva would move to attack position. The probability of attack success (PCA, [0,1]) could be simplified and represented as a linear function: Lprey 3 PCA ¼ max 0:0; minð1:0; 16:7 þ Þ Lfish 2
ð3Þ
Here, Lprey (mm) was the length of the prey, and Lfish (mm) was the length of the larval or juvenile fish (Lfish 2 [3,20] (mm)). As the larva grew and developed, its ability and success in capturing larger prey items increased.
Growth and metabolism module. The feeding module added the total ingested biomass captured during one time-step to the food biomass already stored in the larval gut, where it was used to estimate growth rate. If the larva captured enough prey biomass required to grow at the physiological maximum, the growth was restricted by temperature alone (Folkvord, 2005). If the larva captured less than required, the growth rate was determined by the total energy available in the gut after all other costs (routine and active metabolism, digestion) were accounted for (Kristiansen et al., 2007, 2009b). The temperature- and weight-dependent maximum specific growth rate (Folkvord, 2005) was: SGRðw; TÞ ¼ a b lnðwÞ c lnðwÞ2 þ d lnðwÞ3
ð4Þ
Here, SGR was the specific growth rate (% day1), as a function of temperature (T) and dry mass of the larva (w, mg), parameterized by a = 1.08 + 1.79T, b = 0.074T, c = 0.0965T, d = 0.0112T. The instantaneous growth rate g (day1) was g¼ðlnðSGRdtÞ 7 0.9 100þ1Þ=dt . Routine metabolism R = 57.12e w exp (0.088T) (R, mg day1) of larval cod was parameterized from Finn et al. (2002) as a function of larval weight (dry-weight, w, mg) and temperature (T). Metabolism increased when the larvae were active, and we defined active metabolism as Ractive = 2.5R for larval body length SL ≥5.5 mm and Ractive = 1.4R for larval body length SL
Mechanistic insights into the effects of climate change on larval cod T R O N D K R I S T I A N S E N 1 * , C H A R L E S S T O C K 2 , K E N N E T H F . D R I N K W A T E R 1 and ENRIQUE N. CURCHITSER3 1 Institute of Marine Research, Bergen, 5817 Norway, 2NOAA Geophysical Fluid Dynamics Laboratory, Princeton University, Forrestal Campus, 201 Forrestal Road, Princeton, NJ 08540-6649, USA, 3Department of Environmental Sciences/Institute of Marine and Coastal Sciences, Rutgers University, New Brunswick, NJ 08901 USA
Abstract Understanding the biophysical mechanisms that shape variability in fisheries recruitment is critical for estimating the effects of climate change on fisheries. In this study, we used an Earth System Model (ESM) and a mechanistic individual-based model (IBM) for larval fish to analyze how climate change may impact the growth and survival of larval cod in the North Atlantic. We focused our analysis on five regions that span the current geographical range of cod and are known to contain important spawning populations. Under the SRES A2 (high emissions) scenario, the ESMprojected surface ocean temperatures are expected to increase by >1 °C for 3 of the 5 regions, and stratification is expected to increase at all sites between 1950–1999 and 2050–2099. This enhanced stratification is projected to decrease large (>5 lm ESD) phytoplankton productivity and mesozooplankton biomass at all 5 sites. Higher temperatures are projected to increase larval metabolic costs, which combined with decreased food resources will reduce larval weight, increase the probability of larvae dying from starvation and increase larval exposure to visual and invertebrate predators at most sites. If current concentrations of piscivore and invertebrate predators are maintained, larval survival is projected to decrease at all five sites by 2050–2099. In contrast to past observed responses to climate variability in which warm anomalies led to better recruitment in cold-water stocks, our simulations indicated that reduced prey availability under climate change may cause a reduction in larval survival despite higher temperatures in these regions. In the lower prey environment projected under climate change, higher metabolic costs due to higher temperatures outweigh the advantages of higher growth potential, leading to negative effects on northern cod stocks. Our results provide an important first large-scale assessment of the impacts of climate change on larval cod in the North Atlantic. Keywords: climate change, cod, earth system modeling, fish, Gadus morhua, individual-based modeling, mechanistic modeling, phytoplankton production, recruitment, zooplankton Received 30 August 2012 and accepted 19 November 2013
Introduction Global climate change is expected to alter the physical environment as well as the biological structure and functioning of the world’s ocean ecosystems (Rosenzweig et al., 2007; Doney et al., 2012). Changes in physical features, such as ocean temperature, stratification, and currents, will have considerable impacts on marine ecosystems (Hays et al., 2005; Doney et al., 2012). For example, projected reductions in primary productivity (Steinacher et al., 2010) could have effects throughout the foodweb including the recruitment of fish stocks (Brander, 2007a). Climate change is also expected to *Present address: National Oceanic and Atmospheric Administration, Silver Spring, Maryland, USA Correspondence: Trond Kristiansen, tel. +4797701109, fax +4755238531, e-mail: [email protected]
© 2013 John Wiley & Sons Ltd
affect fish stocks by causing major geographic shifts in species distribution and abundance over the next 50–100 years (Barker & Knorr, 2007; Brander, 2007b; Cheung et al., 2009b). Together these changes could have major direct and indirect impacts on fisheries recruitment, although considerable uncertainty exists with respect to the magnitude of these changes and the mechanisms that underlie them (Brander, 2007a). Observations and climate model simulations suggest that climate change impacts in the North Atlantic will include several physical changes that could affect biological productivity. Observed sea surface temperatures in the North Atlantic over the twentieth century indicate a warming trend everywhere with exception of the far northwestern Atlantic (Deser et al., 2010b). Climate model simulations from the fourth IPCC assessment report project these trends to continue under high CO2 emissions scenarios (Steinacher et al., 2010; Capotondi et al., 2012). Warming and decreased surface salinity 1559
1560 T . K R I S T I A N S E N et al. are projected to increase surface ocean stratification over most of the North Atlantic during the second half of the twenty-first century (Capotondi et al., 2012). Enhanced stratification limits the vertical exchange of cold nutrient-rich water to the surface layer and is projected to decrease primary productivity in most areas of the North Atlantic (Steinacher et al., 2010). The prominent exceptions to this trend occur in extremely highlatitude ecosystems, where increased stratification and warming is projected to increase primary productivity by alleviating strong light limitation and reducing ice cover (Bopp et al., 2001; Steinacher et al., 2010; Doney et al., 2012). These physical changes can have diverse impacts on cod recruitment across the North Atlantic basin. Atlantic cod (Gadus morhua) is a benthic fish species with a pelagic larval stage that is distributed across the North Atlantic in more than 25 different cod populations (Drinkwater, 2005). Past natural climate variability in the North Atlantic basin has caused oscillating cold and warm conditions in cod spawning grounds (Drinkwater, 2006; Sundby & Nakken, 2008). Warm conditions have been found to have a positive effect on survival and recruitment of Atlantic cod through enhanced zooplankton productivity, faster growth rates, and elevated survival (Ottersen et al., 2001). However, the response to ocean warming may depend on the geographical location of the population, the availability of prey resources, and the species involved (Drinkwater, 2005; Cheung et al., 2009a, 2013; Kristiansen et al., 2011). For example, increased ocean temperatures have resulted in stronger recruitment of Barents Sea stocks of cod, haddock, and herring (Ottersen & Loeng, 2000), while increased temperatures in the North Sea have had negative effects on cod recruitment through changes in the prey composition, distribution, and abundance (Beaugrand et al., 2003). In the Barents Sea, warmer temperatures were associated with an increased inflow of water rich in zooplankton, which increased growth and survival through the larval phase (Sundby, 2000). However, in the North Sea, warmer temperatures caused a shift in zooplankton species composition that was lower quality for larval fish, leading to starvation (Beaugrand et al., 2003). Although past conditions in the North Atlantic provide key insights into interactions between cod and climate, we still have limited knowledge about how the physical and biological changes that result from climate change may affect cod recruitment. Cod production, or growth and recruitment, is strongly linked to the mortality rate during larval and juvenile stages (Sundby et al., 1989). Successful feeding, growth, and survival during this period depends on the combination of physical and biological environmental conditions
such as turbulence, prey abundance, light level, and predator abundance (Leggett & Deblois, 1994) among others. For example, although warmer temperatures increase fish metabolic rates (e.g., K€ oster et al., 2003; O’Connor et al., 2007) enabling higher and faster growth potential, they also require higher prey abundance to sustain higher metabolic rates. In addition, the overlap in time and space between larval fish and phytoplankton and zooplankton (Cushing, 1990; Edwards & Richardson, 2004; Durant et al., 2007) can have a significant impact on fisheries recruitment patterns. Because no single variable determines the survival of larval fish, understanding future dynamics will require considering several mechanisms that all have an effect on fish growth, feeding, and survival (Clark et al., 2003; Beaugrand & Kirby, 2010; Kristiansen et al., 2011). While empirical relationships between fish populations and environmental variables (e.g., bioclimate envelope models or niche models) provide the basis for initial assessments of future responses (e.g., Cheung et al., 2009a), such approaches allow for only limited exploration of the underlying interactions between processes. For example, the behavioral response of larval fish to vertical variations in food abundance can only be considered using models that account for the dynamic mechanistic interactions between larvae and the biophysical conditions of the ocean (DeAngelis & Gross, 1992). In addition, the static parameterization of niche models may limit their robustness for climate change applications as ecosystems evolve into novel states. The continued development and application of ecosystem models that can account for ecology and physiology is essential for improving confidence in projections of climate change impacts on living marine resources (Stock et al., 2011). In this study, we forced a highly mechanistic individual-based cod larvae model (Fiksen & MacKenzie, 2002; Kristiansen et al., 2007, 2009b) with climate change projections from a global Earth System Model (ESM). This combination of models provided a powerful tool for elucidating the mechanisms underlying the response of cod larvae to climate changes. Our modeling framework quantified physical changes, the planktonic ecosystem response, and the larval behavioral response within five habitats, or areas, in the North Atlantic that are historically important spawning and nursery grounds for larval and juvenile stages of Atlantic cod and span the current geographical range of the species. Using our modeling framework, we analyzed how survival through the larval stages may be affected by large-scale climate change from 1950 to 2099 and how changes in survival may differ across habitats in the North Atlantic. We also assessed the relative importance of temperature changes and productivity changes in driving the predicted responses. © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
C L I M A T E C H A N G E E F F E C T S O N F I S H E R I E S 1561 depth. An overview of each of the model components is provided below and the relationships between the components are summarized in Fig. 1. This is followed by a description of the implementation of the modeling framework for each study site and the experiments conducted.
Materials and methods This study combined projected changes in physical and planktonic ecosystem dynamics from an ESM with an IBM for cod larvae to project how changes in climate will affect survival and growth rates of larval cod. A zooplankton model was also used to translate changes in primary production projected by the ESM into changes in the larval cod prey field. The IBM is a two-dimensional (time and depth) vertical model and uses physical and plankton forcing specific to five different regional domains. Each of these domains encompasses large internal spatial variability that we did not capture with this analysis where we focused only on variation with time and
Wind
Chl-a
This study used the Geophysical Fluid Dynamics Laboratory (GFDL) prototype Earth System Model (ESM2.1). This model is built on the physical components of GFDL’s CM2.1 physical climate model (Delworth et al., 2006), with added marine biogeochemistry and land vegetation dynamics. The ocean
Phytoplankton production rates
Earth System Model GFDL ESM v.2.1
Temp.
The earth system model (NOAA GFDL ESM2.1)
Zooplankton model
Light
Prey
Forcing: ESM + zooplankton models
Forcing data for IBM
Ingest
Encounter IBM: Mechanistic feeding module
Capture
Approach
- Temperaturedependent growth Stomach - Temperature and food dependent growth
Update: - weight - stomach
IBM: Bio-energetics module
- Metabolism
Calculate mortality from starvation, invertebrates, and visual predators
IBM: Predation module
Calculate vertical behavior based on ingestion and mortality rates
IBM: Vertical behavior module
Fig. 1 A schematic of the modeling setup consisting of the Earth System Model (ESM), a prey model, and the individual-based model. Temperature, wind (used to calculate turbulence), and current velocity data were used to force the individual-based model directly. The phytoplankton values estimated by the ESM model were used as input to a prey model to estimate abundance of zooplankton at 12 different prey stages. The prey items were used as input to a mechanistic individual-based model that simulated the feeding ecology and bioenergetics of larval cod. Growth was either food-limited or temperature-limited (growth saturated) depending on how much prey the larva consumed in one time-step and the amount of food in the stomach. © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
1562 T . K R I S T I A N S E N et al. component of ESM2.1 is the Modular Ocean Model (Griffies et al., 2005) with 1° 9 1° longitude/latitude horizontal resolution except along the equator where resolution is refined to 1/3°. Fifty vertical levels are resolved, with 10 m resolution over the top 200 m. The TOPAZ (Tracers of Ocean Phytoplankton with Allometric Zooplankton) model provided the biogeochemical dynamics. TOPAZ includes major macro- and micronutrients (N, P, Si, and Fe) and three types of phytoplankton (5 lm ESD large phytoplankton, and diazotrophs). Phytoplankton growth and chlorophyll to carbon ratios are based on the algorithm of Geider et al. (1997). Particle remineralization and sinking is modulated by biogenic and lithogenic minerals (Dunne et al., 2005, 2007) and zooplankton grazing is modeled implicitly (Dunne et al., 2005). In the North Atlantic, historical ocean-ice hindcasts using the MOM4p1-TOPAZ coupling have successfully captured North Atlantic bloom dynamics as observed by SeaWiFS (Henson et al., 2009a) and interannual variations in bloom dynamics observed over the past 50 years in the Continuous Plankton Recorder (CPR) time series (Henson et al., 2009b). However, these simulations did not fully capture plankton regime shifts observed in the CPR record (Henson et al., 2009b). The ESM simulations used herein are identical to earlier work (Rykaczewski & Dunne, 2010; Polovina et al., 2011). The simulations were initialized using data from the World Ocean Atlas 2001 (Conkright et al., 2002) and Carbon Dioxide Information Analysis Center (Key et al., 2004). The model was spun up for 1000 years with a pre-industrial CO2 concentration of 286 ppm (Delworth et al., 2006). Simulations were then made based on the A2 scenario of the IPCC Special Report on Emissions Scenarios (Nakicenovic et al., 2000). This scenario calls for CO2 to increase to 850 ppm in 2100 and is considered a high carbon emission scenario. Focus on a single scenario from a single model prevented a full exploration of response uncertainty across models and scenarios. The intent of this study, however, was to elucidate the potential mechanistic responses of larval cod to characteristic climate-change-driven changes in stratification, temperature, and productivity across the North Atlantic shared by most climate projections (Steinacher et al., 2010; Capotondi et al., 2012). The projected trends will be compared against an ensemble of earth system projections from the fifth IPCC assessment. Six output variables were taken from ESM2.1 to force the zooplankton and larval cod modules: temperature, primary production from large and small phytoplankton, chlorophyll, and wind, and surface light. Profiles of temperature and phytoplankton conditions in the upper 60 m of the water column were used to simulate a dynamically changing environment where the larval fish feeds and grows. The wind conditions at the surface were used to estimate the turbulence level of the water column (MacKenzie & Leggett, 1993), while chlorophyll values were used to determine the attenuation of light with depth (Riley, 1956). For temperature, the difference between the modeled (1950–1999) climatology and the World Ocean Atlas (Locarnini et al., 2010) was used to bias correct ESM projections (Table 1). The 1950–1999 average ocean temperature
estimated by the ESM model for the upper 30 m of the water column was within 1 °C of WOA values for all sites except the Georges Bank/Gulf of Maine (GB/GOM) region, where it is 3.2 °C too warm (Table 1). This warm bias is shared by many coarse resolution climate models (Randall et al., 2007) and is generally attributed to the separation of the Gulf Stream in the models being too far to the north. The magnitude of seasonal changes in temperature was reasonably captured by the model (Table 1), so no bias correction for seasonal temperature changes was applied. The modeled average (2002–2011) primary productivity (g C m2 d1) was within the range of observed 14C-based methods or satellite-based estimates (see ESM vs. Eppley-VGPM; Table 2, Fig. S4), so no bias correction was used.
Zooplankton model A zooplankton model was created to simulate a generic species with six nauplii (I–VI) and six copepodite (I–VI) stages, where the abundance and distribution was calculated based on the primary production estimates obtained from the ESM. The zooplankton model estimated the mesozooplankton biomass (Z, mg C m3) based on estimates of mesozooplankton production and mortality. The mesozooplankton production estimate was derived from small and large phytoplankton primary productivity (PPSP and PPLP, converted to mg C m3 day1) taken from the ESM. Consumption by microzooplankton was assumed to be the dominant loss mechanism for small phytoplankton (Calbet & Landry, 2004) and consumption by mesozooplankton was assumed to be the dominant loss mechanism for large phytoplankton and microzooplankton. Therefore, microzooplankton was not considered a prey for larval fish, but only regarded as a prey resource for mesozooplankton. Assuming a zooplankton gross growth efficiency (gge) of 0.25 (Straile, 1997) then yielded the following estimate of mesozooplankton production (PZ, mg C m3 day1): PZ ¼ gge2 PPSP þ gge PPLP
ð1Þ
PZ was added to the total mesozooplankton biomass in the water column (Z). Mesozooplankton mortality (mZ, mg C m3 day1) was assumed to be temperature and density dependent (Ohman & Hirche, 2001): ðTT0 Þ=10
mZ ¼ Q10
mZ0 Z2
ð2Þ
Here, T was the temperature, Q10 = 2.0 was the factor by which the mortality rate change for a 10 °C rise in water temperature, T0 = 6.4 (Ohman & Hirche, 2001), and mZ0 was set to be consistent with observed mesozooplankton mortality rates (0.2 day1, Ohman et al. (2004)) at typical mesozooplankton concentrations in the North Atlantic. The density-dependent (i.e., quadratic) formulation was assumed to arise from a combination of mesozooplankton predators aggregating over areas of increased mesozooplankton biomass and adult mesozooplanton consuming earlier life stages (Ohman & Hirche, 2001). This formulation produced average zooplankton biomass that was comparable with observed values (see Table 2 for details; Moriarty & O’Brien, 2012). The zooplankton
© 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
C L I M A T E C H A N G E E F F E C T S O N F I S H E R I E S 1563 Table 1 Column two shows the climatological temperature from the World Ocean Atlas 2009 (Locarnini et al., 2010), while column three shows the calculated climatology (1950–1999) from Earth System Model (ESM) for the five areas Georges Bank, West Greenland, Iceland, North Sea, and Lofoten. Column four shows the temperature difference (DT, °C) between the model (ESM) and the observed (WOA). DT was used to bias correct the ESM temperature data that forced the individual-based model (IBM) used in this study. Data for both WOA and ESM were depth averaged for the upper 30 m of the water column over the domain (Fig. 1). Columns 6 to 9 show the seasonal and annual temperature (°C) and salinity (psu) from the World Ocean Atlas (climatology) and the ESM averaged between 1950 and 1999. The World Ocean Atlas shows standard deviations between stations found within the regions and therefore indicates a different standard deviation compared to the standard deviation show for the ESM values, which depicts variability between years.
Location Georges Bank
West Greenland
Iceland
North Sea
Lofoten
Clim (WOA, °C)
Clim (ESM, °C)
10.2
13.39
3.9
8.72
9.88
7.53
4.4
7.87
9.28
6.70
DT (ESMWOA, °C) 3.2
0.5
0.86
0.6
0.83
Period Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND
biomass was divided into zooplankton length intervals of 133 lm ranging from 100 to 1600 lm according to the algorithm described in Daewel et al. (2008). The result is a generic prey item consisting of 12 different size categories from nauplii to copepodite.
The individual-based model The IBM was configured for Atlantic cod and has been thoroughly tested against several observational datasets in prior studies. For example, the IBM was able to reproduce growth and feeding patterns as observed in a macrocosm over the course of the first 47 days posthatching. The IBM adequately © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
WOA temp (°C)
WOA salt (psu)
ESM temp (°C; bias corrected)
ESM salt (psu)
10.24 0.58
32.84 0.09
10.62 3.52
35.67 0.25
5.85 8.25 14.98 11.87 3.90
0.05 0.56 1.72 0.10 0.02
32.88 32.73 32.80 32.97 34.31
0.04 0.11 0.14 0.07 0.06
6.42 13.47 13.34 9.38 1.42
1.16 2.59 1.95 1.75 2.70
35.68 35.52 35.70 35.78 34.15
0.22 0.24 0.25 0.21 0.31
2.89 2.99 5.92 3.80 8.73
0.03 0.03 0.16 0.05 0.1
34.60 34.49 34.09 34.06 35.18
0.02 0.04 0.13 0.05 0.0
0.99 4.34 2.79 0.46 7.10
1.08 1.81 1.87 1.21 1.96
34.39 34.25 33.95 34.02 35.36
0.23 0.26 0.27 0.24 0.08
7.29 8.14 10.89 8.59 9.88
0.02 0.12 0.22 0.02 0.33
35.17 35.18 35.21 35.16 34.60
0.0 0.0 0.01 0.0 0.10
5.61 9.58 7.52 5.69 9.24
0.50 1.40 1.58 0.35 3.29
35.41 35.35 35.32 35.38 27.50
0.03 0.09 0.10 0.06 0.29
5.87 8.72 14.23 10.71 7.53
0.02 0.46 0.86 0.01 0.20
34.64 34.58 34.55 34.61 34.42
0.09 0.12 0.10 0.10 0.09
5.22 11.66 12.14 7.94 8.35
1.21 1.86 1.59 1.90 2.33
27.48 27.27 27.53 27.74 34.59
0.24 0.19 0.25 0.26 0.20
5.7 6.43 10.14 7.87
0.07 0.16 0.79 0.08
34.50 34.54 34.37 34.27
0.04 0.06 0.16 0.09
6.01 10.76 9.61 7.01
0.56 1.90 1.52 0.88
34.73 34.66 34.49 34.49
0.15 0.20 0.19 0.16
reproduced observed growth patterns when the observed environment (prey resources, temperature, and turbulence) was used as input to the IBM (Kristiansen et al., 2007). The IBM has also been tested against a very detailed observational dataset on GB and proved to be able to simulate growth, feeding, behavior, and prey selection comparable with observations for two separate years 1993 and 1994 (Kristiansen et al., 2009b). The IBM consisted of 4 modules: (1) a feeding module, (2) a growth module, (3) a predator module, and (4) a behavioral module. These modules consisted of a number of functions that were estimated sequentially (Fig. 1). All processes were estimated for each time-step as responses may vary with time of day, the depth position in the water column, and with
1564 T . K R I S T I A N S E N et al. Table 2 Table showing average (variable period see below) observed vs. average (1958–1999) modeled zooplankton biomass log10 (mg C m3), and annual and seasonal average modeled phytoplankton productivity (g C m2 d1)
Location
Period
Georges Bank
Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND Annual average JFM AMJ JAS OND
West Greenland
Iceland
North Sea
Lofoten
Zooplankton biomass log10 (mg C m3)
Phytoplankton production (g C m2 d1)
ESM
COPEPOD
ESM
Eppley-VGPM
0.85 0.22 (n = 5628)
0.79 0.034 (n = 144)
0.73 0.046 (n = 114)
0.65 0.29 (n = 72)
0.24 1.19 1.25 0.47 0.40
0.06 0.16 0.09 0.18 0.014 (n = 144)
0.27 0.73 1.04 0.61 0.35
0.12 0.26 0.36 0.26 0.029 (n = 85)
0.45 0.078 (n = 43)
0.04 0.75 0.54 0.16 0.40
0.02 0.31 0.08 0.08 0.007 (n = 144)
0.04 0.47 0.42 0.10 0.54
0.04 0.26 0.16 0.11 0.077 (n = 85)
n/a (n = 0)
0.01 0.86 0.69 0.12 0.54
0.00 0.5 0.06 0.08 0.010 (n = 144)
0.05 0.43 0.91 0.10 0.58
0.05 0.27 0.41 0.11 0.028 (n = 95)
0.37 0.26 (n = 165)
0.14 0.91 0.89 0.26 0.43
0.06 0.08 0.05 0.12 0.011 (n = 144)
0.11 0.68 0.63 0.36 0.64
0.13 0.24 0.22 0.15 0.055 (n = 76)
0.01 0.89 0.75 0.08
0.01 0.42 0.14 0.07
0.00 0.67 0.72 0.18
0.00 0.30 0.27 0.23
1.1 0.13 (n = 470)
0.68 0.14 (n = 470)
1.1 0.15 (n = 470)
0.74 0.06 (n = 470)
0.80 0.11 (n = 470)
For both the modeled and observed phytoplankton, the values were averaged for the period 2002–2011. Both observed and modeled average and standard deviation values (n equals number of samples) were calculated based on annual average values for the described time periods over the geographical regions show in Fig. 2. The seasonal averages of observed and modeled phytoplankton productivity were estimated by averaging over the seasons January, February, March (JFM), April, May, June (AMJ), July, August, September (JAS), and October, November, December (OND). All modeled values were taken from the ESM model. Observed zooplankton within each of the 5 regions were extracted from the COPEPOD database (Moriarty & O’Brien, 2012; http:// www.st.nmfs.noaa.gov/plankton/) and averaged over the entire period when observations were available: 1958–2001 for Georges Bank, 1963 for West Greenland and Iceland, North Sea, and 1958–2005 for Lofoten. Eppley-VGPM (Eppley, 1972; Behrenfeld & Falkowski, 1997) primary production estimates were extracted for each region, annually averaged, before average and standard deviation values for 2002–2011 were calculated. The Eppley-VGPM estimates were calculated using observed chlorophyll and temperature data (MODIS), SeaWiFS PAR, and estimates of euphotic zone depth were available from the Oregon State University (www.oregonstate.edu/ocean.productivity/index.php). Also see Fig. S4 for a seasonal comparison between phytoplankton production estimates from ESM and Eppley-VGPM. larval ontogeny. The mechanistic approach relied on a realistic representation of the physical and biological environment, a fundamental understanding of the biology of the prey and the larvae, and the interaction between these components.
Feeding module. The feeding module calculated the number of prey items encountered, approached, captured, and ingested by each larva within one time-step. The encounter
rate was estimated based on the larva’s ability to visually perceive the prey and the prey density (Aksnes & Giske, 1993). The prey differed in length and weight from smaller nauplii to larger copepodites where the smallest prey items were visually hard to detect but easy to catch, and the larger prey were more visible but harder to catch (Kristiansen et al., 2009b). These prey items were visually detected by the larva using a functional response that depended on (1) the eye sensitivity of © 2013 John Wiley & Sons Ltd, Global Change Biology, 20, 1559–1584
C L I M A T E C H A N G E E F F E C T S O N F I S H E R I E S 1565 the larva that improved through ontogeny, (2) the image size of the prey, and (3) the light level (Aksnes & Utne, 1997). The light conditions at the depth of the larva were estimated based on how the surface light attenuated with depth, as well as beam attenuation (loss of light) between the predator and the prey. The attenuation coefficient was calculated assuming a background attenuation of 0.1 m1, which increased depending on the chlorophyll level in the surface layer (Riley, 1956). The chlorophyll and light values were taken from the ESM model. The biological characteristics of the prey, such as its contrast to the background and its image area, affected the visibility to the larva (Fiksen & MacKenzie, 2002). The movement of the larva during feeding was modeled as a pause-travel pattern, where the larva searched for prey in the visible half-sphere in front of its snout during pause (O’Brien et al., 1989). If a prey was located within the field of perception, or a prey item was moved into the field of perception due to turbulence, the larva would move to attack position. The probability of attack success (PCA, [0,1]) could be simplified and represented as a linear function: Lprey 3 PCA ¼ max 0:0; minð1:0; 16:7 þ Þ Lfish 2
ð3Þ
Here, Lprey (mm) was the length of the prey, and Lfish (mm) was the length of the larval or juvenile fish (Lfish 2 [3,20] (mm)). As the larva grew and developed, its ability and success in capturing larger prey items increased.
Growth and metabolism module. The feeding module added the total ingested biomass captured during one time-step to the food biomass already stored in the larval gut, where it was used to estimate growth rate. If the larva captured enough prey biomass required to grow at the physiological maximum, the growth was restricted by temperature alone (Folkvord, 2005). If the larva captured less than required, the growth rate was determined by the total energy available in the gut after all other costs (routine and active metabolism, digestion) were accounted for (Kristiansen et al., 2007, 2009b). The temperature- and weight-dependent maximum specific growth rate (Folkvord, 2005) was: SGRðw; TÞ ¼ a b lnðwÞ c lnðwÞ2 þ d lnðwÞ3
ð4Þ
Here, SGR was the specific growth rate (% day1), as a function of temperature (T) and dry mass of the larva (w, mg), parameterized by a = 1.08 + 1.79T, b = 0.074T, c = 0.0965T, d = 0.0112T. The instantaneous growth rate g (day1) was g¼ðlnðSGRdtÞ 7 0.9 100þ1Þ=dt . Routine metabolism R = 57.12e w exp (0.088T) (R, mg day1) of larval cod was parameterized from Finn et al. (2002) as a function of larval weight (dry-weight, w, mg) and temperature (T). Metabolism increased when the larvae were active, and we defined active metabolism as Ractive = 2.5R for larval body length SL ≥5.5 mm and Ractive = 1.4R for larval body length SL
