---
title: "Mesocosm example analyses"
author: "Jessica Melbourne-Thomas, Ben Raymond, Andrew Constable, Simon Wotherspoon"
date: "2013"
output:
rmarkdown::html_vignette:
fig_caption: yes
link-citations: yes
bibliography: QPress.bib
vignette: >
%\VignetteIndexEntry{Mesocosm example analyses}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
---
```{r setup, include = FALSE}
knitr::opts_chunk$set(collapse = TRUE, warning = FALSE, message = FALSE, comment = "#>",
fig.height = 6, fig.width = 6)
figrefs <- c()
figref <- function(z) paste0("Figure ", ifelse(z %in% figrefs, which(figrefs == z), { warning ("unknown figure reference: ", z); "??" }))
tablerefs <- c()
tableref <- function(z) paste0("Table ", ifelse(z %in% tablerefs, which(tablerefs == z), { warning ("unknown table reference: ", z); "??" }))
```
The `QPress` package provides facilities for qualitatively modelling the impact of a press perturbation upon a network model. This document illustrates
the comparison of alternative models in terms of their responses to perturbation.
Two examples are provided for experimental mesocosms in (I) a lake environment
[@Hulot2000], and (II) Antarctica (Davidson et al. unpublished).
# Introduction
An experimental mesocosm is an ecological tool that enables the experimenter to assess assemblage-level responses to environmental change. Mesocosm experiments are performed in the "natural environment" rather than in a laboratory, but in an enclosure that is small enough that key variables can be directly controlled. Because the responses of component species to perturbation can be observed directly in a mesocosm setting, the results from such studies lend themselves to analysis using qualitative network modelling. In particular, these analyses can be used to make inferences regarding the trophic and competitive interactions that drive assemblage-level responses to perturbation in a mesocosm setting and potentially more broadly.
In this document we use the `QPress` package to analyse two sets of models corresponding to two separate mesocosm experiments. Example I is for experimental mesocosms in a lake in France that were subjected to high and low nutrient treatments [@Hulot2000]. Example II examines the effects of elevated concentrations of CO$_2$ on the Southern Ocean microbial loop, studied in environmentally controlled tanks (referred to as "minicosms" by the authors) in Antarctica (Davidson et al. unpublished).
# Example I: Lake mesocosm
`r figrefs <- c(figrefs, "Lake")`
Hulot et al. [-@Hulot2000] describe an eight-variable lake mesocosm model comprising phosphorous, three algal groups, two size-classes of herbivores, invertebrate
carnivores and carnivorous fish (`r figref("Lake")`). Unambiguous, directional responses of model components to increased phosphorous observed in the mesocosm experiment were increases in each of: phosphorous (Phos), large herbivores (Herb2), fish (Carn2) and periphyton (AlgP). Melbourne-Thomas et al. [-@JMT2012] compare the ability of Hulot et al.'s original model (`r figref("Lake")`a) and a variant of the model with one edge added and one edge deleted (`r figref("Lake")`b) to reproduce these changes in response to a positive perturbation to phosphorous. Specifically, the variant model (`r figref("Lake")`b) includes a self-limitation effect for invertebrate carnivores, and the predator-prey edge from fish to invertebrate carnivores has been deleted.
The models shown in `r figref("Lake")` were created using Dia and so are read using `model.dia`. Here we assign the two model objects to `lake.a` and `lake.b`.
```{r init}
library(QPress)
lake.a <- model.dia("Lake-a.dia")
lake.b <- model.dia("Lake-b.dia")
```
`r figrefs <- c(figrefs, "Selector1", "Perturb1")`
We then simulate for model (a) and use `impact.barplot` to view simulation outcomes corresponding with observations from the mesocosm experiment (`r figref("Selector1")`).
```{r sima, eval = FALSE}
simlake.a <- system.simulate(1000, lake.a)
impact.barplot(simlake.a)
```
Doing the same for model (b) gives `r figref("Perturb1")`.
```{r outc1, echo = FALSE, fig.height = 4, fig.width = 8, fig.cap = paste0("**", figref("Perturb1"), ":** Simulation outcomes from the two lake mesocosm models in response to a positive press perturbation to phosphorous. Each plot shows the fraction of positive (orange), and negative (blue) outcomes at each node. The x-axis indicates the number of simulations (out of 1000 total) that match with the selection for \"Monitor\" shown in ", figref("Selector1"), ".")}
simlake.a <- system.simulate(1000, lake.a)
simlake.b <- system.simulate(1000, lake.b)
opar <- par(mfrow = c(1, 2))
impact.barplot0(simlake.a, c(0, 0, 0, 0, 0, 0, 0, 1), c(NA, NA, 1, NA, 1, NA, 1, 1), main = "Model (a)")
impact.barplot0(simlake.b, c(0, 0, 0, 0, 0, 0, 0, 1), c(NA, NA, 1, NA, 1, NA, 1, 1), main = "Model (b)")
par(opar)
```
Approximately 50% of the simulations match with observations for model (a), whereas for model (b) this increases to over 80%. This suggests that model (b) provides a more parsimonious representation of the lake mesocosm system. Melbourne-Thomas et al. [-@JMT2012] provide a direct comparison of the two models using Bayes factors.
# Example II: Antarctic mesocosm
`r figrefs <- c(figrefs, "Antarctic-a")`
This example is based on experiments conducted at Australia's Davis Station in Antarctica, which assessed the responses of a natural community of Antarctic marine microbes from near-shore waters to elevated concentrations of CO$_2$ using environmentally controlled tanks (mesocosms). Davidson et al. (in prep) describe the main components of the Southern Ocean microbial loop represented in these mesocosms and the ecological interactions between them. Specifically:
- Small and large phytoplankton cells consume macronutrients and iron (Fe), although iron is potentially less important for growth of small phytoplankton;
- Small phytoplankton are consumed by heterotrophic nanoflagellates and microzooplankton;
- Large phytoplankton (diatoms) are consumed by microzooplankton;
- Bacteria break down DOC (dissolved organic carbon) and are consumed by heterotrophic nanoflagellates;
- Phytoplankton, microzooplankton and heterotrophic nanoflagellates supply the DOC pool.
These interactions are represented in a model that was built using Dia (black-coloured edges in `r figref("Antarctic-a")`).
When this system was subjected to elevated CO$_2$ concentrations in experimental mesocosms, the following responses were observed:
- Small phytoplankton increased;
- Large phytoplankton decreased;
- Heterotrophic nanoflagellates decreased; and
- Bacteria increased.
In the model shown in `r figref("Antarctic-a")` we assume that the direct effects of CO$_2$ correspond with these observed changes.
We then analyse the response this model to a positive press perturbation to CO$_2$. In this case we enforce a specific ordering of the nodes, and also enforce self-limitation for each model variable.
`r figrefs <- c(figrefs, "Selector2", "Perturb2")`
```{r anta}
labels <- c("CO2", "Fe", "Macronutrients", "Small phytoplankton", "Large phytoplankton",
"Heterotrophic nanoflagellates", "Microzooplankton", "Bacteria", "DOC")
antarctic.a <- model.dia("Antarctic-a.dia", labels = labels)
antarctic.a <- enforce.limitation(antarctic.a)
```
As for Example I, we then simulate and use `impact.barplot` (with the selections shown in `r figref("Selector2")`) to produce `r figref("Perturb2")`.
```{r simanta, eval = FALSE}
simantarctic.a <- system.simulate(1000, antarctic.a)
impact.barplot(simantarctic.a)
```
```{r simanta2, echo = FALSE, fig.height = 4, fig.cap = paste0("**", figref("Perturb2"), "**: Simulation outcomes in response to a positive press perturbation to CO$_2$ for Antarctic mesocosm model (a). Each plot shows the fraction of positive (orange), and negative (blue) outcomes at each node.")}
simantarctic.a <- system.simulate(1000, antarctic.a)
impact.barplot0(simantarctic.a, c(1, 0, 0, 0, 0, 0, 0, 0, 0), c(NA, NA, NA, 1, -1, -1, NA, 1, NA), main = "Model (a)")
```
`r figrefs <- c(figrefs, "Antarctic-b")`
For this version of the model, less than 20% of the simulations match with our "Monitor" criteria, suggesting that it does not provide a particularly good representation of the Antarctic microbial loop system or its response to increased CO$_2$ concentrations. Given that increases in dissolved CO$_2$ are generally expected to have negative impacts of the growth of marine organisms, we could assume a direct negative effect of CO$_2$ on each of large phytoplankton, small phytoplankton, heterotrophic nanoflagellates, and also the pool of available iron (`r figref("Antarctic-b")`).
Reading this alternative representation from a Dia file
```{r modb}
antarctic.b <- model.dia("Antarctic-b.dia", labels = labels)
antarctic.b <- enforce.limitation(antarctic.b)
```
`r figrefs <- c(figrefs, "Perturb3")`
and exploring simulation outcomes with `impact.barplot` (with the same "Monitor" criteria as before) produces `r figref("Perturb3")`.
```{r simantb, eval = FALSE}
simantarctic.b <- system.simulate(1000, antarctic.b)
impact.barplot(simantarctic.b)
```
`r figrefs <- c(figrefs, "Antarctic-c")`
In this case, we still only see approximately 20% of simulations that meet our criteria, suggesting that we need to re-think the structure of our model. Davidson et al. (in prep) indicate that bacteria "transform dissolved organic carbon (DOC)... thereby repackaging carbon and making it available to higher trophic levels". Looking at our original model, we haven't quite captured this in our representation of interactions between bacteria, DOC and macronutrients. In a third version of our microbial loop model (model (c) --- `r figref("Antarctic-c")`) we therefore assume a negative effect of bacteria on DOC and a positive effect of bacteria on macronutrients (which includes carbon).
`r figrefs <- c(figrefs, "Perturb4")`
Reading a third Dia version of the model
```{r modc}
antarctic.c <- model.dia("Antarctic-c.dia", labels = labels)
antarctic.c <- enforce.limitation(antarctic.c)
```
and repeating the simulation and `impact.barplot` steps from above gives `r figref("Perturb4")`.
```{r simantb2, echo = FALSE, fig.height = 4, fig.cap = paste0("**", figref("Perturb3"), "**: Simulation outcomes in response to a positive press perturbation to the CO$_2$ for the alternative microbial loop model (b) shown in ", figref("Antarctic-b"), ".")}
simantarctic.b <- system.simulate(1000, antarctic.b)
impact.barplot0(simantarctic.b, c(1, 0, 0, 0, 0, 0, 0, 0, 0), c(NA, NA, NA, 1, -1, -1, NA, 1, NA), main = "Model (b)")
```
```{r simantc2, echo = FALSE, fig.height = 4, fig.cap = paste0("**", figref("Perturb4"), "**: Simulation outcomes in response to a positive press perturbation to the CO$_2$ for the alternative microbial loop model (c) shown in ", figref("Antarctic-c"), ".")}
simantarctic.c <- system.simulate(1000, antarctic.c)
impact.barplot0(simantarctic.c, c(1, 0, 0, 0, 0, 0, 0, 0, 0), c(NA, NA, NA, 1, -1, -1, NA, 1, NA), main = "Model (c)")
```
## References
