Constraining the response of phytoplankton to zooplankton grazing and photo-acclimation in a temperate shelf sea with a 1-D model – towards S2P3 v8.0

Abstract

An established one-dimensional Shelf Sea Physics and Primary Production (S2P3) model has been developed into three different new models: S2P3-NPZ which includes a nutrient–phytoplankton–zooplankton (NPZ) framework, where the grazing rate is no longer fixed but instead varies over time depending on different functions chosen to represent the predator–prey relationship between zooplankton and phytoplankton; S2P3-Photoacclim which includes a representation of the process of photo-acclimation and flexible stoichiometry in phytoplankton; and S2P3 v8.0 which combines the NPZ framework and the variable stoichiometry of phytoplankton at the same time. These model formulations are compared to buoy and conductivity–temperature–depth (CTD) observations, as well as zooplankton biomass and in situ phytoplankton physiological parameters obtained in the central Celtic Sea (CCS). Models were calibrated by comparison to observations of the timing and magnitude of the spring phytoplankton bloom, magnitude of the spring zooplankton bloom, and phytoplankton physiological parameters obtained throughout the water column. A sensitivity study was also performed for each model to understand the effects of individual parameters on model dynamics. Results demonstrate that better agreement with biological observations can be obtained through the addition of representations of photo-acclimation, flexible stoichiometry, and grazing provided these can be adequately constrained.

1 Introduction

Shelf seas are ocean regions where water depth is less than a few hundred metres ($\sim \mathrm{200}$ m) and represent only $\sim \mathrm{10}$ % by area of the global ocean. However, these systems have a disproportionate importance because of their exceptionally high biological productivity , being responsible for 15 % to 30 % of the total oceanic primary production (PP) . Research vessels and remote, autonomous vehicles have been used to study shelf sea regions such as in the Shelf Sea Biogeochemistry (SSB) research programme (https://www.uk-ssb.org/, last access: December 2019), whose aim was to increase the understanding of how physical, chemical, and biological processes interact on UK and European shelf seas, collecting observations throughout 2014 and 2015 in different regions of the UK shelf sea, although these data are not synoptic (i.e. they are not sampled at different locations simultaneously). To complement the available data from research vessels, ocean models have been used to study and understand marine biogeochemistry, including a variety of high-spatial-resolution models to represent the biogeochemistry of shelf seas with high complexity and horizontal spatial resolution .

Different models have been developed to study plankton communities, ranging from very simple ones, e.g. the Lotka–Volterra competition model , to more sophisticated ones, adding more degrees of complexity by including representation of the physical processes of advection and diffusion, or more complexity in ecosystem functions through representation of different groups of organisms and/or size structure. For example, coupled models such as the Nucleus for European Modelling of the Ocean (NEMO) and European Regional Seas Ecosystem Model (ERSEM) , Regional Oceanic Modelling System (ROMS) , and Finite Volume Community Ocean Model (FVCOM) . Although complexity can be useful for describing the interacting behaviour of multiple system components, incomplete understanding of the ecology and key processes of the organisms, and the lack of data for validation can reduce the reliability of predictions. Moreover, simulations of models like NEMO-ERSEM, ROMS, and FVCOM rely on high-performance computing resources, and running multiple sensitivity analyses and experiments is difficult. In contrast, simpler models like the Shelf Sea Physics and Primary Production (S2P3) model have been used to study the dynamics of shelf seas and to simulate seasonal stratification with greater computational efficiency by using a 1-D nutrient–phytoplankton (NP) model to represent physical and biogeochemical processes in the water column. In temperate shelf seas, away from advective sources such as the shelf break or plumes from rivers, horizontal processes can be neglected in comparison to vertical processes; thus, although S2P3 does not consider advective fluxes, it can make a good representation of the dynamics in the water column for temperate shelf seas . Understanding and management of shelf sea ecosystems depend on the level of understanding of factors that influence the communities of resident organisms, and ecosystem models can provide a useful tool to explore these processes. Computationally inexpensive 1-D models like S2P3 can provide useful tools for investigating how different drivers, like changes in the physical environment , can influence shelf sea biogeochemistry. The simplicity of such models combined with the ability to run multiple experiments and long time series can facilitate wider adoption of such models, including by graduate and undergraduate students .

The S2P3 v7.0 model was introduced and developed in the work of , where it was outlined that further development of that model would include resolving phytoplankton physiology. In this study, the S2P3 v7.0 model is developed by allowing variations in light intensity to produce phenotypic adjustments in the phytoplankton cells by changing the chlorophyll content of the phytoplankton and therefore the cellular absorption cross section . This phenotypic change in response to variations in the photon flux density is called photo-acclimation . The main property of photo-acclimation is the reduction of photosynthetic pigment content in response to increased irradiance . Moreover, changes in nutrient availability can further alter cellular chlorophyll and nitrogen quotas , incorporating a combined representation of these two processes . This new version of S2P3 v7.0 is termed S2P3-Photoacclim and relates phytoplankton growth rates to cell quota , describing the light, nutrient, and temperature dependencies of phytoplankton growth rate to varying ratios of N : C : Chl . On the other hand, simpler models using an NPZ or NPZD framework, with the use of nutrients, phytoplankton, zooplankton, and detritus as the main model components e.g., have shown good agreement with observations in terms of chlorophyll and PP, by simulating the timing and magnitude of the spring phytoplankton bloom in different regions of the ocean. Despite their relative simplicity, NPZ models can be a better option to approach an understanding of the physics and biology of an ecosystem, which leads to a further development of the S2P3 v7.0 model where the simplest assumption of a fixed proportion of phytoplankton being grazed and remineralised into the DIN pool (grazing rate) is developed into an NPZ framework (S2P3-NPZ model), using a Holling type 2 or Ivlev grazing functional response of zooplankton grazing on phytoplankton , which shows a saturating response to increasing food.

A combination of photo-acclimation, flexible stoichiometry (S2P3-Photoacclim model), and NPZ framework (S2P3-NPZ model) is then performed to produce a newly developed model called S2P3 v8.0 (Fig. 1). This paper presents a thorough analysis in terms of sensitivity to biological parameter values in each new developed model, resulting in differences of the model structure. A comparison between each model demonstrates how structural differences influence the representation of the spring phytoplankton bloom and annual PP. The aim of this paper is to provide a better understanding of the predator–prey relationship between zooplankton and phytoplankton, the effects of photo-acclimation, and flexible stoichiometry in a simple 1-D model. Model outputs are compared with observations of phytoplankton responses to physical forcing to illustrate the importance of different processes' representation.

2 Study region and model setup

This study is focused on the central Celtic Sea (CCS), a region located in the northwestern (NW) European shelf, which is characterised by its tidally dynamic environment and summer stratification . Daily meteorological data, available from the National Centers for Environmental Prediction (NCEP) reanalysis (http://www.esrl.noaa.gov/psd, last access: December 2019), are used to force the model at the CCS site, located at 49.4${}^{\circ }$ N, 8.6${}^{\circ }$ W. Wind speed (m s${}^{-\mathrm{1}}$), cloud coverage (%), air temperature (${}^{\circ }$C), and relative humidity (%) variables from this dataset are all used to force each model version.

The following description of model setup is applied to each model structure developed in this work. Tidal components consist of the $u$ component (semi-major axis) and the $v$ component (semi-minor axis) for the $M_{\mathrm{2}}$, $S_{\mathrm{2}}$, and $N_{\mathrm{2}}$ tidal constituents. Tidal data are obtained from a fine mesh (12 km resolution) covering the UK shelf. Tidal currents are predicted using the Proudman Oceanographic Laboratory Coastal Ocean Modelling Systems (POLCOMS) 3-D shelf model with an output extracted for the CCS location. Moreover, each model is initialised on 1 January of the first year of simulation with a temperature of 10.10 ${}^{\circ }$C at all depths and water column presumed mixed throughout, including a vertical resolution of 1 m (i.e. 140 vertical levels). Initial values of physical variables are consistent with former studies , whereas initial values of biological variables are based on observations of zooplankton biomass at the CCS location over winter months of 0.02 mmol N m${}^{-\mathrm{3}}$ ; phytoplankton chlorophyll correspond to a typical winter value of 0.2 mg Chl m${}^{-\mathrm{3}}$ for the CCS location, and the DIN initial value is 7 mmol N m${}^{-\mathrm{3}}$. The initialised variables are only set up at the start of each simulation and do not reset in between years.

Figure 1

Structure of the phytoplankton growth formulations: (a) S2P3 v7.0 model, with constant Phyto${}_{\text{N}}$ : Phyto${}_{\text{Chl}}$ ratio (Q); (b) S2P3-NPZ model, including explicit zooplankton with an associated ingestion rate of phytoplankton (I) based on the Ivlev grazing type; (c) S2P3-Photoacclim model, with varying ratios of N : C : Chl, with phytoplankton chlorophyll content regulated by a coefficient of chlorophyll synthesis (${\mathit{\rho }}_{\text{Chl}}$), which reflects the ratio of energy assimilated to energy absorbed , and with phytoplankton carbon regulated through a cost associated with biosynthesis ($u\mathit{\zeta }$) , respiration, and grazing ($G$); (d) the S2P3 v8.0 model, including two different varying quotas for Phyto${}_{\text{N}}$ : Phyto${}_{\text{c}}$ and Phyto${}_{\text{N}}$ : Phyto${}_{\text{Chl}}$.

[Figure omitted. See PDF]

3 Model development

The S2P3 v7.0 model can be divided into two different components: a physical part and a biological part. For this research, the model is an improved version of the original described in to be compiled and executed in a Unix environment , allowing a 1-D representation of physical and biological processes in shelf seas by simulating the seasonal cycle of phytoplankton, water column stratification, and PP at a selected location defined by water depth and tidal current amplitude. The physical part of the model has been greatly described in many other studies and is not described in this section again. This model uses the turbulence closure scheme based on . Likewise, the biological part of the S2P3 v7.0 model is described in (Fig. 1a).

In order to explicitly account for the influence of zooplankton grazing, and hence predator–prey dynamics, the S2P3 v7.0 model is developed into an NPZ framework. This new version of the model (S2P3-NPZ) includes zooplankton as a state variable, contrary to the S2P3 v7.0 model where grazing ($G$) is calculated as a fixed seasonal cycle represented as a sink term in the phytoplankton tendency equation. Addition of an explicit grazer also allows comparison to zooplankton biomass observations. As within S2P3 v7.0, the biological part of the S2P3-NPZ model calculates phytoplankton biomass in chlorophyll currency (Phyto${}_{\text{Chl}}$; mg Chl m${}^{-\mathrm{3}}$). Similar to the S2P3 v7.0 model, the growth of phytoplankton biomass ($\mathit{\mu }$) can be either nutrient limited or temperature limited, with the maximum growth rate of phytoplankton being related to temperature through an Eppley function and being modified by a nutrient quota ($Q$), which corresponds to a varying ratio between phytoplankton nitrogen and phytoplankton chlorophyll biomass ($Q={\text{Phyto}}_{\text{N}}/{\text{Phyto}}_{\text{Chl}}$). Additionally, phytoplankton are modelled in terms of nitrogen (mmol N m${}^{-\mathrm{3}}$) represented by Phyto${}_{\text{N}}$. Zooplankton biomass and external DIN are likewise modelled in terms of nitrogen (mmol N m${}^{-\mathrm{3}}$; Fig. 1b).

The S2P3-Photoacclim model is a new version of the S2P3 v7.0 model, incorporating a more complete representation of phytoplankton physiology . This model allows phytoplankton to acclimate to changes in light and nutrients; therefore, the ratios of N : C : Chl and characteristics of phytoplankton physiology can vary, allowing direct comparison to physiological data. The biological part of the S2P3-Photoacclim model uses three currencies of phytoplankton biomass: carbon (Phyto${}_{\text{C}}$; mg C m${}^{-\mathrm{3}}$), nitrogen (Phyto${}_{\text{N}}$; mg N m${}^{-\mathrm{3}}$), and chlorophyll (Phyto${}_{\text{Chl}}$; mg Chl m${}^{-\mathrm{3}}$). The S2P3-Photoacclim model calculates phytoplankton growth as a function of both nitrogen assimilation and carbon fixation (i.e. variable Chl : N and Chl : C ratios). It is assumed that respiration ($R$) is equal for all cellular components as a function of temperature: $R_{\text{C}}=R_{\text{n}}=R_{\text{Chl}}=R_{\text{ref}}{T}_{\text{function}}$, where $R_{\text{ref}}$ (d${}^{-\mathrm{1}}$) is a degradation rate constant at a reference temperature (Fig. 1c).

The S2P3 v8.0 model describes a combination of zooplankton and physiological acclimation components in order to provide a more realistic representation of the ecosystem dynamics (Fig. 1d). This model can be divided into two different components: a physical part and a biological part. Full details and model equations of the biological part of the model are described below. Details of the variables and parameters are listed in Appendix A.

The biological part of the S2P3 v8.0 model calculates changes in phytoplankton carbon biomass (Phyto${}_{\text{C}}$) over time as

1$$\begin{array}{rl}& {\displaystyle \frac{\partial {\text{Phyto}}_{\text{C}}}{\partial t}}=\\ & {\displaystyle \frac{\partial }{\partial z}}\left(K_Z{\displaystyle \frac{\partial {\text{Phyto}}_{\text{C}}}{\partial z}}\right)+{\text{Phyto}}_{\text{C}}(\mathit{\mu }-R_{\text{C}}{T}_{\text{func}}-u\mathit{\zeta })-I{\displaystyle \frac{Z}{Q_P}},\end{array}$$ where $\mathit{\mu }$ is the growth rate of phytoplankton, $R_{\text{C}}$ is a respiration rate constant of phytoplankton, $T_{\text{func}}$ is a temperature-response function of phytoplankton, $u$ is the phytoplankton carbon-specific nitrate uptake rate, $\mathit{\zeta }$ is the cost of biosynthesis, $I$ is the ingestion rate of zooplankton, and $Q_P$ is the cellular nutrient : carbon quota.

Phytoplankton biomass is also modelled in terms of internal nitrogen, Phyto${}_{\text{N}}$. S2P3 v8.0 calculates Phyto${}_{\text{N}}$ as 2$$\begin{array}{rl}& {\displaystyle \frac{\partial {\text{Phyto}}_{\text{N}}}{\partial t}}=\\ & {\displaystyle \frac{\partial }{\partial z}}\left(K_Z{\displaystyle \frac{\partial {\text{Phyto}}_{\text{N}}}{\partial z}}\right)+u{\text{Phyto}}_{\text{C}}-{\text{Phyto}}_{\text{N}}\left(R_{\text{n}}{T}_{\text{func}}\right)-IZ,\end{array}$$ where $R_{\text{n}}$ is the nitrate remineralisation constant.

The rate of change of phytoplankton biomass in terms of chlorophyll (Phyto${}_{\text{Chl}}$) is described by 3$$\begin{array}{rl}{\displaystyle \frac{\partial {\text{Phyto}}_{\text{Chl}}}{\partial t}}& ={\displaystyle \frac{\partial }{\partial z}}\left(K_Z{\displaystyle \frac{\partial {\text{Phyto}}_{\text{Chl}}}{\partial z}}\right)+u{\mathit{\rho }}_{\text{Chl}}{\text{Phyto}}_{\text{C}}\\ & -{\text{Phyto}}_{\text{Chl}}\left(R_{\text{Chl}}{T}_{\text{func}}\right)-I{\displaystyle \frac{Z}{Q}},\end{array}$$ where ${\mathit{\rho }}_{\text{Chl}}$ is a chlorophyll synthesis regulation term, $R_{\text{Chl}}$ is the chlorophyll degradation rate constant, and $Q$ is the cellular nutrient quota (N : Chl).

The change in time of external DIN is calculated as 4$$\begin{array}{rl}{\displaystyle \frac{\partial N}{\partial t}}& ={\displaystyle \frac{\partial }{\partial z}}\left(K_Z{\displaystyle \frac{\partial N}{\partial z}}\right)+{\mathit{\gamma }}_{\mathrm{1}}IZ+{\mathit{\gamma }}_{\mathrm{2}}mZ\\ & +{\text{Phyto}}_{\text{N}}\left(R_{\text{n}}{T}_{\text{func}}\right)-u{\text{Phyto}}_{\text{C}},\end{array}$$ where ${\mathit{\gamma }}_{\mathrm{1}}$ is the grazing inefficiency or “messy feeding” that returns a fraction of grazed material back into the DIN pool, ${\mathit{\gamma }}_{\mathrm{2}}$ is the fraction of dead zooplankton that goes into the sediments, and $m$ is the loss rate of zooplankton due to predation and physiological death.

Zooplankton grazing on phytoplankton response depends on a Holling type 2 or Ivlev grazing , with the ingestion rate of zooplankton ($I$) described as 5$$I=R_{\text{max}}\left(\mathrm{1}-e^{\left(-\mathit{\lambda }{\text{Phyto}}_{\text{N}}\right)}\right),$$ where $R_{\text{max}}$ is the zooplankton maximal grazing rate (d${}^{-\mathrm{1}}$) and $\mathit{\lambda }$ is the rate at which saturation is achieved with increasing food levels (mmol N m${}^{-\mathrm{3}}$)${}^{-\mathrm{1}}$.

Zooplankton biomass is therefore modelled as 6$${\displaystyle \frac{\partial Z}{\partial t}}={\displaystyle \frac{\partial }{\partial z}}\left(K_Z{\displaystyle \frac{\partial Z}{\partial z}}\right)+(\mathrm{1}-{\mathit{\gamma }}_{\mathrm{1}})IZ-mZ.$$

The carbon-specific, light-saturated photosynthetic rate depends on the internal nitrogen of phytoplankton described as a factor $f$ based on the work of : 7$$P_{\text{max}}=P_{\text{max}}^{\text{C}}f,$$ where $P_{\text{max}}^{\text{C}}$ is the maximum value of the carbon-specific rate of photosynthesis and the factor $f$ is defined as 8$$f={\displaystyle \frac{Q-Q_{\text{min}}}{Q_{\text{max}}-Q_{\text{min}}}}.$$ Nitrogen assimilation is calculated as a saturating function of external nutrient content, the internal nutrient quota and the maximum nitrogen assimilation rate ($u_{\text{max}}$) according to 9$$u=u_{\text{max}}\left({\displaystyle \frac{\mathrm{1}-f}{\mathrm{1.015}-f}}\right)\left({\displaystyle \frac{N}{k_{\text{n}}+N}}\right),$$ where $k_{\text{n}}$ is the half-saturation constant for nitrate uptake.

The carbon-specific photosynthesis is a saturating function of irradiance and it is calculated as 10$$\mathit{\mu }=P_{\text{max}}\left(\mathrm{1}-e^{-\left(\frac{{\mathit{\alpha }}^{\text{Chl}}{I}_{\text{PAR}}\mathit{\theta }}{P_{\text{max}}}\right)}\right),$$ where ${\mathit{\alpha }}^{\text{Chl}}$ is the chlorophyll-specific initial slope of the photosynthesis–light curve, $\mathit{\theta }$ is the cellular chlorophyll : phytoplankton carbon ratio, and $I_{\text{PAR}}$ is the photosynthetically available radiation.

Finally, chlorophyll a synthesis depends on the rates of photosynthesis and light absorption: 11$${\mathit{\rho }}_{\text{Chl}}={\mathit{\theta }}_{\text{max}}^{\text{N}}\left({\displaystyle \frac{\mathit{\mu }}{{\mathit{\alpha }}^{\text{Chl}}{I}_{\text{PAR}}\mathit{\theta }}}\right),$$ where ${\mathit{\theta }}_{\text{max}}^{\text{N}}$ is the maximum value of the cellular chlorophyll : phytoplankton nitrogen ratio.

4 Validation of the models: observations

To calibrate or tune each model, they were adjusted on a trial-and-error basis until disagreement with the in situ observations was minimised, allowing investigation of the sensitivity of each model to changes in the parameters listed in Table 1.

4.1 UK SSB programme

Time series of surface chlorophyll a concentrations (mg Chl m${}^{-\mathrm{3}}$) from long-term mooring deployments including the Carbon and Nutrient Dynamics and Fluxes over Shelf Systems (CaNDyFloSS) Smartbuoy were collected at the CCS location (49.4${}^{\circ }$ N, 8.6${}^{\circ }$ W; depth of 145.8 m), gathering data for 5 min every 30 min during the years 2014 and 2015 as part of the research cruise expeditions DY029 and DY033. The phytoplankton community fluorescence from the water samples was calculated as a proxy for chlorophyll a and calibrated taking into account daytime fluorescence quenching which results in a reduction of fluorescence per unit chlorophyll. For this study, daytime data were removed.

Conductivity–temperature–depth (CTD) casts were performed in different locations on the NW European shelf from the CCS location to the shelf break, with discrete samples of temperature, DIN, and chlorophyll a collected using Niskin bottles as part of the research cruise expeditions DY029 and DY033. At the CCS location, the CTD samples were collected from pre-dawn to midday with a 1 m vertical resolution over the whole water column (140 m depth) during the year 2015. CTD casts for the CCS location were chosen to validate the model during spring and summer. Relevant information about dates and positions from other CTD casts taken during spring and summer of the year 2015 is listed in Table A1. CTD casts chosen to compare the each model are marked in bold in Table B1.

4.2 Zooplankton biomass

Zooplankton biomass samples were collected at the CCS location (49${}^{\circ }$25 N, 8${}^{\circ }$35 W; $\sim \mathrm{150}$ m water depth) during four periods: 5–12 August, 10–29 November 2014, 3–28 April, and 13–31 July 2015 for the cruises DY026, DY018, DY029, and DY033, respectively . Zooplankton were fractionated into microplankton, small mesozooplankton, and large mesozooplankton by using different mesh sizes. For zooplankton biomass samples, net rings of 57 cm diameter were used and fitted with two different mesh sizes of 63 and 200 $\mathrm{\mu }$m. The nets had a closing mechanism when deployed, sampling zooplankton biomass during daytime and nighttime at different depths: above and below the thermocline, and, when present, across the deep chlorophyll maximum (DCM; determined based on fluorescence measurements). The thermocline and DCM were determined from CTD casts immediately prior to the net deployments. The 63 and 200 $\mathrm{\mu }$m mesh nets were hauled at 0.2 and 0.5 m s${}^{-\mathrm{1}}$, respectively.

For the S2P3-NPZ and the S2P3 v8.0 models, only mesozooplankton biomass were considered, with a community composition that included amphipods, appendicularian, Chaetognatha, copepods, Euphausiacea, Polychaeta, and others (e.g. cladocerans, dinoflagellates, echinoderm, eggs, foraminifera, Gymnosomata, unidentified larvae, nauplii, ostracods, and radiolarian, all of which contributed $<\mathrm{3}$ % in all samples). A FlowCam (Fluid Imaging Technologies Inc.) and a ZooScan were used to scan zooplankton individuals with images processed using ZooProcess 7.19 and Plankton Identifier 1.3.4 software (Gorsky et al., 2010). From these images, biovolume spectra were calculated and converted into image-derived dry weight (DW). The total 246 net hauls collected for biomass samples provided 44 vertical depth profiles, integrating zooplankton biomass typically between 0 and 120 m at the CCS location. The complete dataset can be obtained from the British Oceanographic Data Centre (BODC; http://www.bodc.ac.uk/data, last access: December 2019.) as reported in .

4.3 Physiological observations

Water samples were collected in the Celtic Sea from cruises JR98 and CD173 using 10 L Niskin bottles for four of the sites (CS1, CS2, CS3, and IS1) (Fig. 2) during 24 h periods on 31 July, 29 July, 5 August, and 2 August, respectively. Multiple samples were collected at the surface and in deeper layers (in the surface mixed layer (SML) and the DCM) to obtain different phytoplankton populations throughout the photoperiod. The JR98 cruise was undertaken from 24 July to 14 August of 2003. Stations, from the Irish Sea to the Celtic Sea shelf break, ranged in characteristics from very strong, narrow thermoclines in the southern Celtic Sea (CS1), to the weak, deep surface layer associated with internal wave mixing at the shelf edge (CS2). On the other hand, the CD173 cruise was undertaken from 15 July to 6 August of 2005, from the stratified region of the Celtic Sea shelf (stations D2, CS1, CS3, U2) and shelf break (stations CS2, N1). For this work, observations from both cruises were used, considering only the stations from the seasonally stratified sites (B2, CS1, CS3, D2, JB1, OB, P1, U2, and ctd16) and excluding stations CS2 and N1 which are close to the shelf edge where advective fluxes are more relevant than in the stations closer to the CCS location and none of the models used here consider advective fluxes.

Photosynthesis vs. irradiance ($P$ vs. $E$) experiments were conducted in short-term incubations (2–4 h) using a photosynthetron . From these $P$ vs. $E$ experiments, chlorophyll a normalised PP was derived from ${}^{\mathrm{14}}$C uptake to obtain the chlorophyll a specific maximum light-saturated photosynthesis rate $P_{\text{max}}^{\text{Chl}}$ (mg C (mg Chl $a$)${}^{-\mathrm{1}}$ h${}^{-\mathrm{1}}$) and the maximum light utilisation coefficient, ${\mathit{\alpha }}^{\text{Chl}}$ (mg C (mg Chl $a$)${}^{-\mathrm{1}}$ h${}^{-\mathrm{1}}$ ($\mathrm{\mu }$$E$ m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}$)${}^{-\mathrm{1}})$ . Values of ${\mathit{\alpha }}^{\text{Chl}}$ and the light saturation parameter, $E_k$ ($\mathrm{\mu }$$E$ m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}$) (given by $E_k=P_{\text{max}}^{\text{Chl}}/{\mathit{\alpha }}^{\text{Chl}}$) were spectrally corrected to the in situ irradiance at the sample depth according to the phytoplankton light absorption . The maximum light utilisation coefficient (${\mathit{\alpha }}^{\text{Chl}}$) was constrained for the S2P3-Photoacclim and the S2P3 v8.0 models by finding the mean from all observations: ${\mathit{\alpha }}^{\text{Chl}}=\mathrm{9.16}\times {\mathrm{10}}^{-\mathrm{6}}$ mg C (mg Chl $a$)${}^{-\mathrm{1}}$ h${}^{-\mathrm{1}}$ ($\mathrm{\mu }$$E$ m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}{)}^{-\mathrm{1}}$ (ranges of this parameter accounted from $\mathrm{3.58}\times {\mathrm{10}}^{-\mathrm{6}}$ to $\mathrm{3.59}\times {\mathrm{10}}^{-\mathrm{5}}$ mg C (mg Chl $a$)${}^{-\mathrm{1}}$ h${}^{-\mathrm{1}}$ ($\mathrm{\mu }$$E$ m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}{)}^{-\mathrm{1}}$; SD $=\mathrm{4.38}\times {\mathrm{10}}^{-\mathrm{6}}$), and therefore the values of $P_{\text{max}}^{\text{Chl}}$ and $E_k$ were used as variables for comparison with equivalent modelled values.

Figure 2

Map of study area and stations for the JR98 and CD173 cruises, including the CCS location (in red colour). The image was created with MATLAB using the repository data for gridded bathymetry provided by General Bathymetric Chart of the Oceans (GEBCO). Bathymetric data are only considered for the shelf sea region (0 to 300 m depth) with open ocean depth neglected (deeper than 300 m). Continents considered in black colour (over 0 m elevation).

[Figure omitted. See PDF]

Table 1

List of parameter values, including units and definitions for the calibrated S2P3-NPZ, S2P3-Photoacclim, and S2P3 v8.0 models.

Parameters (units)	Definition	S2P3-NPZ	S2P3-Photoacclim	S2P3 v8.0
${\mathit{\gamma }}_{\mathrm{1}}$ (dimensionless)	Grazing inefficiency or “messy feeding” (0.0–1.0); returns a fraction of grazed material backinto the DIN pool	0.2		0.1
${\mathit{\gamma }}_{\mathrm{2}}$ (dimensionless)	Fraction of dead zooplankton (0.0–1.0) that goes into the sediments or higher trophic levels	0.5		0.4
$\mathit{\lambda }$ (mmol N m${}^{-\mathrm{3}}{)}^{-\mathrm{1}}$	Rate at which saturation is achieved with increasing food levels	0.053		0.014
$R_{\text{max}}$ (d${}^{-\mathrm{1}}$)	Zooplankton maximal grazing rate	2.5		3.5
$m$ (d${}^{-\mathrm{1}}$)	Loss rate of zooplankton due to predation andphysiological death	0.05		0.02
$P_{\text{max}}^{\text{C}}$ (d${}^{-\mathrm{1}}$)	Maximum value of the carbon-specific rate ofphotosynthesis		2.0	3.5
$Q_{\text{max}}$ (mmol N (mg C)${}^{-\mathrm{1}}$)	Maximum value of the cellular nutrient quota		0.028	0.032
${\mathit{\theta }}_{\text{max}}^{\text{N}}$ (mg Chl (mmol N)${}^{-\mathrm{1}}$)	Maximum value of the chlorophyll : phytoplankton nitrogen ratio		4.2	2.1
$R_{\text{C}}=R_{\text{n}}=R_{\text{Chl}}$ (d${}^{-\mathrm{1}}$)	Respiration rates		0.02	0.02

4.4 Model calibration

For this study, approximately 2 years of phytoplankton chlorophyll data were available for the CCS location, while in the case of the zooplankton biomass observations, these were collected only during certain days per year allowing only a discrete representation of the seasonal cycle of zooplankton. Finally, profiles of physiological data were only collected during summertime of the years 2003 and 2005. The calibrated version of the S2P3-NPZ model shows differences in the timing of the spring phytoplankton bloom for the year 2015 in comparison to observations of surface chlorophyll (Fig. 3a), with a later bloom from the S2P3-NPZ model, reaching a peak bloom about a month later. Additionally, the magnitude of the spring phytoplankton bloom is also higher in the model in comparison to observations. Phytoplankton are able to escape grazing control in April and early May, with the spring zooplankton bloom occurring about a month later (Fig. 3b). Similar differences can be observed between the calibrated S2P3-Photoacclim model, demonstrating that constraining the timing of the spring phytoplankton bloom is a complex process. Tuning model parameters to generate earlier blooms modifies the magnitude of the spring phytoplankton bloom by increasing it to unrealistic levels. Figure 3 shows that between the three models the best agreement found in comparison to buoy observations correspond to the S2P3 v8.0 model, although the timing and magnitude of the spring phytoplankton bloom show some remaining small differences, with higher concentrations of surface chlorophyll a during spring ($\sim \mathrm{10}$ mg Chl m${}^{-\mathrm{3}}$). Moreover, the timing of the spring phytoplankton bloom during the year 2014 matches the observations but a delayed bloom is shown during the year 2015. On the other hand, zooplankton biomass is higher than in the S2P3-NPZ and the predator–prey relationship is well represented, with the spring zooplankton bloom happening approximately a month later than the spring phytoplankton bloom during the year 2014, but this difference is lower during the year 2015, with the start of the zooplankton bloom happening about half a month after the start of the spring phytoplankton bloom. Remaining mismatches between the calibrated models and observations may be driven by water column processes including advection and diffusion that were not considered in these 1-D models but affects the real water column where the observations were taken. Furthermore, photo-acclimation and grazing depend on changes in temperature and other ecosystem processes which are not explicitly represented, potentially explaining remaining mismatches between the observations and the calibrated S2P3 v8.0 model.

To provide a quantitative index of bloom timing, we consider a threshold criterion . In this study, the spring phytoplankton bloom is defined as when surface chlorophyll reaches more than 1.5 mg Chl m${}^{-\mathrm{3}}$ (see Table B2). For the S2P3-Photoacclim model in comparison to CTD cast during spring (Fig. 4b), the model has not yet reached the spring phytoplankton bloom, whereas observations indicate the spring phytoplankton bloom has already started; the S2P3-Photoacclim model shows low chlorophyll a concentrations at the surface, with DIN concentrations not being depleted at this stage (Fig. 4c). Similar results can be seen for the S2P3-NPZ, although Table B2 shows that the spring phytoplankton bloom has already started in this model, but it is still later than in the observations. The S2P3 v8.0 model, on the contrary, shows a better agreement in terms of the timing of the spring phytoplankton bloom with the CTD observations. During summer months, the models are able to reproduce the subsurface mixed layer observed in the CTD profile (Fig. 4e) with a similar magnitude but shallower by approximately 20 m. On the other hand, the physical structure of the model shows good agreement with the observations during spring (Fig. 4a), but there are differences during summer (Fig. 4d). The lack of a marked mixed layer depth in all the models is likely related to short-term meteorology forcing: during 24 July 2015, the air temperature was high ($\sim \mathrm{20}$ ${}^{\circ }$C) and wind speed was low ($\sim \mathrm{5}$ m s${}^{-\mathrm{1}}$). However, the thermocline shows a sharp development in the CTD observations that cannot be constrained better in the models by parameterising values of the turbulent closure scheme, light attenuation in the water column, and mixing control parameters (data not shown).

The S2P3-Photoacclim and the S2P3 v8.0 models were further compared with phytoplankton physiological variability observations (Fig. 5). Figure 5 shows that near the sea surface, where light levels are high, photo-acclimation of phytoplankton, values of $P_{\text{max}}^{\text{Chl}}$ and $E_k$ are higher than in deeper layers of the water column. The observed $P_{\text{max}}^{\text{Chl}}$ ranges approximately 0.5–2.5 $\times$ ${\mathrm{10}}^{-\mathrm{3}}$ (mg C (mg Chl $a$)${}^{-\mathrm{1}}$ s${}^{-\mathrm{1}}$) in the surface waters (first 5 m), while this range is smaller in deeper layers (0–1.0 $\times {\mathrm{10}}^{-\mathrm{3}}$ mg C (mg Chl $a$)${}^{-\mathrm{1}}$ s${}^{-\mathrm{1}}$), e.g. at 40 m depth. Similar variability can be observed with $E_k$, having lower values in deeper layers of the water column, but the largest variability occurs in the surface layer ($\sim$ 100–250 $\mathrm{\mu }$$E$ m${}^{-\mathrm{2}}$ s${}^{-\mathrm{1}}$). The S2P3-Photoacclim and S2P3 v8.0 models show good agreement with the observations, with plausible values of $P_{\text{max}}^{\text{Chl}}$ and $E_k$ found through the water column particularly in terms of the magnitude of the vertical gradients. However, the version of the S2P3 v8.0 model which has the best parameterisations when tuned to fit the timing of the spring bloom has remaining discrepancies when compared to the physiological observations, showing an overestimation of $P_{\text{max}}^{\text{Chl}}$ and $E_k$ at depth. Calibration of these models against physiological data is relatively novel as such comparisons remain rare. Greater complexity allowed the S2P3 v8.0 model to resolve a more diverse range of biogeochemical dynamics, explicitly accounting for zooplankton biomass and for the dynamics of internal quotas of phytoplanktonic cells, with phytoplankton biomass being in carbon, nitrogen, and chlorophyll currencies, allowing the decoupling of nutrient uptake from carbon fixation . Including additional parameters in models can add more unconstrained degrees of freedom but also allows for more parameter combinations and therefore more flexibility to constrain S2P3 v8.0 in order to reproduce observations. The new added parameters and variables in the S2P3 v8.0 model, were carefully chosen to allow the model to be constrained against additional data, specifically zooplankton abundances and photosynthetic physiological measurements; therefore, a higher complexity allowed better representations of the temporal dynamics. Additionally, despite having more sophisticated formulations of the ecosystem, the S2P3 v8.0 model continues to be a 1-D model, allowing multiple experiments to be run at the same time with relatively low computational cost.

Figure 3

(a) SSB observations of surface chlorophyll a (black line), along with the modelled surface chlorophyll a for the S2P3-NPZ (red line), S2P3-Photoacclim (cyan line), and S2P3 v8.0 (green line) calibrated models. (b) Observations of zooplankton biomass presented as discrete points for day nets (red dots) and night nets (blue dots) taken during the cruises DY026, DY018, DY029, and DY033; modelled zooplankton biomass from the S2P3-NPZ (red line) and S2P3 v8.0 (green line) calibrated models.

[Figure omitted. See PDF]

Figure 4

CTD observations from the SSB programme (black lines) including data for springtime (20 April 2015) (a) temperature, (b) chlorophyll a, and (c) DIN (black dots); for summertime (24 July 2015) (d) temperature, (e) chlorophyll a, and (f) DIN (black dots) along the S2P3-NPZ (red line), S2P3-Photoacclim (cyan line), and S2P3 v8.0 (green line) calibrated models.

[Figure omitted. See PDF]

Figure 5

Observations from the cruises CD173 and JR98 for (a) chlorophyll $a$ specific maximum light-saturated photosynthesis rate $\left(P_{\text{max}}^{\text{Chl}}\right)$ in different locations in the Celtic Sea and for the calibrated S2P3-Photoacclim model (cyan lines) and the S2P3 v8.0 model (green lines); (b) observations of the light saturation parameter ($E_k$) for different stations across the Celtic Sea and for the calibrated S2P3-Photoacclim model (cyan lines) and the S2P3 v8.0 model (green lines). The data from both models were plotted for the same days that the observations were collected.

[Figure omitted. See PDF]

The behaviour of the S2P3 v8.0 model calibrated for the CCS location is displayed in Figs. 6 and 7. Figure 6 shows contour plots from daily profiles of the S2P3 v8.0 model for temperature (Fig. 6a), phytoplankton chlorophyll $a$ (Fig. 6b), zooplankton biomass (Fig. 6c), phytoplankton chlorophyll : phytoplankton carbon ratio (Fig. 6d), and DIN (Fig. 6e) for the years 2014 and 2015 which correspond to the observation period of the SSB programme. Figure 6 allows a more detailed overview of the model dynamics. Water column temperature increases from April of each year, reaching a maximum value at the surface during summer months. At the same time, the spring phytoplankton bloom can be observed during April, reaching $\sim \mathrm{10}$ mg Chl m${}^{-\mathrm{3}}$. Additionally, the spring zooplankton bloom can be observed approximately a month after the spring phytoplankton bloom is developed, with zooplankton being able to grow during summer months and decreasing until a minimum value during winter. These spring blooms also mark the start of DIN depletion at the surface, a state that lasts until the end of summer. Finally, the phytoplankton chlorophyll : phytoplankton carbon ratio shows the highest values during winter months when irradiance levels are low, with the Chl : N ratio decreasing until the end of summer due to the lower concentrations of chlorophyll in the cell to avoid internal damage due to high irradiance during this period, highlighting the photo-acclimation of phytoplankton and flexible stoichiometry of the S2P3 v8.0 model.

Figure 7 provides a general overview of the dynamics of the calibrated S2P3 v8.0 model, representing the interannual variability of each variable using the median (black lines), and lower and upper quantiles (red lines; 95 % of the data distribution) during 1965 to 2015. Figure 7a shows the start of thermal stratification during early spring with observable interannual variability in the extent of stratification. Once the water column is stratified, the spring phytoplankton bloom can be observed (Fig. 7b), followed by the start of the spring zooplankton bloom (Fig. 7c). As phytoplankton grow, DIN concentrations at the surface start to deplete reaching a minimum value during spring and summer months, and increasing when thermal stratification breaks down during winter months (Fig. 7d). Finally, net primary production (NPP) time series show seasonal and interannual phytoplankton dynamics (Fig. 7e). All variables of the S2P3 v8.0 model shown in Fig. 7 provide the spectrum of interannual variability of 95 % of the data using the upper and lower quantile values (red lines), demonstrating that the interannual variability of thermal stratification provides variability for the timing and magnitude of the spring phytoplankton bloom as well as summer growth.

Figure 6

Contoured daily vertical profiles for the start of 2014 to the end of 2015 for the calibrated version of the S2P3 v8.0 model including (a) temperature (${}^{\circ }$C), (b) phytoplankton chlorophyll $a$ (mg Chl m${}^{-\mathrm{3}}$), (c) zooplankton biomass (mmol N m${}^{-\mathrm{3}}$), (d) phytoplankton chlorophyll : phytoplankton carbon ratio (mg Chl (mg C)${}^{-\mathrm{1}}$), and (e) DIN (mmol N m${}^{-\mathrm{3}}$).

[Figure omitted. See PDF]

Figure 7

Annual representation of the median (black lines) calculated from 1965 to 2015 for the S2P3 v8.0 model calibrated for the CCS location and forced with all meteorological components (i.e. wind speed, cloud coverage, air temperature, and relative humidity). Red lines represent the annual lower and upper quantiles for each variable of the model (95 % distribution of the data) over 1965–2015. (a) Surface temperature minus bottom temperature (${}^{\circ }$C), (b) surface chlorophyll $a$ (mg Chl m${}^{-\mathrm{3}}$), (c) surface zooplankton biomass (mmol N m${}^{-\mathrm{3}}$), (d) surface DIN (mmol N m${}^{-\mathrm{3}}$), and (e) net primary production (NPP) (mg C m${}^{-\mathrm{2}}$ d${}^{-\mathrm{1}}$).

[Figure omitted. See PDF]

Table 2

List of the most sensitive experiments run for the S2P3 v8.0 model calibrated for the CCS location, including the year of observations, timing and magnitude of the spring phytoplankton bloom, and total annual zooplankton biomass values.

		Timing spring	Magnitude spring	Total annual
		phytoplankton	phytoplankton	zooplankton
Experiments	Years	bloom (date)	bloom (mg Chl m${}^{-\mathrm{3}}$)	biomass (g DW m${}^{-\mathrm{2}}$)
Calibrated S2P3 v8.0	2014	3 April	137.7	6110
	2015	9 April	172.6	6094
${\mathit{\theta }}_{\text{max}}^{\text{N}}\downarrow$	2014	19 May	103.5	3055
	2015	21 May	103.1	3514
${\mathit{\theta }}_{\text{max}}^{\text{N}}\uparrow$	2014	25 March	110.2	6733
	2015	27 March	211.5	6797
$R_{\text{max}}\downarrow$	2014	29 March	292.5	6840
	2015	3 April	359.2	6812
$R_{\text{max}}\uparrow$	2014	14 April	61.4	5283
	2015	25 April	159.5	4992
$m\downarrow$	2014	3 May	238	8362
	2015	14 May	203.5	8369
$m\uparrow$	2014	31 March	104.4	4311
	2015	6 April	134.5	4368

5 Sensitivity studies

An analysis was performed to assess the sensitivity of the model to selected parameter values for the S2P3 v8.0 model. Each parameter listed in Table 1 was varied in turn in order to understand how sensitive the model is to those changes and the effect that they have on the modelled ecosystem dynamics at the CCS location. In each case, parameters were varied from the best calibrated value by $+\mathrm{50}$ % and $-\mathrm{50}$ %. Sensitivity studies are important tools to improve the accuracy of shelf sea models , but developing these analyses has to be done carefully in order to identify which processes are responsible for the observed model behaviour . A direct comparison was calculated in terms of the S2P3 v8.0 model attributes considering the timing and magnitude of the spring phytoplankton bloom, and the total annual zooplankton biomass for the calibrated version of the model and each experiment, providing better insights on the effects that each parameter produces in the behaviour of each model. Table 2 only shows the experiments involving the zooplankton maximal grazing rate ($R_{\text{max}}$), mortality of zooplankton (m), and the maximum value of the Chl : N ratio (${\mathit{\theta }}_{\text{max}}^{\text{N}}$), as they were demonstrated to be the most significant parameters in terms of the sensitivity of the model according to the attributes calculated, with the rest of the experiments omitted for this discussion. The S2P3 v8.0 model is strongly influenced by NPZ parameters, with the zooplankton maximal grazing rate ($R_{\text{max}}$) and zooplankton mortality rate ($m$) having the largest effect in the magnitude of the spring phytoplankton bloom (Figs. C1a, C4a) and in the total annual zooplankton biomass (Figs. C1b, C4b). The effect of $R_{\text{max}}$ implies that lower values in the maximum ingestion rate of phytoplankton can produce earlier and larger spring phytoplankton blooms compared to the calibrated S2P3 v8.0 model. On the other hand, a higher value of $R_{\text{max}}$ shows a delayed spring phytoplankton bloom (Fig. C2) compared to the CTD observations. Additionally, zooplankton mortality produced differences in the timing of the spring phytoplankton bloom, with delays of 30 d (year 2014) and 35 d (year 2015) when there is less zooplankton mortality, affecting the timing and magnitude of the spring zooplankton bloom, therefore generating low values of surface chlorophyll a during spring (Fig. C5). It is well known that zooplankton are key players in the biogeochemical cycling of carbon and nutrients in marine ecosystems , influencing the export of organic matter to the deep ocean . Additionally, grazing responses comprise the dominant losses for phytoplankton in the ocean , influencing plankton stocks and primary production .

Table 2 presents differences for the S2P3 v8.0 model calibrated for the CCS location and for selected sensitivity experiments in terms of the timing of the spring phytoplankton bloom (days), defined as in Table B2, using a threshold for phytoplankton biomass ($>$ 1.5 mg Chl m${}^{-\mathrm{3}}$); the magnitude of the spring phytoplankton bloom (mg Chl m${}^{-\mathrm{3}}$); and the total annual zooplankton biomass (g DW m${}^{-\mathrm{2}}$). Representation of phytoplankton physiology had an important influence on the timing of the spring phytoplankton bloom (Table 2), with ${\mathit{\theta }}_{\text{max}}^{\text{N}}$ affecting the S2P3 v8.0 model the most in terms of this attribute of the model structure, showing less productive and delayed spring blooms when ${\mathit{\theta }}_{\text{max}}^{\text{N}}$ is lower (Figs. C7a, C8b). Changes in the timing and magnitude of the spring zooplankton bloom coincide with the changes of the timing and magnitude of the phytoplankton blooms (Fig. C7b). These changes in the plankton communities over the year due to different values of ${\mathit{\theta }}_{\text{max}}^{\text{N}}$, also have an effect on the values of DIN (Fig. C8c, f), with the largest differences shown during springtime at the surface (Fig. C8c). These differences agree with the results found by , where it was demonstrated that taking into account photo-acclimation and variable stoichiometry of phytoplankton growth in marine ecosystem models produces qualitative and quantitative differences in phytoplankton dynamics. Moreover, these quota formulations in S2P3 v8.0 were compared to the available dataset of physiological observations (Fig. C9). It is interesting to note that the sensitivity analysis of NPZ parameters produced differences in the physiological variables $P_{\text{max}}^{\text{Chl}}$ and $E_k$, specially at the surface (Figs. C3, C6), suggesting that the predator–prey interactions are indirectly influencing phytoplankton physiology, presumably through feedbacks between zooplankton and the nutrient cycling, which subsequently have an effect on phytoplankton physiology due to the dependency of nutrient quotas to the availability of inorganic nutrients. The current study thus demonstrates how a greater variety of data, spanning multiple trophic levels and incorporating information on physiological status as well as standing stocks, provides additional constraints on model validation and hence constraints on parameterisation.

Table 3

Comparison between the S2P3 v7.0, S2P3-NPZ, S2P3-Photoacclim, and S2P3 v8.0 models calibrated for the CCS location, in terms of the total spring NPP, total summer NPP, and total annual NPP calculated from 1965 to 2015, including the mean, maximum, minimum, and SD values.

Model	Characteristic (units)	Mean	Maximum	Minimum	SD
S2P3 v7.0	Total spring NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	39.6	54.7	33.5	5.8
	Total summer NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	17.2	22.9	56.9	5.1
	Total annual NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	57.1	61.1	53.1	1.8
S2P3-NPZ	Total spring NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	21.1	49.3	7.2	9.8
	Total summer NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	28.2	39.7	6.8	7.7
	Total annual NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	55.7	60.0	52.7	1.5
S2P3-Photoacclim	Total spring NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	35.6	40.7	25.7	4.8
	Total summer NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	3.8	11.4	0.6	4.2
	Total annual NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	39.4	42.4	37.0	1.7
S2P3 v8.0	Total spring NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	25.5	39.5	14.5	3.7
	Total summer NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	10.3	12.9	6.8	1.1
	Total annual NPP (g C m${}^{-\mathrm{2}}$ yr${}^{-\mathrm{1}}$)	37.8	47.8	33.4	2.2

6 Comparison of overall model performance

The S2P3 v7.0 , S2P3-NPZ, S2P3-Photoacclim, and S2P3 v8.0 models were calibrated for the CCS location and further analysis was undertaken by running each model for a extended period (1965–2015) to evaluate the statistics of productivity, partitioned between spring and summer. Table 3 shows that for the S2P3 v7.0 model, on average, 69.2 % of the annual phytoplankton production occurs during the spring phytoplankton bloom. On the other hand, for the S2P3-NPZ model, on average, only 37.8 % of the annual production occurs during spring months, showing that the predator–prey relationship has a strong influence on the magnitude of the spring phytoplankton bloom every year. Additionally, the S2P3-Photoacclim model shows a very strong spring phytoplankton bloom, corresponding to 90 % of the total annual NPP. Finally, for the S2P3 v8.0 model, only 67.4 % of the annual production corresponds to the spring bloom period. It is clear that, on average, the least productive model overall was S2P3 v8.0, followed by the S2P3-NPZ, S2P3 v7.0, and S2P3-Photoacclim. This shows the impact and complexity that the predator–prey relationship has on the model dynamics, with the addition of explicit zooplankton and their grazing activity as one of the main losses of phytoplankton . On the other hand, the S2P3 v7.0 model has, on average, more total annual NPP than the S2P3-NPZ model, suggesting that the influence of a constant grazing rate is not as strong in comparison to the one provided by the zooplankton grazing (NPZ framework), because the predator–prey relationship cannot be entirely represented in the S2P3 v7.0 model.

7 Conclusions

This study demonstrates that the combination of an NPZ framework, photo-acclimation, and flexible stoichiometry of phytoplankton in one model produces a better representation of the ecosystem based on the comparison to observations. This combined framework offers an improvement to the S2P3 v7.0 model for application to the CCS location and more broadly within shelf sea systems. The model validation, using both zooplankton biomass and physiological rates of phytoplankton observations, is rarely found in the literature, providing a novel contribution to the marine biogeochemistry modelling field of shelf seas. The development of the S2P3 v8.0 model provides a better fit to observations in comparison to the S2P3 v7.0, S2P3-NPZ, and S2P3-Photoacclim models. Improved confidence in the S2P3 v8.0 model thus suggests improved insights in studies about the effects of physical forcing through tides , intra- and interannual variations in meteorology and other drivers on PP, and phytoplankton dynamics would be possible.

Appropriate parameterisations to represent shelf seas is a subject that should be further supported by fieldwork campaigns and future work should aim to include additional datasets with longer time series . For this study, constraining the seasonal cycle of the phytoplankton physiology is not possible due to the lack of physiological observations during other periods of the year; furthermore, phytoplankton and zooplankton biomass datasets only include the years 2014 and 2015 but longer time series would help to improve the model calibration. Many model parameters quantities are poorly constrained observationally mainly due to the fact that model state variables are highly integrated pools, which are affected by biotic and abiotic factors in the environment, making them difficult to be determined by in situ measurements .

The tuning and sensitivity analysis performed in this work allows a better understanding of the ecosystem dynamics represented in the model and how it is influenced by each parameter. By considering the timing and magnitude of the spring phytoplankton bloom, the annual zooplankton biomass, and summer phytoplankton photosynthetic physiology as features to be calculated, a quantitative and clearer comparison between each model could be developed. Finally, the model calibration will never be in perfect agreement with all the observations, particularly in this case, where only one type of phytoplankton and zooplankton was considered. Thus, responses must represent some typical or average dynamics and cannot represent any effects of competition between types. Additionally, no stages of zooplankton growth were taken into account, which might affect the predator–prey interactions . Despite such potential limitations, the S2P3 v8.0 model is able to reasonably represent the integrated behaviour of the mixture of species that inhabit the NW European shelf sea.

Appendix A Model variables and parameters

Table A1

List of all variables for the biological part of the S2P3 v8.0 model.

Variables	Definition	Units
$I_{\text{PAR}}$	Photosynthetically available radiation	W m${}^{-\mathrm{2}}$
$K_z$	Depth-dependant vertical eddy diffusivity	m${}^{\mathrm{2}}$ s${}^{-\mathrm{1}}$
$N$	External DIN concentration	mmol N m${}^{-\mathrm{3}}$
Phyto${}_{\text{C}}$	Phytoplankton biomass in carbon currency	mg C m${}^{-\mathrm{3}}$
Phyto${}_{\text{Chl}}$	Phytoplankton biomass in chlorophyll currency	mg Chl m${}^{-\mathrm{3}}$
Phyto${}_{\text{N}}$	Phytoplankton biomass in nitrogen currency	mmol N m${}^{-\mathrm{3}}$
$Z$	Zooplankton biomass in nitrogen currency	mmol N m${}^{-\mathrm{3}}$
$I$	Ingestion rate of phytoplankton	d${}^{-\mathrm{1}}$
$P_{\text{max}}$	Carbon-specific rate of photosynthesis	d${}^{-\mathrm{1}}$
$u$	Phytoplankton carbon-specific nitrate uptake rate	mmol N (mg C)${}^{-\mathrm{1}}$ d${}^{-\mathrm{1}}$
$Q_P$	Cellular nutrient quota (N : C)	mmol N (mg C)${}^{-\mathrm{1}}$
$Q$	Internal nitrogen : phytoplankton chlorophyll ratio	mmol N (mg Chl)${}^{-\mathrm{1}}$
$z$	Vertical coordinate (positive upwards)	m
$\mathit{\theta }$	Chlorophyll : phytoplankton carbon ratio	mg Chl (mg C)${}^{-\mathrm{1}}$
$\mathit{\mu }$	Carbon-specific rate of photosynthesis	d${}^{-\mathrm{1}}$
${\mathit{\rho }}_{\text{Chl}}$	Chlorophyll synthesis regulation term	mg Chl (mmol N)${}^{-\mathrm{1}}$

Table A2

List of parameters, with their respective definitions, units, and initialised values for the biological part of the S2P3 v8.0 model.

Parameters	Definition	Units	Value
$k_{\text{n}}$	Half-saturation constant for nitrate uptake	mmol N m${}^{-\mathrm{3}}$	0.014
$P_{\text{max}}^{\text{C}}$	Maximum value of the carbon-specific rate of photosynthesis	d${}^{-\mathrm{1}}$	3.5
$Q_{\text{max}}$	Maximum value of the cellular nutrient quota	mmol N (mg C)${}^{-\mathrm{1}}$	0.032
$Q_{\text{min}}$	Minimum value of the cellular nutrient quota	mmol N (mg C)${}^{-\mathrm{1}}$	0.0028
$R_{\text{C}}$	Respiration rate constant	d${}^{-\mathrm{1}}$	0.02
$R_{\text{Chl}}$	Chlorophyll degradation rate constant	d${}^{-\mathrm{1}}$	0.02
$R_{\text{n}}$	Nitrate remineralisation rate constant	d${}^{-\mathrm{1}}$	0.02
$u_{\text{max}}$	Maximum value of the phytoplankton carbon-specific nitrate uptake rate	mmol N (mg C)${}^{-\mathrm{1}}$ d${}^{-\mathrm{1}}$	0.004
${\mathit{\alpha }}^{\text{Chl}}$	Chlorophyll-specific initial slope of the photosynthesis–light curve	mg C (mg Chl)${}^{-\mathrm{1}}$ d${}^{-\mathrm{1}}$	$\mathrm{1.99}\times {\mathrm{10}}^{-\mathrm{6}}$
		(W m${}^{-\mathrm{2}}{)}^{-\mathrm{1}}$
$\mathit{\zeta }$	Cost of biosynthesis	mg C (mmol N)${}^{-\mathrm{1}}$	0
${\mathit{\theta }}_{\text{max}}^{\text{N}}$	Maximum value of the chlorophyll : phytoplankton nitrogen ratio	mg Chl (mmol N)${}^{-\mathrm{1}}$	2.1
$T_{\text{func}}$	Temperature-response function	dimensionless	1.0
$m$	Loss rate of zooplankton due to predation and physiological death	d${}^{-\mathrm{1}}$	0.02
$R_{\text{max}}$	Zooplankton maximal grazing rate	d${}^{-\mathrm{1}}$	3.5
${\mathit{\gamma }}_{\mathrm{1}}$	Grazing inefficiency or “messy feeding” (0.0–1.0), returns a	dimensionless	0.1
	fraction of grazed material back into the DIN pool
${\mathit{\gamma }}_{\mathrm{2}}$	Fraction of dead zooplankton (0.0–1.0) that goes into the sediments	dimensionless	0.4
$\mathit{\lambda }$	Rate at which saturation is achieved with increasing food levels	(mmol N m${}^{-\mathrm{3}}{)}^{-\mathrm{1}}$	0.014

Appendix B Model calibration

Table B1

List of relevant CTD casts for the CCS location from DY029 and DY033 cruises considering the date, location (latitude and longitude), and depth. CTD casts in bold are the ones chosen in this work to validate each model during spring and summer.

Cruise name	Date	Latitude (${}^{\circ }$ N)	Longitude (${}^{\circ }$ W)	Depth (m)
DY029	3 Apr 2015	49.38	8.59	147
DY029	4 Apr 2015	49.38	8.59	146
DY029	5 Apr 2015	49.38	8.59	146
DY029	6 Apr 2015	49.4	8.58	147
DY029	11 Apr 2015	49.39	8.58	145
DY029	15 Apr 2015	49.4	8.59	147
DY029	20 Apr 2015	49.4	8.6	147
DY029	21 Apr 2015	49.4	8.62	148
DY029	25 Apr 2015	49.4	8.59	148
DY029	26 Apr 2015	49.4	8.58	146
DY029	28 Apr 2015	49.4	8.58	146
DY033	13 Jul 2015	49.43	8.59	144
DY033	14 Jul 2015	49.42	8.54	144
DY033	15 Jul 2015	49.37	8.61	145
DY033	24 Jul 2015	49.36	8.62	145
DY033	25 Jul 2015	49.41	8.59	148
DY033	29 Jul 2015	49.42	8.57	147
DY033	30 Jul 2015	49.4	8.57	148
DY033	1 Aug 2015	49.38	8.58	146

Table B2

Quantitative comparison of the timing of the spring phytoplankton bloom between the buoy observations, S2P3-NPZ, S2P3-Photoacclim, and S2P3 v8.0 models.

Years	Buoy	S2P3-NPZ	S2P3-Photoacclim	S2P3 v8.0
2014	5 April	14 April	25 April	3 April
2015	16 March	24 April	6 May	9 April

Appendix C Sensitivity studies

Figure C1

(a) SSB observations of surface chlorophyll a (red line), along with the modelled surface chlorophyll a for the calibrated S2P3 v8.0 model (green line), and experiments $R_{\text{max}}=\mathrm{1.75}$ (black line) and $R_{\text{max}}=\mathrm{5.25}$ (cyan line). (b) Observations of zooplankton biomass presented as discrete points for day nets (red dots) and night nets (blue dots) taken during the cruises DY026, DY018, DY029, and DY033; modelled zooplankton biomass from the calibrated S2P3 v8.0 model (green line) and experiments $R_{\text{max}}=\mathrm{1.75}$ (black line) and $R_{\text{max}}=\mathrm{5.25}$ (cyan line).

[Figure omitted. See PDF]

Figure C2

CTD observations from the SSB programme (red line) including data for springtime (20 April 2015) (a) temperature, (b) chlorophyll a, and (c) DIN (red dots); for summertime (24 July 2015) for (d) temperature, (e) chlorophyll a, and (f) DIN (red dots) along the calibrated S2P3 v8.0 model (green line), experiments $R_{\text{max}}=\mathrm{1.75}$ (black line) and $R_{\text{max}}=\mathrm{5.25}$ (cyan line).

[Figure omitted. See PDF]

Figure C3

Observations from the cruises CD173 and JR98 in different locations in the Celtic Sea, including the calibrated S2P3 v8.0 model (green lines), experiments $R_{\text{max}}=\mathrm{1.75}$ (black lines) and $R_{\text{max}}=\mathrm{5.25}$ (cyan lines) for (a) chlorophyll a specific maximum light-saturated photosynthesis rate $\left(P_{\text{max}}^{\text{Chl}}\right)$ and (b) light saturation parameter ($E_k$).

[Figure omitted. See PDF]

Figure C4

(a) SSB observations of surface chlorophyll a (red line), along with the modelled surface chlorophyll a for the calibrated S2P3 v8.0 model (green line), and experiments $m=\mathrm{0.01}$ (black line) and $m=\mathrm{0.03}$ (cyan line). (b) Observations of zooplankton biomass presented as discrete points for day nets (red dots) and night nets (blue dots) taken during the cruises DY026, DY018, DY029, and DY033; modelled zooplankton biomass from the calibrated S2P3 v8.0 model (green line) and experiments $m=\mathrm{0.01}$ (black line) and $m=\mathrm{0.03}$ (cyan line).

[Figure omitted. See PDF]

Figure C5

CTD observations from the SSB programme (red line) including data for springtime (20 April 2015) (a) temperature, (b) chlorophyll a, and (c) DIN (red dots); for summertime (24 Jul 2015) for (d) temperature, (e) chlorophyll a, and (f) DIN (red dots) along the calibrated S2P3 v8.0 model (green line), experiments $m=\mathrm{0.01}$ (black line) and $m=\mathrm{0.03}$ (cyan line).

[Figure omitted. See PDF]

Figure C6

Observations from the cruises CD173 and JR98 in different locations in the Celtic Sea, including the calibrated S2P3 v8.0 model (green lines), experiments $m=\mathrm{0.01}$ (black lines) and $m=\mathrm{0.03}$ (cyan lines) for (a) chlorophyll a specific maximum light-saturated photosynthesis rate ($P_{\text{max}}^{\text{Chl}}$) and (b) light saturation parameter ($E_k$).

[Figure omitted. See PDF]

Figure C7

(a) SSB observations of surface chlorophyll a (red line), along with the modelled surface chlorophyll a for the calibrated S2P3 v8.0 model (green line), and experiments ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.075}$ (black line) and ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.225}$ (cyan line). (b) Observations of zooplankton biomass presented as discrete points for day nets (red dots) and night nets (blue dots) taken during the cruises DY026, DY018, DY029, and DY033; modelled zooplankton biomass from the calibrated S2P3 v8.0 model (green line) and experiments ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.075}$ (black line) and ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.225}$ (cyan line).

[Figure omitted. See PDF]

Figure C8

CTD observations from the SSB programme (red line) including data for springtime (20 April 2015) (a) temperature, (b) chlorophyll a, and (c) DIN (red dots); for summertime (24 July 2015) for (d) temperature, (e) chlorophyll a, and (f) DIN (red dots) along the calibrated S2P3 v8.0 model (green line), experiments ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.075}$ (black line) and ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.225}$ (cyan line).

[Figure omitted. See PDF]

Figure C9

Observations from the cruises CD173 and JR98 in different locations in the Celtic Sea, including the calibrated S2P3 v8.0 model (green lines), experiments ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.075}$ (black lines) and ${\mathit{\theta }}_{\text{max}}^{\text{N}}=\mathrm{0.225}$ (cyan lines) form (a) chlorophyll a specific maximum light-saturated photosynthesis rate $\left(P_{\text{max}}^{\text{Chl}}\right)$ and (b) light saturation parameter ($E_k$).

[Figure omitted. See PDF]

Code and data availability

The current version of model is available from https://doi.org/10.5281/zenodo.3600467 under the Creative Commons Attribution 4.0 international license. The exact version of the model used to produce the results used in this paper is archived on Zenodo (https://zenodo.org/record/3600467#.XyqbKhNKhBw, last access: 10 August 2020), as are input data and scripts to run the model and produce the plots for all the simulations presented in this paper.

Unzipped and uncompressed, the directory /s2p3v8.0 contains several subdirectories:

• /main contains the source code, s2p3_v8.f90, which is compiled “stand-alone”, and executed using accompanying scripts.

• /domain contains location (latitude and longitude), tidal components, and bathymetry (total depth) for the central Celtic Sea in the western English Channel (s12_m2_s2_n2_h_tim.dat).

• /met contains 2000–2015 meteorological forcing (Falmouth_met_i143_j17_2000-2015 ASCII file).

• /output contains example output data from experiments that show the calibrated version of each model (S2P3v8, S2P3-Photoacclim, and S2P3-NPZ).

• /plotting contains MATLAB scripts for plotting time series data from experiments S2P3v8, S2P3-Photoacclim, and S2P3-NPZ. This folder also include all the datasets from observations used to validate each model.

The ancillary files needed for simulations are available on request from the author (e-mail: [email protected]).

Author contributions

AABD developed the model code, performed the simulations, and prepared the manuscript with contributions from all co-authors.

Competing interests

The authors declare that they have no conflict of interest.

Acknowledgements

Hydrographic and biological observations used were acquired through the Natural Environment Research Council (NERC) funded UK Shelf Sea Biogeochemistry (SSB) research programme. Data used were requested from and are available from the British Oceanographic Data Centre (BODC). Sari Giering is thanked both for original collection and generation of the zooplankton data and discussion on interpretation of this data. NCEP/NCAR reanalysis 1 data were provided by the NOAA/OAR/ESRL PSD, Boulder, Colorado, USA, from their website at https://www.esrl.noaa.gov/psd/ (last access: December 2019). We thank the topical editor Paul Halloran and two anonymous reviewers for a series of insightful comments that helped to improve the manuscript.

Financial support

This research has been supported by the CONICYT PFCHA/DOCTORADO BECAS CHILE/2015 (grant no. 72160249).

Review statement

This paper was edited by Paul Halloran and reviewed by two anonymous referees.