Overview

StrathE2E2 is an R package for modelling the ‘big-picture’, whole ecosystem effects of hydrodynamics, temperature, nutrient additions, and fishing on continental shelf marine food webs.
StrathE2E2 has two linked parts - a fishing fleet model and an ecology model. The fishing model integrates harvesting, discarding and seabed disturbance rates across a range of gears and passes the results into the ecology model.
The ecology model is a network of coupled ordinary differential equations representing the rates of change in organic detritus, dissolved inorganic nutrient, and coarse guilds of living biomass spanning microbes to megafauna, in two coarse spatial compartments (inshore and offshore) and eight seabed habitats.
The equations which define the model include representations of feeding, metabolism, reproduction, active migrations, advection and mixing.
Environmental driving data include temperature, irradiance, hydrodynamics, and nutrient inputs from rivers, atmosphere and ocean boundaries.
The package includes functions for parameter optimization, global sensitivity analysis, and Monte Carlo estimation of credible intervals for model outputs.
This tutorial, which should take about an hour, provides an introduction to using the package to run ecosystem models, write your own code to develop sceanrios, and work with the model outputs.
Support for the developmemt of this course came from the EU Horizon 2020 projects Mission Atlantic (https://missionatlantic.eu/) and SEAwise (https://seawiseproject.org/)

Requirements

Laptop with R version 4.3 or higher installed
Familiarity with writing scripts in R to manipulate data and generate graphs
Familiarity with the the model concepts and structure - which can be gained from our online course at https://outreach.mathstat.strath.ac.uk/apps/StrathE2EApp/

Installing the package

Install the package from within an R session using:

install.packages("StrathE2E2", 
                 repos="https://www.marineresourcemodelling.maths.strath.ac.uk/sran/")

Alternatively, download the installation zip file for your operating system from https://www.marineresourcemodelling.maths.strath.ac.uk/strathe2e/index.html and install into your R version.

Note that Strathe2E2 requires installation of the following additional packages which are available from any CRAN site:

deSolve
NetIndices

Further information on building the package from the source code is available at the StrathE2E GitLab site and repository

Once installed, load the package before using it by:

library(StrathE2E2)

Installing the supplementary data-package

The model package has a supplementary data-package - Strathe2E2examples - which contains example outputs from the more computationally intensive functions in the package. These are provided to illustrate the performance of various post-processing and plotting functions for these analyses that are available.

Manually install the package using:

install.packages("StrathE2E2examples", 
                 repos="https://www.marineresourcemodelling.maths.strath.ac.uk/sran/")

Once installed, details of the data sets are available by typing:

library(StrathE2E2examples)

Details of the datasets are available by typing:

help(StrathE2E2examples)

Initial tour of the package

To see an overview of the StrathE2E2 package, launch an R-session and then:

library(StrathE2E2)     # Load the package
help(StrathE2E2)        # Launch the overview document in a browser window

The package overview contains links to help and example pages for each of the available functions. These are also acessible by typing, for example:

help(e2e_run)        # Launch the help page for the e2e_run function in a browser window

Cheatsheet - quick-reference guide to StrathE2E2 functions

Reference Manual - compilation of the help pages for all the package functions generated by CRAN-check

User Manual - guide to using the package

Technical Manual - documentation on the structure of input and output R-objects and files

In addition the following supporting documents are available:

Download documentation for the North Sea model which is provided in the package

Publications inventory for StrathE2E models and applications

How to cite the package:

citation("StrathE2E2")
To cite package 'StrathE2E2' in publications use:

  Heath, M.R., Speirs, D.C., Thurlbeck, I. and Wilson, R. (2021).
  StrathE2E2: an R package for modelling the dynamics of marine food
  webs and fisheries. Methods in Ecology and Evolution, 12, 280-287.

A BibTeX entry for LaTeX users is

  @Article{,
    title = {StrathE2E2: an R package for modelling the dynamics of marine food webs 
          and fisheries},
    author = {Michael R. Heath and Douglas C. Speirs and Ian Thurlbeck and Robert J. Wilson},
    journal = {Methods in Ecology and Evolution},
    volume = {12},
    pages = {280-287},
    year = {2021},
    issn = {2041-210X},
    url = {https://besjournals.onlinelibrary.wiley.com/doi/10.1111/2041-210X.13510},
    keywords = {food web, dynamic model, ecosystem, fisheries, 
                 ordinary differential equations, optimization, 
                 sensitivity analysis, Monte-Carlo, R},
  }

EXERCISE 1. Spend a few minutes roaming around the package help, starting from help(StrathE2E2)

Let’s get started

Loading the package

To start a session, begin by loading the package:

library(StrathE2E2)

To check the details of the package you have installed:

e2e_info_pk()
Package: StrathE2E2
Version: 4.0.1
Title: End-to-End Marine Food Web Model
Author: Michael Heath [aut], Ian Thurlbeck [ctb]
Maintainer: Michael Heath <m.heath@strath.ac.uk>
Depends: R (>= 3.6.0),
Imports: deSolve, methods, NetIndices
Description: A dynamic model of the big-picture, whole ecosystem
        effects of hydrodynamics, temperature, nutrients, and fishing
        on continental shelf marine food webs. The package is described
        in: Heath, M.R., Speirs, D.C., Thurlbeck, I. and Wilson, R.J.
        (2021) <doi:10.1111/2041-210X.13510> StrathE2E2: An R package
        for modelling the dynamics of marine food webs and fisheries.
        Methods in Ecology and Evolution 12, 280-287.
License: GPL (>= 2)
NeedsCompilation: yes
Date: 2025-01-04
URL: https://www.marineresourcemodelling.maths.strath.ac.uk/strathe2e/
Suggests: knitr, StrathE2E2examples, testthat (>= 2.1.0)
VignetteBuilder: knitr
Encoding: UTF-8
RoxygenNote: 7.3.2
Additional_repositories:
        https://www.marineresourcemodelling.maths.strath.ac.uk/strathe2e/articles/Pkg_versions.html
Packaged: 2025-01-14 21:55:36 UTC; ais04103
Built: R 4.4.2; x86_64-w64-mingw32; 2025-01-14 21:55:39 UTC; windows
Archs: x64

-- File: C:/Users/ais04103/AppData/Local/R/win-library/4.4/StrathE2E2/Meta/package.rds

Load a model

The package has two variants of a model of the North Sea embedded within it:

e2e_ls()
Current working directory is 'C:/Users/ais04103/OneDrive - University of Strathclyde/Documents/GitLab/Downloads/public/resources/StrathE2E2/documents/training/source_files'
List of package models in system folder "extdata/Models", with examples of how to read them:
 Model: "North_Sea"
  Variant: "1970-1999"  model <- e2e_read("North_Sea", "1970-1999")
  Variant: "2003-2013"  model <- e2e_read("North_Sea", "2003-2013")

You can view details of these models using the function e2e_info_md():

e2e_info_md()
Path to the model is: 'C:/Users/ais04103/AppData/Local/R/win-library/4.4/StrathE2E2/extdata/Models'
 [1] "Model region      : North_Sea"                                                                                                                                       
 [2] "StrathE2E2 version: 4.0.1"                                                                                                                                           
 [3] ""                                                                                                                                                                    
 [4] "Model variant name: 1970-1999"                                                                                                                                       
 [5] "Release date      : 30-05-2022"                                                                                                                                      
 [6] "Implementation originally produced in the NERC Marine Ecosystem Programme (https://marine-ecosystems.org.uk/)"                                                       
 [7] "Driving data represent a climatological annual cycle for the period 1970-1999"                                                                                       
 [8] "Hydrodynamnic driving data from NEMO-ERSEM 7km model (https://www.uk-ssb.org/science_components/work_package_4/)"                                                    
 [9] "Nutrient driving data from ICES and BODC data archives and informed by NEMO-ERSEM"                                                                                   
[10] "Atmospheric driving data from EMEP (https://www.emep.int/)"                                                                                                          
[11] "River driving data from a statistical models of rainfall river flow and nutrent monitoring data (https://strathprints.strath.ac.uk/18588/6/strathprints018588.pdf)"  
[12] "Fishing data from ICES and STECF (https://www.ices.dk/data/dataset-collections/Pages/Fish-catch-and-stock-assessment.aspx and https://stecf.jrc.ec.europa.eu/dd/fdi)"
[13] "Setup parameters modified May 2022 to be consistent with StrathE2E 4.0.0"                                                                                            
[14] "Documentation: https://www.marineresourcemodelling.maths.strath.ac.uk/resources/StrathE2E2/documents/4.0.0/StrathE2E2_North_Sea_model.pdf"                           
[15] ""                                                                                                                                                                    
[16] "Model variant name: 2003-2013"                                                                                                                                       
[17] "Release date      : 30-05-2022"                                                                                                                                      
[18] "Implementation originally produced in the NERC Marine Ecosystem Programme (https://marine-ecosystems.org.uk/)"                                                       
[19] "Driving data represent a climatological annual cycle for the period 2003-2013"                                                                                       
[20] "Hydrodynamnic driving data from NEMO-ERSEM 7km model (https://www.uk-ssb.org/science_components/work_package_4/)"                                                    
[21] "Nutrient driving data from ICES and BODC data archives and informed by NEMO-ERSEM"                                                                                   
[22] "Atmospheric driving data from EMEP (https://www.emep.int/)"                                                                                                          
[23] "River driving data from a statistical models of rainfall river flow and nutrent monitoring data (https://strathprints.strath.ac.uk/18588/6/strathprints018588.pdf)"  
[24] "Fishing data from ICES and STECF (https://www.ices.dk/data/dataset-collections/Pages/Fish-catch-and-stock-assessment.aspx and https://stecf.jrc.ec.europa.eu/dd/fdi)"
[25] "Setup parameters modified May 2022 to be consistent with StrathE2E 4.0.0"                                                                                            
[26] "Documentation: https://www.marineresourcemodelling.maths.strath.ac.uk/resources/StrathE2E2/documents/4.0.0/StrathE2E2_North_Sea_model.pdf"                           
[27] ""

Here is a map of the North Sea model domain:

Figure 1. Domain of the North Sea model include in the package showing the seabed habiats in the inshore and offshore zones.

e2e_read() is a key function in the package which returns a list object containing all of the input data and parameters gathered from a set of csv files that define the model implementation.

To read the 1970-1999 North Sea model:

model <- e2e_read("North_Sea", "1970-1999")   # Read the 1970-1999 North Sea model from within the package
Current working directory is... 
'C:/Users/ais04103/OneDrive - University of Strathclyde/Documents/GitLab/Downloads/public/resources/StrathE2E2/documents/training/source_files'
No 'results.path' specified so any csv data requested
will be directed to/from the temporary directory...
'C:\Users\ais04103\AppData\Local\Temp\RtmpsXcJPr'

Model setup and parameters gathered from ...
StrathE2E2 package folder
Model results will be directed to/from ...
'C:\Users\ais04103\AppData\Local\Temp\RtmpsXcJPr/North_Sea/1970-1999/'

Viewing the model input data

You can view the climatological annual cycles of input environmental driving data in the model object using the function e2e_plot_edrivers(). This displays the following data:

Table 1. Environmental input data dispalayed in each panel of the plot produced by e2e_plot_edrivers()

Data	Units	Comments
Sea surface irradiance	$\mu$E.m^-2.d^-1
Suspended particulate matter	g.m^-3
Temperature	deg-C
Vertical diffusivity gradient	m.d^-1	derived from the vertical diffusivity (m².s^-1) and mixing length scale (m)
External inflows	m³ per m² sea surface of model domain	derived from proportion input per layer volume, layer thicknesses and areas
River discharge	m³ per m² sea surface of model domain	derived from proportion input to inshore volume, and inshore layer thickness and area
Inshore significant wave height	m
Proportion of seabed disturbed	d^-1	aggregated over the three sediment classes in each zone
External boundary nitrate concentration	mMN.m^-3
External boundary ammonia concentration	mMN.m^-3
External boundary phytoplankton concentration	mMN.m^-3
External boundary detritus concentration	mMN.m^-3
River nitrate concentration	mMN.m^-3
River ammonia concentration	mMN.m^-3
Atmospheric nitrate deposition flux	mMN.m^-2.d^-1	area averaged flux density over each zone
Atmospheric ammonia deposition flux	mMN.m^-2.d^-1	area averaged flux density over each zone

e2e_plot_edrivers(model)

Figure 2. Plot of climatological annual cycles of the environmental drivers of the 1970-1999 North Sea model

The other main sets of input data to the model concerns the fishing fleets. The fishing fleet model can represent up to 12 different categories of fishing, each defined by its:

activity rate
proportional spatial distribution across the inshore and offshore zone and seabed habitats
catching power for each of the catchable guilds in the model (i.e. selectivity)
discard rate for each of the catchable guilds
processing-at-sea rate for each of the catchable guilds (which generates offal)
seabed abrasion rate
proportionality coefficient relating effort (activity x power) to fishing mortality rate on each guild

The fishing fleet model integrates across the fleets and generates the overall harvest ratio (equivalent to fishing mortality rate), discard and offal production rates for each guild, and the seabed abrasion rates for each seabed habitat.

You can view the fishing fleet input data using the function e2e_plot_fdrivers(). This function has two arguments:

Table 2. Arguments of the functions e2e_plot_fdrivers()

Argument	Description
model	R-list object defining a model configuration and compiled by the e2e_read() function
selection	Text string from a list identifying the class of fishing-related driving data to be plotted. Select from: “ACTIVITY”, “ABRASION”, “HARVESTR”, “DISCARDS”, “OFFAL”. Remember to include the phrase within “” quotes

Regardless of the variable selected, the function returns both a graphical image and a list object containing the matix of data plotted and vectors of the axis labels, which uses can use to generate your own plots if desired.

Example use of e2e_plot_fdrivers() using the argument selection = “ACTIVITY”. This plots a matrix of seabed habitats vs fishing gears with each cell shaded to indicate $log_e$ transformed activity density (white = 0, purple = high).

The spatial distribution of fishing gear activity in the model is defined by two input data sets: * Vector of whole domain activity density of each fishing gear (s.d^-1 per m² of whole model domain) * Matrix of the proportional distribution of whole domain activity density of each gear (rows) across seabed habitats (columns)

The activity density of each gear in each habitat is then obtained by multiplying the whole domain activity density into the proportional distribution matrix, and dividing by the area-proportions of habitats in the domain. The units of habitat-specific activity density are then (s.d^-1.m^-2).

# Load the 1970-1999 version of the North Sea model supplied with the package, 
# run, and generate a plot:
plotted_data<-e2e_plot_fdrivers(model, selection="ACTIVITY")

Figure 3. Distribution of gear activity rates across habitats in the 1970-1999 North Sea model.

EXERCISE 2. Read the 2003-2013 variant of the North Sea model and have a go at plotting the various aspects the input data.

Run a model and plot basic results

e2e_run() is a key function in the package wihich returns a list object which contains all of the output data from a model run.

e2e_plot_ts() with the arguments model object name, results object name and “ECO” produces a coarse diagnostic summary of model outputs - time series of the daily-interval masses of each of the main state variables in the model.

To load the 1970-1999 North Sea model, run and view the results:

model <- e2e_read("North_Sea", "1970-1999",silent=TRUE) # Read the 1970-1999 North Sea model from within the package
results <- e2e_run(model,nyears=3)                      # Run the model for 3 years
e2e_plot_ts(model,results,"ECO")                        # Plot time series of the output variables

Figure 4. Time series output from the 1970-1999 North Sea model included in the package.

You will note that these model results show exactly repeating annual cycles of outputs. We refer to this as a steady or stationary state. The driving input data to the model are exactly repeating annual cycles of climatologically average environmental conditions in this case for the period 1970-1999 in the North Sea. The model has previously been run until the the annual cycles of output are also exactly repeating. Then, the values of all the variables at the end of a year are saved and used for subsequent re-starts of the model.

Replacing the argument “ECO” with “CATCH” produces a summary of the time series of annual integrated fishery catch outputs from the model brohen down by spatial zone, landings and disacards:

e2e_plot_ts(model, results,selection="CATCH")

Figure 5. Time-series plot of catch data output, units mMN.y^-1 in the model domain as a whole (1m^-2 sea surface).

EXERCISE 3. Run the 2003-2013 variant of the North Sea model and plot the basic results.

More detailed visualization of model outputs

For the final year of a model run we can generate more detailed visualizations of the model results using the functions:

e2e_plot_eco() - Plot annual cycles of ecological ouputs
e2e_plot_catch() - Plot annual catches disaggregared by gear, guild, landings, discards
e2e_plot_biomass() - Plot annual mean biomass of all guilds by inshore/offshore zone
e2e_plot_migration() - Plot annual cycles of active migration rates between zones
e2e_plot_trophic() - Plot annual mean trophic level and omnivory index of each guild

All of these functions have arguments to control the use of credible interval data from the examples package - which we will come to later. Here we just need the arguments to control plotting data from a single run.

Final year time series of ecological data (e2e_plot_eco())

Table 3. Relevant arguments of the functions e2e_plot_eco()

Argument	Description
model	R-list object defining the model configuration compiled by the e2e_read() function
selection	Text string from a list identifying the group of model output variables to be plotted. Select from: “NUT_PHYT”, “SEDIMENT”, “ZOOPLANKTON”, “FISH”, “BENTHOS”, “PREDATORS”, “CORP_DISC”, “MACROPHYTE”, Remember to include the phrase within “” quotes
results	R-list object of model output generated by the e2e_run() function

Table 4. Variables plotted (and units) given values of the argument “selection”

selection	Variables plotted	Units
NUT_PHYT	Water column nitrate, ammonia, detritus and phytoplankton	mMN.m^-3
SEDIMENT	Sediment porewater nitrate and ammonia, detritus and corpses	mMN.m^-3 or gN.gDW^-1
ZOOPLANKTON	Omnivorous and carnivorous zooplankton	mMN.m^-2
FISH	Planktivorous and demersal fish and fish larvae	mMN.m^-2
BENTHOS	Susp/dep and carn/scav benthos and benthos larvae	mMN.m^-2
PREDATORS	Birds, pinnipeds, cetaceans and migratory fish	mMN.m^-2
CORP_DISC	Corpses and discards	mMN.m^-2
MACROPHYTE	Macrphytes and macrophyte debris	mMN.m^-2

Example plot of final year nutrient and phytoplankton time series

# Load the 1970-1999 version of the North Sea model supplied with the package, 
# run, and generate a plot:
model <- e2e_read("North_Sea", "1970-1999",silent=TRUE)
results <- e2e_run(model, nyears=3)
e2e_plot_eco(model, selection="NUT_PHYT",results=results)

Figure 6. Annual cycles of nutrients and phytoplankton in the final year of a 1970-1999 North Sea model.

Plotting distributions of catches by gear and guild e2e_plot_catch()

Table 5. Arguments of the functions e2e_plot_catch()

Argument	Description
model	R-list object defining the model configuration compiled by the e2e_read() function
results	List object of single-run model output generated by running the function e2e_run() function
selection	Text string from a list identifying the group of model output variables to be plotted. Select from: “BY_GUILD”, “BY_GEAR”. With the former, each panel represents a different guild; with the latter, each panel represnts a different gear. Remember to include the phrase within “” quotes

Create stacked barplots of the distribution of offshore and inshore landings and discards across gears by guild, or across guilds by gears, in the final year of a model run.

Data in zonal landings and discards by guild and gear in the final year of a model run are generated as a standard output from the model, and saved both as csv files and in the results object returned by the e2e_run() function. This function organises these data in two different ways to display as barplots.

The first display is a multi-panel plot in which each panel represents a different guild, and the bars show the zonal landings by each gear. The alternative display has each panel as a different gear, and the bars show the landings and discards of each guild.

The unit of the displayed data are mMN.y^-1 from the model domain as a whole, which is taken as being 1 m².

e2e_plot_catch(model, results, selection="BY_GEAR")

Figure 7. Annual landings and discards of each guild, by gear, in the final year of a run of the 1970-1999 North Sea model.

e2e_plot_catch(model, results, selection="BY_GUILD")

Figure 8. Annual landings and discards by each gear, by guild, in the final year of a run of the 1970-1999 North Sea model.

Zonal distributions of annual average biomass e2e_plot_biomass()

Table 6. Relevant arguments of the functions e2e_plot_biomass()

Argument	Description
model	R-list object defining the model configuration compiled by the e2e_read() function
results	R-list object of model output generated by the e2e_run() function

Generate plots showing the annual average biomass densities of each guild in the ecology model in the inshore (blue) and offshore (red) zones during the final year of a run.

Note that in this plot the biomass are expressed as mMN.m^-2, meaning that the mass in each zone has been scaled to the zonal sea surface area. So the data in this plot are area densities and not mass.

e2e_plot_biomass(model, results=results)

Figure 9. Annual average biomass densities for the final year of a single run of the 1970-1999 North Sea model.

Final year time series of migration fluxes (e2e_plot_migration())

Table 7. Relevant arguments of the functions e2e_plot_migration()

Argument	Description
model	R-list object defining the model configuration compiled by the e2e_read() function
results	R-list object of model output generated by the e2e_run() function

Each panel of the plot shows a time-series of the net flux densities (mMN.d^-1 in the model domain as a whole, assumed to be 1 m² sea surface area) of one of the migratory guilds in the model (all three guilds of fish, birds, pinnipeds and cetaceans) between the inshore and offshore zones of the model, over the final year of a run. These migration fluxes are the dynamic product of gradients in the ratio of food concentration to predator concentration across the inshore-offshore boundary. Positive values of the net migration flux indicate net movement from the offshore to inshore zone. Negative values indicate net movement from inshore to offshore.

e2e_plot_migration(model, results=results)

Figure 10. Annual cycles of migration fluxes between the offshore and inshore zones in the final year of a 1970-1999 North Sea model.

Plotting mean trophic level and omnivory indices e2e_plot_trophic()

Table 8. Relevant arguments of the functions e2e_plot_trophic()

Argument	Description
model	R-list object defining the model configuration compiled by the e2e_read() function
results	R-list object of model output generated by the e2e_run() function

Generate a two-panel plot showing: (upper panel) the mean trophic level of each guild in the ecology model, and (lower panel) the omnivory index of each guild. The data are generated by the NetIndices package from a flow matrix of nutrient fluxes through, into and out of the ecosystem during the final year of a run. The data are generated automatically as part of the output from every call of the e2e_run() function.

Note the the NetINdices package assigns trophic level 1 to inorganic nutrients and detritus, to autotrophs and detritivores are ranked as trophic level 2.

e2e_plot_trophic(model, results=results)

Figure 11. Annual mean trophic level and omnivory indices for the final year of a single run of the 1970-1999 North Sea model.

EXERCISE 4. Spend a few minutes exploring the final year plotting functions with output from a North Sea model run. Use help(function name) for information and assistance.

e2e_plot_eco() - Plot annual cycles of ecological ouputs
e2e_plot_catch() - Plot annual catches disaggregared by gear, guild, landings, discards
e2e_plot_biomass() - Plot annual mean biomass of all guilds by inshore/offshore zone
e2e_plot_migration() - Plot annual cycles of active migration rates between zones
e2e_plot_trophic() - Plot annual mean trophic level and omnivory index of each guild

Writing code to create model scenarios

Typical use of the model might involve

comparing a baseline run with a scenario run involving some changes in driving data (e.g. different temperature conditions or different activity rates of selected gears),
conducting an analysis of the sensitivity to systematic changes in driving data (e.g. incrememts in river nutrient inputs over some range of values).

The most efficient way to create model scenarios is to write scripts to modify the R list object created by the e2e_read() function, assigning a unique identifier to the outputs from each run.

The scope for developing fishing and/or environmemtal scenarios by coding to modify a default model object is limitless, but requires an understanding of the structure and content of the list object produced by the e2e_read() function. Use the R str() function to explore the list and its sub-components, and consult the StrathE2E2 Technical Manual.

model<-e2e_read("North_Sea", "2003-2013",silent=TRUE)
str(model,max.level=1)                    # View the top 2 levels of input list object
List of 2
 $ setup:List of 7
 $ data :List of 9

str(model$data,max.level=1)                    # View the contents of the **data** level in the mopdel object
List of 9
 $ fixed.parameters   :List of 51
 $ fitted.parameters  :List of 171
 $ physical.parameters:List of 66
 $ physics.drivers    :'data.frame':    12 obs. of  28 variables:
 $ chemistry.drivers  :'data.frame':    12 obs. of  27 variables:
 $ biological.events  :'data.frame':    26 obs. of  2 variables:
 $ fleet.model        :List of 19
 $ initial.state      :List of 416
 $ biology.scenarios  :List of 2

To mine deeper into the model object, for example, to show the fleet model inputs:

str(model$data$fleet.model,max.level=1)
List of 19
 $ gear_labels                : chr [1:12] "Pelagic_Trawl+Seine" "Sandeel+sprat_trawl(Otter30-70mm+TR3)" "Longline_mackerel" "Beam_Trawl_BT1+BT2" ...
 $ gear_codes                 : chr [1:12] "PTS" "SST" "LLm" "BTf" ...
 $ gear_activity              : num [1:12] 2.17e-06 4.23e-06 1.68e-06 1.15e-05 1.72e-08 2.16e-05 7.92e-06 1.27e-05 1.72e-05 2.40e-05 ...
 $ gear_group_rel_power       :'data.frame':    12 obs. of  10 variables:
 $ gear_group_discard         :'data.frame':    12 obs. of  10 variables:
 $ gear_group_gutting         :'data.frame':    12 obs. of  10 variables:
 $ gear_ploughing_rate        : num [1:12] 0 8.8 0 54.1 22.4 17.1 0 13.5 78.9 0 ...
 $ gear_habitat_activity      :'data.frame':    12 obs. of  9 variables:
 $ HRscale_vector             : Named num [1:10] 0.0697 0.0755 0.3696 12.662 0.6301 ...
  ..- attr(*, "names")= chr [1:10] "PF_HR_scale" "DF_HR_scale" "MF_HR_scale" "SB_HR_scale" ...
 $ HRscale_vector_multiplier  : int [1:10] 1 1 1 1 1 1 1 1 1 1
 $ offal_prop_live_weight     : num 0.1
 $ gear_mult                  : int [1:12] 1 1 1 1 1 1 1 1 1 1 ...
 $ quota_nonquota_parms_vector: num [1:6] 0.16 0.07 0.67 0.02 0.4 0.02
 $ DFsize_SWITCH              : num 0
 $ DFdiscard_SWITCH           : num 1
 $ plough_thickness           : num 0.05
 $ plough_depth_vector        : num [1:9] 0 0.4316 0.0582 0.0455 0 ...
 $ BSmort_gear                : num 0.2
 $ BCmort_gear                : num 0.2

Below is an example of code to configure and run a scenario case of the the 2003-2013 North Sea model. It involves editing the fleet.model data to set the activity of demersal beam trawlers to zero.

The identities of the fishing fleets in the model are listed in:

model$data$fleet.model$gear_labels
 [1] "Pelagic_Trawl+Seine"                  
 [2] "Sandeel+sprat_trawl(Otter30-70mm+TR3)"
 [3] "Longline_mackerel"                    
 [4] "Beam_Trawl_BT1+BT2"                   
 [5] "Demersal_Seine"                       
 [6] "Demersal_Otter_Trawl_TR1"             
 [7] "Gill_Nets+Longline_demersal"          
 [8] "Beam_Trawl_shrimp"                    
 [9] "Nephrops_Trawl_TR2"                   
[10] "Creels"                               
[11] "Mollusc_Dredge"                       
[12] "Whaler"

The activity rates of each gear in the ‘out of the box’ model are listed in:

model$data$fleet.model$gear_activity
 [1] 2.170000e-06 4.230000e-06 1.680000e-06 1.150000e-05 1.720000e-08
 [6] 2.160000e-05 7.920000e-06 1.270000e-05 1.720000e-05 2.400000e-05
[11] 3.110000e-06 1.691535e-07

The vector gear_mult (length 12) is a set of multipliers which are applied to the ‘out of the box’ activity rates (defaults = 1):

model$data$fleet.model$gear_mult
 [1] 1 1 1 1 1 1 1 1 1 1 1 1

Code to build a model scenario

# Example of code to run a baseline case of the North Sea model with 2003-2013 conditions, 
# and then edit the model list object to create a scenario run with beam trawl activity
# reset to zero
#--------------------------------------------------------
# Read the embedded North Sea 2003-2013 model and assign an identifier for the results
base_model<-e2e_read("North_Sea", "2003-2013",model.ident="base",silent=TRUE)
# Here model.ident assigns a text string to add to any output csv files
#
# Run the model for 3 years and save the results to a named list object
base_results<-e2e_run(base_model,nyears=3)
# You could visualize the output from the run but we know this shows a repeating annual cycle with no trend)
# e2e_plot_ts(base_model,base_results,"ECO")

# copy the baseline model list object to a new model list
scenario_model<-base_model
#
# Reset the activity multiplier of demersal beam trawlers to zero
# Beam trawlers are element 4 of the fleet list
scenario_model$data$fleet.model$gear_mult[4] <- 0
# Assign a new identifier for any .csv outputs
scenario_model$setup$model.ident <- "no_beam"
#
# Run the model for 10 years and save the results to a new list object
scenario_results<-e2e_run(scenario_model,nyears=10)
# Visualize the output from the run (should show trends in some outputs due to the change in gear activity)
e2e_plot_ts(scenario_model,scenario_results,"ECO")

Figure 12. Scenario model run for 2003-2013 in the North Sea with no beam trawlers

It’s hard to see the impacts of the scenario in these basic time series plots, but we can compare the final year annual averaged mass of each food web component using the function e2e_compare_runs_bar() to draw a tornado plot:

mdiff_results <- e2e_compare_runs_bar(selection="AAM",
            model1=NA, use.saved1=FALSE, results1=base_results,
            model2=NA, use.saved2=FALSE, results2=scenario_results,
            log.pc="PC", zone="W",
            bpmin=(-25),bpmax=(+25),
            maintitle="Effect of eliminating beam trawlers")
[1] "Using baseline data held in memory from an existing model run"
[1] "Using scenario data held in memory from an existing model run"

# The dataframe mdiff_results contains the data used in the plot so you can make your own graph if you prefer

Figure 13. Tornado plot of the percentage different in annual avareged mass of model components between a baseline 2003-2013 model and a scenario in which beam trawling is eliminated. Green bars to the right indicate higher in the scenarion case, red bars to the left indicate less.

Check help(e2e_compare_runs_bar) for options available with this function.

EXERCISE 5. Try of a small coding project of your choice using the examples above as a template.

Ideas:

Reset seabed ploughing rate of all or selected gears to zero
Reset pelagic fishing to more or less than the baseline
Apply a feeding rate penalty to a selected guild - to represent the impact of litter

Hints: the model object contains some ‘multiplier factors’ to make it a bit easier to create scenarios. These are:

model$data$fleet.model$gear_mult - a vector of 12 values - multipler applied to gear avtivity rates (default = 1)
model$data$fleet.model$HRscale_vector_multiplier - a vector of 12 values - multipler applied to fishing mortlaity rates for each guild (default = 1; order of elements as in model$data$fleet.model$HRscale_vector
model$data$biology.scenarios$uptake.multiplier - 16 x 2 dataframe - multiplier applied to the maximum uptake rates of each guild in the inshore and offshore zones (default = 1) - represents a non-lethal effect eg plastic ingestion
model$data$biology.scenarios$additional.mortality - 16 x 2 dataframe - additional mortality (% per year) applied to each guild in the inshore and offshore zones (default = 0) - represents a lethal effect eg wind turbine collisions

Fishery yield curves - a special type of scenario analysis

Fishery yield curves show the steady state response of biomass and catch of a guild to fishing mortality (harvest ratio).

The catch response should be dome-shaped with the peak representing Maximum Sustainable Yield (MSY). The fishing mortality at MSY is called F_MSY.

It is relatively straightforward to write a script to loop through a sequence of harvest ratio multipliers, run the model and save the relevant results:

Set up a dataframe to hold the results
Set up a loop to cycle through multipliers of harvest ratio (ie fishing mortality rate)
Each loop, run the model to a steady state (minimum 20 years) for the current value of the harvest ratio multiplier
Each loop, accumulate the annual average masses of each food web component and the annual landings and discards in the results dataframe
On completion, plot the annual average mass and catch of planktivorous or demersal fish against the harvest ratio multiplier

model<-e2e_read("North_Sea","1970-1999")

PFHRvector <- seq(0,3,by=0.5)  # set up the vector of harvest ratios
Nsteps     <- length(PFHRvector)

#setup data store to hold the results
datastore<-data.frame(PFHR<-PFHRvector,PlankFishbiom<-rep(0,Nsteps),PlankFishland<-rep(0,Nsteps),PlankFishdisc<-rep(0,Nsteps))

nyr<-20  # number of years to run each harvest ratio - needs to be enough to achieve steady state

for(kkk in 1:(Nsteps)){

message("Run ",kkk,"; pf harvest ratio multiplier = ",PFHRvector[kkk])

model$data$fleet.model$HRscale_vector_multiplier[1]<-PFHRvector[kkk]

results<-e2e_run(model,nyears=nyr,csv.output=FALSE)

#Identify the rows of data to be pulled out of the relevent output object
PF_browid<-which(results$final.year.outputs$mass_results_wholedomain$Description=="Planktivorous_fish")
PF_lrowid<-which(results$final.year.outputs$annual_flux_results_wholedomain$Description=="Plank.fish_landings_live_weight")
PF_drowid<-which(results$final.year.outputs$annual_flux_results_wholedomain$Description=="Plank.fish_discards")

datastore$PlankFishbiom[kkk]<-results$final.year.outputs$mass_results_wholedomain[PF_browid,1]
datastore$PlankFishland[kkk]<-results$final.year.outputs$annual_flux_results_wholedomain[PF_lrowid,1]
datastore$PlankFishdisc[kkk]<-results$final.year.outputs$annual_flux_results_wholedomain[PF_drowid,1]

}

#Plot the results
par(mfrow=c(2,1))
plot(datastore$PFHR,datastore$PlankFishbiom,type="l",
      xlab="Harvest ratio mutiplier",ylab="Planktivorous fish biomass")
plot(datastore$PFHR,datastore$PlankFishland+datastore$PlankFishdisc,type="l",
      xlab="Harvest ratio mutiplier",ylab="Planktivorous fish catch")

This procedure is coded into the functions e2e_run_ycurve() and e2e_plot_ycurve() with some additional details. Check help(e2e_plot_ycurve). Here is a visualization of example data generated with e2e_run_ycurve():

# Plot example data for one of the North Sea model versions internal to the package
model <- e2e_read("North_Sea", "1970-1999",silent=TRUE)
pf_yield_data<-e2e_plot_ycurve(model, selection="PLANKTIV", use.example=TRUE)
Reading example results from StrathE2E2examples data package for the North_Sea 1970-1999 model

Figure 14. Yield data plot for planktivoirous fish in the 1970-1999 North Sea model with baseline demersal harvest ratios. Upper panel, annual average biomass of planktivorous fish (mMN.m^-2 in the whole modle domain) as a function of harvest ratio. Lower panel: catch divided into landings and discards of planktivorous fish (mMN.m^-2.y^-1) as a function of harvest ratio.

In the above example the dataframe pf_yield_data holds the masses and catches of all the other components of the food web for each increment of harvest ratio multiplier. Below is some code to visualize these data and reveal the indirect food web connectivity impacts of the gradient of planktivorous fish harvesting rates.

#First have a look at the dataframe pf_yield_data
head(pf_yield_data)
  PlankFishHRmult DemFishHRmult PlankFishbiom PlankFishland PlankFishdisc
1             0.0             1      49.60182      0.000000   0.000000000
2             0.2             1      45.53717      2.350435   0.007812564
3             0.4             1      40.46815      4.175492   0.013881520
4             0.5             1      37.56036      4.842612   0.016101612
5             0.6             1      34.50485      5.336611   0.017745678
6             0.7             1      31.40302      5.665270   0.018837143
  DemFishbiom DemFishland DemFishdisc MigFishbiom MigFishland MigFishdisc
1   0.5429663  0.08745704  0.06770612    1.854822   0.7798520   0.1897934
2   1.4744002  0.27704674  0.20522092    1.855700   0.7802267   0.1898805
3   3.5514836  0.72413974  0.48864064    1.854904   0.7798640   0.1897890
4   4.7243989  0.99029530  0.63566380    1.854077   0.7794935   0.1896979
5   5.7373503  1.22757005  0.75568829    1.853198   0.7791021   0.1896025
6   6.5250663  1.41658564  0.84484067    1.852407   0.7787506   0.1895179
  Benthsuspdepbiom Benthsuspdepland Benthsuspdepdisc Benthcarnscavbiom
1         50.38394       0.09458175     0.0005483064          9.770246
2         50.27922       0.09417077     0.0005459316          9.050089
3         50.17636       0.09385823     0.0005441198          8.672532
4         50.14555       0.09380304     0.0005437963          8.687213
5         50.13035       0.09381763     0.0005438747          8.857696
6         50.12300       0.09387306     0.0005441875          9.155849
  Benthcarnscavland Benthcarnscavdisc CarnivZoobiom CarnivZooland CarnivZoodisc
1        0.05070891        0.02912916      2.812613             0             0
2        0.04569241        0.02558717      3.390310             0             0
3        0.04294070        0.02356553      4.395327             0             0
4        0.04293470        0.02349139      5.131129             0             0
5        0.04401586        0.02418609      6.018549             0             0
6        0.04601601        0.02554877      7.007124             0             0
     Birdbiom Birdland     Birddisc Pinnipedbiom Pinnipedland Pinnipeddisc
1 0.012773516        0 1.691732e-06  0.004178594            0 6.719167e-06
2 0.011901478        0 1.571796e-06  0.006705981            0 1.365352e-05
3 0.011125727        0 1.466808e-06  0.013439402            0 3.078911e-05
4 0.010669995        0 1.405013e-06  0.016852473            0 3.961488e-05
5 0.010110695        0 1.328970e-06  0.019257713            0 4.624084e-05
6 0.009446687        0 1.238947e-06  0.020515792            0 5.033734e-05
  Cetaceanbiom Cetaceanland Cetaceandisc Macrophytebiom Macrophyteland
1   0.08222495 1.573626e-04 0.0015347621       19.80999              0
2   0.07604189 1.406299e-04 0.0013715675       19.84066              0
3   0.06946032 1.227186e-04 0.0011968782       19.84287              0
4   0.06554912 1.120964e-04 0.0010932799       19.83227              0
5   0.06104986 9.988807e-05 0.0009742112       19.81608              0
6   0.05605551 8.630867e-05 0.0008417710       19.79696              0
  Macrophytedisc Phytoplanktonbiom OmnivZoobiom
1              0          6.779332     31.83528
2              0          6.681376     32.63153
3              0          6.639763     32.97345
4              0          6.649560     32.88339
5              0          6.678847     32.61636
6              0          6.721173     32.22405

#Then make a plot of the dependence of cetaceans on planktivorous fish harvesting rate
par(mfrow=c(2,1))
par(mar=c(3.2,5,2,0.8))
ym<-1.1*max(pf_yield_data$Cetaceanbiom)
plot(pf_yield_data$PlankFishHRmult,pf_yield_data$Cetaceanbiom,ylim=c(0,ym),type="l",
     lwd=3,yaxt="n",xaxt="n",ann=FALSE)
   abline(v=1,lty="dashed")
   axis(side=1,las=1,cex.axis=0.9)
   axis(side=2,las=1,cex.axis=0.9)
   mtext("Planktiv. fish harvest ratio multiplier",cex=1,side=1,line=2)
   mtext("Cetacean biomass",cex=1,side=2,line=3.5)
   mtext(bquote("mMN.m"^-2),cex=0.7,side=3,line=-0.05,adj=-0.18)
ym<-1.1*max(pf_yield_data$Cetaceandisc)
plot(pf_yield_data$PlankFishHRmult,pf_yield_data$Cetaceandisc,ylim=c(0,ym),type="l",
     lwd=3,yaxt="n",xaxt="n",ann=FALSE)
   abline(v=1,lty="dashed")
   axis(side=1,las=1,cex.axis=0.9)
   axis(side=2,las=1,cex.axis=0.9)
   mtext("Planktiv. fish harvest ratio multiplier",cex=1,side=1,line=2)
   mtext("Cetacean by-catch",cex=1,side=2,line=3.5)
   mtext(bquote("mMN.m"^-2 ~ ".y"^-1),cex=0.7,side=3,line=-0.05,adj=-0.18)

Figure 15. Example of a possible user-generated yield data plot for cetaceans in the 1970-1999 North Sea model as a function of planktivorous harvest ratio, with baseline demersal harvest ratios. Upper panel, annual average biomass of cetaceans (mMN.m^-2 in the whole modle domain) as a function of planktivorous fish harvest ratio. Lower panel: catch divided into targeted landings and by-catch (discards) of cetaceans (mMN.m^-2.y^-1) as a function of harvest ratio.

To repeat for demersal fish, replace selection=“PLANKTIV” with selection=“DEMERSAL”

Brief tour of optimization, sensitivity analysis and MCMC functions

See the package website for full details of parameter optimization, sensitivity and Monte Carlo methodology

Parameter optimization

Parameter optimization is a key part of developing a model implementation and assessing its credibility. The package includes a simulated annealing scheme for maximizing the likelihood of a set of empirical data on the state of the ecosystem given proposed parameters:

e2e_optimize_eco() - Maximize the likelihood of observed ecosystem data by optimizing ecology model parameters
e2e_optimize_hr() - Maximize the likelihood of observed ecosystem data by optimizing harvest ratio scaling parameters
e2e_optimize_act() - Maximize the likelihood of observed ecosystem data or known harvest ratios by optimizing fishing activity parameters

The help pages for these functions provide a full explanation of the control options and example code for running them.

For the e2e_optimize_eco() process, details of the empirical target data to which the parameters are optimized is held in the model results object at resultsobject$final.year.outputs$annual.target.data, (where “resultsobject” is the name assigned to a given results object).

The correspondence between a given model results output and the target data is visualised by:

# Load the 1970-1999 version of the North Sea model supplied with the package, run, and 
# generate a plot. Note that these are the observational data that were used as the target 
# for optimizing the model parameters
model <- e2e_read("North_Sea", "1970-1999",silent=TRUE)
results <- e2e_run(model, nyears=3)
e2e_compare_obs(selection="ANNUAL", model, results=results)

Figure 16. Comparison of 1970-1999 observational data from the North Sea (black box-plots) with corresponding values derived from a single run of the model (red symbols). For the observations, the box spans 50% of the distribution of observations and the whiskers 99%. The median is shown by a black tick mark.

Sensitivity analysis

Sensitivity analysis helps to highlight parameters which are:

the most influence in the simulated annealing optimization process
the most sensitive for any particular guild or process in the model

The e2e_run_sens() function performs a global analysis of parameter sensitivity using the Morris Method (Morris, 1991). This involves a factorial sampling scheme to generate parameter values for replicate model runs and generate values for the ‘Elementary Effect’ (EE) of each parameter. The mean and standard deviation of EE values over many runs allows ranking of the parameters in terms of the sensitivity of their effects on the model, and the strength of interactions with other parameters. For a review of global sensitivity analysis methods, including the Morris Method, see Wu et al. (2013).

Morris, M.D. (1991). Factorial sampling plans for preliminary computational experiments. Technometrics, 33, 161-174.
Wu, J., Dhingra, R., Gambhir, M. & Remais, J.V. (2013). Sensitivity analysis of infectious disease models: methods, advances and their application. Journal of the Royal Society Interface, 10: 20121018, 14pp.

The sensitivity analysis involves many thousands of model runs taking several days on a single processor. Procedures for parallel processing of are included in the package. Results for the North Sea model are included in the examples data package.

Sensitivity analysis results are visualized with the function e2e_plot_sens_mc()

# Plot the results of a parameter sensitivity analysis from the examples data package.
model<-e2e_read("North_Sea","1970-1999",silent=TRUE)
e2e_plot_sens_mc(model, selection="SENS", use.example=TRUE)  
Reading example results from StrathE2E2examples data package for the North_Sea 1970-1999 model

Figure 17. Sensitivity of the likelihood of observed data to each of the parameters for the 1970-1999 North Sea model from the examples data package. Each symbol in the plot represents a single parameter in the model. The x-axis represents the sensitivity of the overall likelihood of the observed data to each parameter. The y-axis represents the extent to which each parameter interacts with others. The parameters are colour-coded to indicate 6 different types - fitted and parameters of the ecology model, fishing fleet model parameters, fishery harvest ratios, environmental drivers, and physical configuration parameters. The wedge formed by the two dashed lines corresponds to $\pm$ 2 standard errors of the mean, so for points falling outside of the wedge there is a significant expectation that the distribution of elementary effects is non-zero. For points falling within the wedge the distribution of elementary effects is not significantly different from zero.

You can examine the list of individual model parameters from the example dataset in descending order of Elementary Effect using:

head(StrathE2E2examples::example.results$North_Sea$variant_1970.1999$sorted_parameter_elementary_effects[,c(5,6,7,9,10,11)])

Monte Carlo analysis

The Monte Carlo scheme provided with the StrathE2E2 package (e2e_run_mc()) generates an ensemble of model runs by sampling the ecology model parameters from uniform distibutions centred on an initial (baseline) parameter set. Ideally, this should be the maximum likelihood parameter set arrived at by prior application the various optimization functions in the package to ‘fit’ the model to a suite of observational data on the state of ecosystem. Each iteration in the Monte Carlo scheme then generates a unique time series of model outputs at daily intervals, together with the overall likelihood of the obervational target data. From these data, the function then derives credible intervals of each output by weighting the values from each parameter set by the associated likelihood.

Running a Monte Carlo analysis is time consuming so the examples data package contains results for the North Sea model. Most of the results visualization functions that we explored earlier can also plot credible intervals of model output from these MC outputs. A few of examples are shown here:

# Plot credible intervals of annual trophic indices derived from a Monte Carlo
# analysis for the 1970-1999 North Sea model:
model <- e2e_read("North_Sea", "1970-1999",silent=TRUE)
e2e_plot_trophic(model,ci.data=TRUE,use.example=TRUE)
Reading example results from StrathE2E2examples data package for the North_Sea 1970-1999 model

Figure 18. Credible intervals of annual mean trophic level and omnivory indices for the final year the 1970-1999 North Sea model. The box spans 50% of the distribution of credible outputs and the whiskers 99%. The median is shown by a black tick mark. The maximum likelihood value is indicated by a red tick mark.

# Plot credible intervals around the annual cycle of nutreint and phytoplankton
# outputs for the 1970-1999 North Sea model:
model <- e2e_read("North_Sea", "1970-1999",silent=TRUE)
e2e_plot_eco(model, selection="NUT_PHYT",ci.data=TRUE,use.example=TRUE)
Reading example results from StrathE2E2examples data package for the North_Sea 1970-1999 model
PlotOffshore zone upper layerandSuspended bact & detritus
PlotOffshore zone upper layerandAmmonia
PlotOffshore zone upper layerandNitrate
PlotOffshore zone upper layerandPhytoplankton
PlotOffshore zone lower layerandSuspended bact & detritus
PlotOffshore zone lower layerandAmmonia
PlotOffshore zone lower layerandNitrate
PlotOffshore zone lower layerandPhytoplankton
PlotInshore zoneandSuspended bact & detritus
PlotInshore zoneand Ammonia
PlotInshore zoneandNitrate
PlotInshore zoneandPhytoplankton

Figure 19. Credible intervals of daily interval nutrient, detritus and phytoplankton outputs for the final year the 1970-1999 North Sea model. The grey shading spans 50% of the distribution of credible outputs and the dashed lines 99%. The median is shown by a black line. The maximum likelihood value is indicated by a red line.

The full list of functions able to display either single run or credible intervals of model outputs is:

e2e_compare_obs()
e2e_compare_runs_box()
e2e_plot_eco()
e2e_plot_migration()
e2e_plot_trophic()
e2e_plot_biomass()

EXERCISE 6. Spend a few minutes exploring the credible interval results from the examples package

For example: use the e2e_plot_eco() function to plot credible intervals of different food web components

Installing and using other model implementations

If you want to use the package to work with Case Study model implementations other than the North Sea you wil have to set up a folder on your local disk and then point e2e_read() to read from there.

Using your file manager, make a new folder - eg called “SE2E_workspace” or whatever
Make a sub-folder called “Models”, and a sub-folder called “Results”

In a web-browser:

Navigate to: https://www.marineresourcemodelling.maths.strath.ac.uk/strathe2e/articles/Implementations.html

Click on the model region of your choice and if available download the implentation zipfile into your newly created “Models” folder
Extract the zipfile. The resulting path should look like e.g. SE2E_workspace/Models/modelname/modelvariant/ Be careful when extracting not to introduce and additional level at “Models”.

In R

“Change dir” to make “SE2E_workspace” or whatever you called it, the home directory for R
Read your new model (eg Celtic Sea):

model<-e2e_read(model.name="Celtic-Sea",
        model.variant="2003-2013",
        models.path="Models",
        results.path="Results",
        model.ident="test")

You can also copy the North Sea model from within the package unto your new Models folder using eg :

e2e_copy("North_Sea", "2003-2013",dest.path="Models")
# Note that the "Models" address is relative to the current home folder for R.

EXERCISE 7. Download a Case Study model into a new workspace and get it running. Please bear in mind that some of them are working drafts.

That’s the end of today’s training session on using the StrathE2E R package. Hope you found it useful. Please contact us for help and to report any bugs that you find when you start using it.

StrathE2E2 (v 4.0.1) tutorial

Michael Heath (m.heath@strath.ac.uk)
Jack Laverick (jack.laverick@strath.ac.uk)
Douglas Speirs (d.c.speirs@strath.ac.uk)

University of Strathclyde, Department of Mathematics and Statistics,
Glasgow, Scotland, UK

September 08, 2025

Overview

Requirements

Installing the package

Installing the supplementary data-package

Initial tour of the package

EXERCISE 1. Spend a few minutes roaming around the package help, starting from help(StrathE2E2)

Let’s get started

Loading the package

Load a model

Viewing the model input data

EXERCISE 2. Read the 2003-2013 variant of the North Sea model and have a go at plotting the various aspects the input data.

Run a model and plot basic results

EXERCISE 3. Run the 2003-2013 variant of the North Sea model and plot the basic results.

More detailed visualization of model outputs

Final year time series of ecological data (e2e_plot_eco())

Plotting distributions of catches by gear and guild e2e_plot_catch()

Zonal distributions of annual average biomass e2e_plot_biomass()

Final year time series of migration fluxes (e2e_plot_migration())

Plotting mean trophic level and omnivory indices e2e_plot_trophic()

EXERCISE 4. Spend a few minutes exploring the final year plotting functions with output from a North Sea model run. Use help(function name) for information and assistance.

Writing code to create model scenarios

EXERCISE 5. Try of a small coding project of your choice using the examples above as a template.

Fishery yield curves - a special type of scenario analysis

Brief tour of optimization, sensitivity analysis and MCMC functions

Parameter optimization

Sensitivity analysis

Monte Carlo analysis

EXERCISE 6. Spend a few minutes exploring the credible interval results from the examples package

Installing and using other model implementations

EXERCISE 7. Download a Case Study model into a new workspace and get it running. Please bear in mind that some of them are working drafts.

Data	Units	Comments
Sea surface irradiance	\(\mu\)E.m^-2.d^-1
Suspended particulate matter	g.m^-3
Temperature	deg-C
Vertical diffusivity gradient	m.d^-1	derived from the vertical diffusivity (m².s^-1) and mixing length scale (m)
External inflows	m³ per m² sea surface of model domain	derived from proportion input per layer volume, layer thicknesses and areas
River discharge	m³ per m² sea surface of model domain	derived from proportion input to inshore volume, and inshore layer thickness and area
Inshore significant wave height	m
Proportion of seabed disturbed	d^-1	aggregated over the three sediment classes in each zone
External boundary nitrate concentration	mMN.m^-3
External boundary ammonia concentration	mMN.m^-3
External boundary phytoplankton concentration	mMN.m^-3
External boundary detritus concentration	mMN.m^-3
River nitrate concentration	mMN.m^-3
River ammonia concentration	mMN.m^-3
Atmospheric nitrate deposition flux	mMN.m^-2.d^-1	area averaged flux density over each zone
Atmospheric ammonia deposition flux	mMN.m^-2.d^-1	area averaged flux density over each zone

StrathE2E2 (v 4.0.1) tutorial

Michael Heath (m.heath@strath.ac.uk) Jack Laverick (jack.laverick@strath.ac.uk) Douglas Speirs (d.c.speirs@strath.ac.uk) University of Strathclyde, Department of Mathematics and Statistics, Glasgow, Scotland, UK

September 08, 2025

Overview

Requirements

Installing the package

Installing the supplementary data-package

Initial tour of the package

EXERCISE 1. Spend a few minutes roaming around the package help, starting from help(StrathE2E2)

Let’s get started

Loading the package

Load a model

Viewing the model input data

EXERCISE 2. Read the 2003-2013 variant of the North Sea model and have a go at plotting the various aspects the input data.

Run a model and plot basic results

EXERCISE 3. Run the 2003-2013 variant of the North Sea model and plot the basic results.

More detailed visualization of model outputs

Final year time series of ecological data (e2e_plot_eco())

Plotting distributions of catches by gear and guild e2e_plot_catch()

Zonal distributions of annual average biomass e2e_plot_biomass()

Final year time series of migration fluxes (e2e_plot_migration())

Plotting mean trophic level and omnivory indices e2e_plot_trophic()

EXERCISE 4. Spend a few minutes exploring the final year plotting functions with output from a North Sea model run. Use help(function name) for information and assistance.

Writing code to create model scenarios

EXERCISE 5. Try of a small coding project of your choice using the examples above as a template.

Fishery yield curves - a special type of scenario analysis

Brief tour of optimization, sensitivity analysis and MCMC functions

Parameter optimization

Sensitivity analysis

Monte Carlo analysis

EXERCISE 6. Spend a few minutes exploring the credible interval results from the examples package

Installing and using other model implementations

EXERCISE 7. Download a Case Study model into a new workspace and get it running. Please bear in mind that some of them are working drafts.

Michael Heath (m.heath@strath.ac.uk)
Jack Laverick (jack.laverick@strath.ac.uk)
Douglas Speirs (d.c.speirs@strath.ac.uk)

University of Strathclyde, Department of Mathematics and Statistics,
Glasgow, Scotland, UK