\documentclass[11pt]{article} \usepackage[top=2cm, bottom=3cm, left=2cm, right=2cm]{geometry} % Page margins \usepackage[utf8]{inputenc} \usepackage{amsmath} % /eqref \usepackage[authoryear,round]{natbib} \usepackage{booktabs} % Some macros to improve tables \usepackage{amsfonts} % /checkmark \usepackage{url} \usepackage[none]{hyphenat} % No hyphens %\VignetteIndexEntry{Cheddar plots and statistics} %\VignetteKeyword{food web} %\VignetteKeyword{body mass} %\VignetteKeyword{numerical abundance} %\VignetteKeyword{community} %\VignetteKeyword{allometry} \newcommand{\code}[1]{\texttt{#1}} \newcommand{\R}{\textsf{R} } \begin{document} \sloppy % Prevent hyphenated words running into margins \title{Cheddar plots and statistics (\Sexpr{packageDescription('cheddar', fields='Version')})} \author{Lawrence Hudson} \date{\Sexpr{packageDescription('cheddar', fields='Date')}} \maketitle \tableofcontents <>= library(cheddar) # Makes copy-paste much less painful options(continue=' ') options(width=90) options(prompt='> ') options(SweaveHooks = list(fig=function() par(mgp=c(2.5,1,0), mar=c(4,4,2,1), oma=c(0,0,1,0), cex.main=0.8))) @ \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.5\textwidth} \section{Introduction} This vignette describes Cheddar's flexible system for plotting properties of ecological communities. You should read the `Community' vignette before reading this one. Cheddar's plotting system is built upon two general-purpose functions: \code{PlotNPS} (for \textbf{Plot} \textbf{N}ode \textbf{P}ropertie\textbf{S}) and \code{PlotTLPS} (\textbf{Plot} \textbf{T}rophic-\textbf{L}ink \textbf{P}ropertie\textbf{S}), which allow a range of properties, first-class and computed, Cheddar-provided and user-defined, to be plotted in a variety of ways. \code{PlotNPS} and \code{PlotTLPS} are used by a large number of convenience `wrapper' functions that provide different views of a community, focussing on food web data that are augmented with estimates of body mass, $M$, and/or numerical abundance, $N$, following \cite{CohenEtAl2003PNAS} (Table \ref{tab:plot_functions}). \begin{table}[h!] \begin{center} \begin{tabular}{cccll} \toprule Trophic links & $M$ & $N$ & Plot function & Description \\ \midrule & \checkmark & & \code{PlotMvRankM} & $\log_{10} (M)$ versus rank$(M)$ \\ & \checkmark & & \code{PlotMDistribution} & Distribution of $\log_{10} (M)$ \\ & & \checkmark & \code{PlotNvRankN} & $\log_{10} (N)$ versus rank$(N)$ \\ & & \checkmark & \code{PlotNDistribution} & Distribution of $\log_{10} (N)$ \\ & \checkmark & \checkmark & \code{PlotBvRankB} & $\log_{10} (B)$ versus rank$(B)$ \\ & \checkmark & \checkmark & \code{PlotBDistribution} & Distribution of $\log_{10} (B)$ \\ & \checkmark & \checkmark & \code{PlotNSpectrum} & Approximate numerical abundance spectrum \\ & \checkmark & \checkmark & \code{PlotNvM} * & $\log_{10} (N)$ versus $\log_{10} (M)$ \\ & \checkmark & \checkmark & \code{PlotBSpectrum} & Approximate biomass abundance spectrum \\ & \checkmark & \checkmark & \code{PlotBvM} * & $\log_{10} (B)$ versus $\log_{10} (M)$ \\ \checkmark & & & \code{PlotPredationMatrix} & Predation matrix \\ \checkmark & & & \code{PlotWebByLevel} & Food web plotted vertically by trophic level \\ \checkmark & & & \code{PlotDegreeDistribution} & Histogram of node degree \\ \checkmark & & & \code{PlotCircularWeb} & An overview of complexity \\ \checkmark & & & \code{PlotWagonWheel} & Concentric circles around a focal node \\ \checkmark & \checkmark & & \code{PlotMRvMC} * & $\log_{10} (M)$ of resource versus $\log_{10} (M)$ of consumer \\ \checkmark & & \checkmark & \code{PlotNRvNC} * & $\log_{10} (N)$ of resource versus $\log_{10} (N)$ of consumer \\ \checkmark & & \checkmark & \code{PlotNPyramid} & $\log_{10} (\sum(N))$ by trophic level \\ \checkmark & \checkmark & \checkmark & \code{PlotBRvBC} * & $\log_{10} (B)$ of resource versus $\log_{10} (B)$ of consumer \\ \checkmark & \checkmark & \checkmark & \code{PlotBPyramid} & $\log_{10} (\sum(B))$ by trophic level \\ \checkmark & \checkmark & \checkmark & \code{PlotAlowervAupper} & Upper link-angle against lower link angle \\ \bottomrule \end{tabular} \caption{Cheddar's high-level plot functions.} \label{tab:plot_functions} \end{center} \end{table} You can use whatever mixture of \code{PlotNPS}, \code{PlotTLPS} and the functions in Table \ref{tab:plot_functions} that meets your needs and fits your data. The examples below use the datasets from the two pelagic epilimnion communities of Tuesday Lake, Michigan, USA \citep{CarpenterAndKitchell1996, CohenEtAl2003PNAS, JonssonEtAl2005AER} and the community of Ythan Estuary in north Scotland \citep{HallAndRaffaelli1991JAnimEcol, EmmersonAndRaffaelli2004JAnimEcol}. \section{Plotting node properties} \subsection{\code{PlotNPS}} \code{PlotNPS} plots one node property against another. Every Cheddar function that plots one point per node delegates the job of plotting to \code{PlotNPS} and all of its power and flexibility are available to all of the functions discussed in this section. \code{PlotNPS} function accepts the names of properties to plot on the x and y axes. These can be either `first-class' properties, such as body mass, $M$, and numerical abundance, $N$, or they can be the names of functions that compute and return node properties. The relationship between species-level log-transformed $N$ and $M$ has implications for community- \citep{ReumanEtAl2008EcolLett} and ecosystem-level \citep{BrownEtAl2004Ecology} theories. We can plot $\log_{10}(N)$ against $\log_{10}(M)$ by using \code{PlotNPS} together with the \code{Log10M} and \code{Log10N} functions. \begin{center} <>= data(TL84) PlotNPS(TL84, 'Log10M', 'Log10N') @ \end{center} \code{PlotNPS} uses the `category' node property, if present, to decide plotting colours and symbols. Producers are shown by green circles, invertebrates by blue squares and vertebrate ectotherms by purple diamonds. Many of Cheddar's plotting an analysis functions make use of the `category' node property; this property is optional but, if included in a community, it should contain one of `producer', `invertebrate', `vert.ecto', `vert.endo' or should be empty. This is explored more in Section \ref{sec:nps_colours and symbols}. \newpage \noindent The \code{TL84} dataset includes trophic links so \code{PlotNPS} will show the food web and highlights cannibals by plotting them in a light-coloured filled circle with a dark border, as shown in the previous example. We can turn off these behaviours. \begin{center} <>= PlotNPS(TL84, 'Log10M', 'Log10N', show.web=FALSE, highlight.nodes=NULL) @ \end{center} We can plots labels rather than points. When plotting labels, nodes given by `highlight.nodes' (by default, those nodes that are cannibals) are shown in a darker colour. \begin{center} <>= PlotNPS(TL84, 'Log10M', 'Log10N', show.nodes.as='labels', show.web=FALSE) @ \end{center} In the above example, labels are taken from the order of the nodes in the \code{TL84} community. We can plot a different label by setting the `node.labels' property. The example below plots node names. \begin{center} <>= PlotNPS(TL84, 'Log10M', 'Log10N', show.nodes.as='labels', show.web=FALSE, node.labels='node', cex=0.5) @ \end{center} The parameter \code{node.labels} also accepts a vector of labels, so we could use letters to label the nodes. \begin{center} <>= lots.of.letters <- c(letters, LETTERS, paste(LETTERS,letters,sep='')) PlotNPS(TL84, 'Log10M', 'Log10N', show.nodes.as='labels', show.web=FALSE, node.labels=lots.of.letters[1:NumberOfNodes(TL84)]) @ \end{center} \code{PlotNPS} can plot both symbols and labels. \begin{center} <>= PlotNPS(TL84, 'Log10M', 'Log10N', show.nodes.as='both', show.web=FALSE, cex=2) @ \end{center} When \code{show.nodes.as} is `both', the size and colours of labels are goverened by \code{label.cex} and \code{label.colour}, the latter discussed further in Section \ref{sec:nps_colours and symbols}. The axes titles in the above examples are taken from the names of the properties. Cheddar provides some functions for producing nicely-formatted labels, shown below. \begin{center} <>= PlotNPS(TL84, 'Log10M', 'Log10N', xlab=Log10MLabel(TL84), ylab=Log10NLabel(TL84)) @ \end{center} \code{PlotNPS} can plot any computed node property. For example, we can examine allometric degree distributions \citep{JonssonEtAl2005AER, OttoEtAl2007Nature, DigelEtAl2011Oikos}, which describe how species' numbers of trophic links scale with their log-transformed body masses. \begin{center} \SweaveOpts{width=8,height=2.7} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(1,3)) PlotNPS(TL84, 'Log10M', 'OutDegree', show.web=FALSE) abline(lm(OutDegree(TL84) ~ Log10M(TL84))) PlotNPS(TL84, 'Log10M', 'InDegree', show.web=FALSE) abline(lm(InDegree(TL84) ~ Log10M(TL84))) PlotNPS(TL84, 'Log10M', 'Degree', show.web=FALSE) abline(lm(Degree(TL84) ~ Log10M(TL84))) @ \end{center} To examine the relationship between trophic structure and body mass \citep{CohenEtAl2003PNAS, JonssonEtAl2005AER, GilljamEtAl2011AER, JacobEtAl2011AER, RiedeEtAl2011EcolLett, deVisserEtAl2011JAnimEcol} we can plot trophic level against log-transformed body mass. The \code{PreyAveragedTrophicLevel} uses the matrix-inversion method of \cite{Levine1980JTheorBiol} to compute trophic levels. \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.5\textwidth} <>= PlotNPS(TL84, 'Log10M', 'PreyAveragedTrophicLevel') @ \end{center} An alternative measure of trophic level is provided by \code{ChainAveragedTrophicLevel}, which computes the mean position of each node in every chain in the food web. \begin{center} <>= PlotNPS(TL84, 'Log10M', 'ChainAveragedTrophicLevel') @ \end{center} Cheddar offers the six different measures of trophic level described by \cite{WilliamsAndMartinez2004AmNat}; see the `Communities' vignette and the help page for \code{TrophicLevels} for details. \newpage \noindent \code{PlotNPS} accepts the range of graphical parameters provided by \R's plotting system, such as \code{ylim} and \code{main}. To directly compare prey-averaged and chain-averaged trophic level we can produce side-by-side plots with the same y-axis limits. \begin{center} \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(1,2)) PlotNPS(TL84, 'Log10M', 'PreyAveragedTrophicLevel', ylim=c(1, 6), main='Prey-averaged') PlotNPS(TL84, 'Log10M', 'ChainAveragedTrophicLevel', ylim=c(1, 6), main='Chain-averaged') @ \end{center} \newpage \noindent It is common to plot some combination of $\log_{10}$-transformed body mass, $M$, numerical abundance, $N$, or biomass, $B$ (e.g.\ \citealp{CohenEtAl2003PNAS, JonssonEtAl2005AER, WoodwardEtAl2005AER}). The convenience functions \code{PlotMvN}, \code{PlotNvM}, \code{PlotBvM} and \code{PlotMvB} are `wrappers' around \code{PlotNPS} that create these plots. \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(2,2)) PlotMvN(TL84) PlotNvM(TL84) PlotBvM(TL84) PlotMvB(TL84) @ \end{center} All of the modifying parameters that can be sent to \code{PlotNPS} to do things like changing the highlighting etc., can also be sent to to these convenience functions. \newpage \subsection{Rank of properties} \code{PlotRankNPS} plots the rank of a property. As with \code{PlotNPS}, properties can be either first-class or computed. To plot $\log_{10}(N)$ against the rank of $N$ \citep{JonssonEtAl2005AER}: \begin{center} <>= PlotRankNPS(TL84, 'Log10N') @ \end{center} To plot $\log_{10}(N)$ against the rank of $M$: \begin{center} <>= PlotRankNPS(TL84, 'Log10N', rank.by='M') @ \end{center} \code{PlotRankNPS} uses \code{PlotNPS} to do the actual plotting so all of \code{PlotNPS}'s power and parameters available. For example, the default behaviour of \code{PlotRankNPS} is to not show the food web. We can use \code{PlotNPS}'s `show.web' parameter to show the food web. \begin{center} <>= PlotRankNPS(TL84, 'Log10N', rank.by='M', show.web=TRUE) @ \end{center} To plot trophic level against rank of $M$ \citep{JonssonEtAl2005AER}: \begin{center} <>= PlotRankNPS(TL84, 'PreyAveragedTrophicLevel', rank.by='M') @ \end{center} \newpage \noindent The \code{PlotRankNPS} function can log-transform the rank values shown on the x axis \citep{JonssonEtAl2005AER}: \begin{center} <>= PlotRankNPS(TL84, 'PreyAveragedTrophicLevel', rank.by='M', log10.rank=TRUE) @ \end{center} The convenience functions \code{PlotMvRankM}, \code{PlotNvRankN} and \code{PlotBvRankB} are `wrappers' around \code{PlotRankNPS} that plot rank $\log_{10}$-transformed body mass, $M$, numerical abundance, $N$, or biomass, $B$ \citep{JonssonEtAl2005AER}. \begin{center} \SweaveOpts{width=8,height=2.7} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(1,3)) PlotMvRankM(TL84) PlotNvRankN(TL84) PlotBvRankB(TL84) @ \end{center} All of the modifying parameters that can be sent to \code{PlotNPS} to do things like changing the highlighting etc., can also be sent to to these convenience functions. \newpage \subsection{Distribution of a property} The \code{PlotNPSDistribution} function plots the distribution of a node property. \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.5\textwidth} <>= PlotNPSDistribution(TL84, 'Log10M') @ \end{center} You can pass arguments to \R's \code{density} function using \code{density.args}, to change the bandwidth, for example: \begin{center} <>= PlotNPSDistribution(TL84, 'Log10M', density.args=list(bw=3)) @ \end{center} The convenience `wrapper' functions \code{PlotMDistribution}, \code{PlotNDistribution} and \code{PlotBDistribution} plot $\log_{10}$-transformed body mass, numerical abundance and biomass abundance respectively. \subsection{Colours and symbols} \label{sec:nps_colours and symbols} Cheddar builds on \R's highly flexible plotting mechanisms. Symbol, outline colour, fill colour and symbol size can controlled either by a first-class property, such as `category' or taxonomic resolution, or by a computed property such as trophic species numbers. As with any \R plot, you can directly set colours, symbols, a graph title etc. If the community being plotted has a node property called `category' the function uses `category' to decide plotting colours and symbols using the specifications given in \code{DefaultCategoryColours} and \code{DefaultCategorySymbols}. \begin{center} <>= PlotNvM(TL84, col=1, pch=19, highlight.nodes=NULL) @ \end{center} \newpage \noindent Plot each node in different colour by providing a vector of colours. \begin{center} <>= PlotNvM(TL84, col=1:56, pch=19, highlight.nodes=NULL) @ \end{center} Plot each level of taxonomic resolution in a different colour. In this example, a value of `colour.spec' is not given so \code{PlotNPS} converts the values of the first-class node property `resolved.to' to a factor and uses the factor levels to colour each node. This is a `quick and dirty' way to set node colours. \begin{center} <>= PlotNvM(TL84, colour.by='resolved.to', pch=19, highlight.nodes=NULL) @ \end{center} \newpage \noindent We can improve on this by showing each level of taxonomic resolution in a specific colour by providing the `colour.spec' parameter. \begin{center} <>= colour.spec <- c(Species='purple3', Genus='green3', 'red3') PlotNvM(TL84, colour.by='resolved.to', colour.spec=colour.spec, pch=19, highlight.nodes=NULL) legend("topright", legend=names(colour.spec), pch=19, col=colour.spec) @ \end{center} The `Unclassified flagellates' node does not have a taxonomic resolution. This node is matched by the unnamed value in \code{colour.spec} and so is shown in red. \newpage \noindent The example below uses the `kingdom' level of classification to determine node colours and symbols. \begin{center} <>= symbol.spec = c(Bacteria=21, Plantae=22, Chromista=23, Protozoa=24, Animalia=25, 19) colour.spec = c(Bacteria='purple3', Plantae='green3', Chromista='blue3', Protozoa='orange3', Animalia='red3', 'black') PlotNvM(TL84, symbol.by='kingdom', symbol.spec=symbol.spec, bg.by='kingdom', bg.spec=colour.spec, colour.by='kingdom', colour.spec=colour.spec, highlight.nodes=NULL) legend("topright", legend=names(colour.spec), pch=symbol.spec, col=colour.spec, pt.bg=colour.spec) @ \end{center} The `Unclassified flagellates' node does not have a classification so is matched by the unnamed values in \code{symbol.spec} and \code{colour.spec} and is shown by a black circle. \newpage \noindent We can fit and plot linear regressions for each kingdom by using the \code{NvMLinearRegressions} and \code{PlotLinearModels} functions. The food web is not shown in this plot to make the regression lines clearer. \begin{center} <>= symbol.spec = c(Bacteria=21, Plantae=22, Chromista=23, Protozoa=24, Animalia=25, 19) colour.spec = c(Bacteria='purple3', Plantae='green3', Chromista='blue3', Protozoa='orange3', Animalia='red3', 'black') PlotNvM(TL84, symbol.by='kingdom', symbol.spec=symbol.spec, bg.by='kingdom', bg.spec=colour.spec, colour.by='kingdom', colour.spec=colour.spec, highlight.nodes=NULL, show.web=FALSE) legend("topright", legend=names(colour.spec), pch=symbol.spec, col=colour.spec, pt.bg=colour.spec) models <- NvMLinearRegressions(TL84, class='kingdom') colours <- PlotLinearModels(models, colour.spec=colour.spec) @ \end{center} The model fitted through all data points is matched by the empty value in \code{colour.spec} so is drawn in black. \code{NvMLinearRegressions} is described in more detail in Section \ref{sec:n_v_m}. \newpage \noindent You can use \code{pch=NA} to turn off symbols and just show the food web. \begin{center} <>= PlotNvM(TL84, pch=NA, highlight.nodes=NULL) @ \end{center} See the help page for the \R \code{par} function for more information about colours and symbols. \newpage \noindent All of Cheddar's plot functions add axis ticks to the top and right of the plot by default. This can be turned off by using the \code{cheddarTopAndRightTicks} option. \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} \begin{center} <>= par(mfrow=c(1,2)) # Don't add ticks options(cheddarTopAndRightTicks=FALSE) PlotNvM(TL84) # Add ticks options(cheddarTopAndRightTicks=TRUE) PlotNvM(TL84) @ \end{center} \newpage \subsection{Highlights and lowlights} \code{PlotNPS} (and by extension all of the other node plotting functions) allow individual nodes to be highlighted or lowlighted using the \code{highlight.nodes} and \code{lowlight.nodes} parameters. These parameters accept functions. Functions should take a single \code{Community} object and return node names. Cheddar provides many such functions, such as \code{Cannibals}, which returns the names of all nodes that consume themselves. \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.5\textwidth} \begin{center} <>= PlotNvM(TL84, highlight.nodes=Cannibals) @ \end{center} To highlight nodes without trophic-links to any other nodes we can use the \code{IsolatedNodes} function. \begin{center} <>= PlotNvM(TL84, highlight.nodes=IsolatedNodes) @ \end{center} These parameters also accept the names or IDs of nodes. We can highlight an individual node. \begin{center} <>= PlotNvM(TL84, highlight.nodes='Chaoborus punctipennis') @ \end{center} \code{PlotNPS} also allows trophic links to be highlighted. The \code{highlight.links} parameter accepts either a \code{data.frame} with the columns `resource' and `consumer' or a function that accepts a single community as its only parameters and returns such a \code{data.frame}. One such function is \code{ResourceLargerThanConsumer}, which returns a \code{data.frame} containing trophic links in which the resource has a larger body mass, $M$, than its consumer. To highlight these links we would simply pass this function as the \code{highlight.links} parameter. \begin{center} <>= PlotNvM(TL84, highlight.links=ResourceLargerThanConsumer) @ \end{center} The example below uses the helper function \code{TrophicLinksForNodes} to highlight all links in to and out of \textit{Chaoborus punctipennis}. \begin{center} <>= PlotNvM(TL84, highlight.nodes='Chaoborus punctipennis', highlight.links=TrophicLinksForNodes(TL84, 'Chaoborus punctipennis')) @ \end{center} \newpage \subsection{Nodes that lack properties} Communities are never completely sampled and some species may be represented by too few individuals for properties such as $M$ and $N$ to be accurately estimated. Other resources such as detritus, commonly eaten by primary consumers in many communities, have no clearly-defined $M$ or $N$. Many datasets will therefore contain food-web nodes that lack one or more properties, for example the Ythan Estuary dataset has a large number of primary consumers that depend upon detritus, which has no $M$ or $N$. \code{PlotNPS}, and hence the higher-level functions that use it, will place nodes without a property at the lowest extent of the relevant axis if the `show.na' parameter is \code{TRUE}. By default, `show.na' is \code{FALSE}. The second of the $\log_{10}(N)$-versus-$\log_{10}(M)$ plots below shows detritus at the bottom-left of the plot. \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} \begin{center} <>= data(YthanEstuary) par(mfrow=c(1,2)) PlotNvM(YthanEstuary) PlotNvM(YthanEstuary, show.na=TRUE) @ \end{center} \newpage \noindent The location of the detritus node is determined by the extent of the data, or by the axes limits, if given. \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.5\textwidth} \begin{center} <>= PlotNvM(YthanEstuary, xlim=c(-10, 4), ylim=c(-10, 13), show.na=TRUE) @ \end{center} The example below shows how \code{PlotNPS} behaves in four different cases - node lacks x property, node lacks y property, node lacks both x and y properties and many nodes lack both x and y properties. \code{lowlight.nodes} defaults to nodes that lack either one of the properties being plotted, so the plots below lowlight the node(s) lacking $M$ and/or $N$. \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(2,2)) np <- NPS(TL84) np[1,'M'] <- NA PlotNvM(Community(nodes=np, trophic.links=TLPS(TL84), properties=CPS(TL84)), main='Node 1 M=NA', show.nodes.as='both', cex=2, show.na=TRUE) np <- NPS(TL84) np[1,'N'] <- NA PlotNvM(Community(nodes=np, trophic.links=TLPS(TL84), properties=CPS(TL84)), main='Node 1 N=NA', show.nodes.as='both', cex=2, show.na=TRUE) np <- NPS(TL84) np[1,'M'] <- NA np[1,'N'] <- NA PlotNvM(Community(nodes=np, trophic.links=TLPS(TL84), properties=CPS(TL84)), main='Node 1 M=NA and N=NA', show.nodes.as='both', cex=2, show.na=TRUE) np <- NPS(TL84) np[c(10, 20, 30, 40),'M'] <- NA np[c(10, 20, 30, 40),'N'] <- NA PlotNvM(Community(nodes=np, trophic.links=TLPS(TL84), properties=CPS(TL84)), main='Nodes 10, 20, 30 and 40 M=NA and N=NA', show.nodes.as='both', cex=2, show.na=TRUE) @ \end{center} \newpage \noindent The \code{PlotMvRankM} function delegates work to \code{PlotNPS}. The behaviour of `show.na' for \code{PlotMvRankM} is shown below. The left-hand plot does not show the detritus node whereas the right-hand plot shows it as a grey circle. \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} \begin{center} <>= par(mfrow=c(1,2)) PlotMvRankM(YthanEstuary) PlotMvRankM(YthanEstuary, show.na=TRUE) @ \end{center} \newpage \section{Plotting trophic-link properties} \subsection{\code{PlotTLPS}} \code{PlotTLPS} (for \textbf{Plot} \textbf{T}rophic-\textbf{L}ink \textbf{P}ropertie\textbf{S}) plots one trophic-link property against another. Every Cheddar function that plots a point per trophic link delegates the job of plotting to \code{PlotTLPS} and all of its power and flexibility are available to all of the functions discussed in this section. \code{PlotTLPS} accepts the names of properties to plot on the x and y axes. These can be both `first-class' and computed properties of trophic links or nodes (where node properties are here viewed as properties of the trophic link of which the node is a part). If a property name begins with either `resource.' or `consumer.' then the remainder of the name is assumed to be a node property. For example, to plot $\log_{10}(M)$ of consumers against $\log_{10}(M)$ of consumers we can use the helper function \code{Log10M}. \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.5\textwidth} \begin{center} <>= PlotTLPS(TL84, 'resource.Log10M', 'consumer.Log10M') @ \end{center} \newpage \noindent We can make the plot more informative if the the limits of both the x and y axes are the same. Setting the `axes.limits.equal' parameter to \code{TRUE} instructs \code{PlotTLPS} to do just that, as shown below. \begin{center} <>= PlotTLPS(TL84, 'resource.Log10M', 'consumer.Log10M', axes.limits.equal=TRUE) @ \end{center} \code{PlotTLPS} uses the \code{TLPS} function to assemble data for plotting. \code{PlotTLPS} uses the \code{category} node property of resources, if present, to decide plotting colours and symbols. Resources that are producers are shown by green circles, invertebrates by blue squares and vertebrate ectotherms by purple diamonds. This is explored more in Section \ref{sec:nps_colours and symbols}. \newpage \noindent The convenience functions \code{PlotPredationMatrix}, \code{PlotMRvMC}, \code{PlotMCvMR}, \code{PlotNRvNC}, \code{PlotNCvNR}, \code{PlotBRvBC}, \code{PlotBCvBR} are `wrappers' around \code{PlotTLPS} that plot a predation matrix (a binary matrix with species shown in node order, starting at the top-left, points on the dashed line are cannibals) or $\log_{10}$-transformed body mass, $M$, numerical abundance, $N$, or biomass, $B$. \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.7\textwidth} <>= par(mfrow=c(2,2)) PlotPredationMatrix(TL84) PlotMRvMC(TL84) PlotNCvNR(TL84) PlotBRvBC(TL84) @ \end{center} \newpage \subsection{Colours and symbols} \label{sec:tlps_colours and symbols} \code{PlotTLPS} uses a similar mechanism to \code{PlotNPS}. If the community being plotted has a node property called `category' the function uses `resource.category' to decide plotting colours and symbols using the specifications given in \code{DefaultCategoryColours} and \code{DefaultCategorySymbols}, as shown by this plot of consumer-versus-resource $\log_{10}$-transormed body mass. \begin{center} <>= PlotMRvMC(TL84) @ \end{center} \newpage \noindent To plot the same data using the consumer's category to decide colours and symbols. \begin{center} <>= PlotMRvMC(TL84, colour.by='consumer.category', bg.by='consumer.category', symbol.by='consumer.category') @ \end{center} The \code{bg.by} parameter is used to determine the background colour of each point. See the section for `bg' in the help page for the \R \code{par} function for more information. \newpage \section{Values by class} Cheddar provides functions that compute values over classes of nodes. For example, we can examine the total $M$, $N$ and biomass in each category. <<>>= SumMByClass(TL84) SumNByClass(TL84) SumBiomassByClass(TL84) @ These functions can show the sums by any property of the community so we could look at the same totals by kingdom. <<>>= SumMByClass(TL84, 'kingdom') SumNByClass(TL84, 'kingdom') SumBiomassByClass(TL84, 'kingdom') @ The `Unclassified flagellates' node does not have a kingdom, and the resulting value is labelled `'. These convenience functions are `wrappers' around the more general-purpose \code{ApplyByClass} function. The following two calls produce the same result. <<>>= SumBiomassByClass(TL84) ApplyByClass(TL84, 'Biomass', 'category', sum) @ See the help page for \code{ApplyByClass} for more examples. \newpage \section{$\log_{10}(N)$-versus-$\log_{10}(M)$ statistics} \label{sec:n_v_m} Cheddar contains some helper functions that assist in the commonly-performed $\log_{10}(N)$-versus-$\log_{10}(M)$ statistics. \code{NvMLinearRegressions} returns a list of linear regressions through $\log_{10}(N)$-versus-$\log_{10}(M)$ node data. By default it fits a regression through all of the nodes (called `all') and a separate regression through each node category. <<>>= models <- NvMLinearRegressions(TL84) names(models) @ You can extract the slopes and intercepts for each of these regression models. <<>>= sapply(models, 'coef') @ The `class' parameter defines the sets of nodes to which \code{NvMLinearRegressions} fits regressions. This defaults to `category', if this is present in the community. You can set the `class' parameter to fit regressions through different sets of nodes, such as through each phylum as shown below. <<>>= models <- NvMLinearRegressions(TL84, class='phylum') names(models) @ The `' label refers to a model fitted to nodes that do not have a phylum - just the `unclassifified flagellates' node in this case. It is not possible to fit a linear regression to a single data point and \code{NvMLinearRegressions} returns \code{NULL} for classes that contain just a single node. The data contain one node in phylum Euglenozoa and one node without a phylum, so these entries are \code{NULL}. <<>>= sapply(models, is.null) @ \code{NvMLinearRegressions} also returns \code{NULL} for classes where all or all bar one of the nodes have $M$ and/or $M$ of \code{NA}. The Broadstone stream dataset contains detrital and algal nodes that all lack both $M$ and $N$; these nodes are in the category. <<>>= data(BroadstoneStream) models <- NvMLinearRegressions(BroadstoneStream) sapply(models, is.null) @ Three helper functions return the slope and/or intercept of an ordinary linear regression through all data. <<>>= NvMSlope(TL84) NvMIntercept(TL84) NvMSlopeAndIntercept(TL84) @ Three more helper functions extract these coefficients for each class. <<>>= NvMSlopeByClass(TL84) NvMInterceptByClass(TL84) NvMSlopeAndInterceptByClass(TL84) @ You can fit regressions to whatever class you like. <<>>= NvMSlopeByClass(TL84, class='kingdom') NvMInterceptByClass(TL84, class='kingdom') NvMSlopeAndInterceptByClass(TL84, class='kingdom') @ \newpage \section{Tri-trophic statistics} \cite{CohenEtAl2009PNAS} proposed a novel analysis that simultaneously considers body mass, $M$, numerical abundance, $N$, and the community's food web. The \code{NvMTriTrophicStatistics} function provides this analysis. It takes a single community as its only parameter. Cannibalistic links and all nodes with an $N$ or $M$ of \code{NA} are removed before statistics are computed. It returns a list containing three matrices: `links', `three.node.chains' and `trophic.chains'. The \code{NvMTriTrophicTable} function returns a table of statistics for a collection of communities. The examples in Section \ref{sec:tri_trophic_examples} show how to use this function and the three matrices that it returns to reconstruct the important figures and tables from \cite{CohenEtAl2009PNAS}. \newpage \section{Quantitative descriptors} \cite{BersierEtAl2002Ecology} described food-web descriptors based upon the magnitude of the trophic interactions between taxa. The examples in Section \ref{sec:quantitative_descriptors_examples} show how to use these functions. \newpage \section{Examples through recreating figures and tables from published articles} The following sections show how to use Cheddar to recreate some published figures. \subsection{\cite{JonssonEtAl2005AER}} \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 3 (p 30)} \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} \begin{center} <>= data(TL86) par(mfrow=c(1,2)) PlotMvN(TL84, show.nodes.as='both', cex=2, xlim=c(-2, 10), ylim=c(-14, 0), highlight.nodes=NULL, highlight.links=NULL, main='') PlotMvN(TL86, show.nodes.as='both', cex=2, xlim=c(-2, 10), ylim=c(-14, 0), highlight.nodes=NULL, highlight.links=NULL, main='') title(main='Jonsson et al. (2005) AER, Fig. 3 (p 30)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 4 (p 33)} \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} \begin{center} <>= par(mfrow=c(1,2)) PlotMCvMR(TL84, xlim=c(-14, 0), ylim=c(-14, 0), main='') abline(a=0, b=1, lty=2) PlotMCvMR(TL86, xlim=c(-14, 0), ylim=c(-14, 0), main='') abline(a=0, b=1, lty=2) title(main='Jonsson et al. (2005) AER, Fig. 4 (p 33)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 5 (p 37)} \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(2,2)) PlotNvM(TL84, xlim=c(-14, 0), ylim=c(-2,10), show.web=FALSE, main='') PlotNvM(TL86, xlim=c(-14, 0), ylim=c(-2,10), show.web=FALSE, main='') PlotBvM(TL84, xlim=c(-14, 0), ylim=c(-8,2), show.web=FALSE, main='') PlotBvM(TL86, xlim=c(-14, 0), ylim=c(-8,2), show.web=FALSE, main='') title(main='Jonsson et al. (2005) AER, Fig. 5 (p 37)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 7 (p 47)} \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.8\textwidth} <>= par(mfrow=c(2,2)) PlotNCvNR(TL84, xlim=c(0, 10), ylim=c(-2,10), main='') abline(a=0, b=1, lty=2) PlotNCvNR(TL86, xlim=c(0, 10), ylim=c(-2,10), main='') abline(a=0, b=1, lty=2) PlotBCvBR(TL84, xlim=c(-8, -2), ylim=c(-8, -2), main='') abline(a=0, b=1, lty=2) PlotBCvBR(TL86, xlim=c(-8, -2), ylim=c(-8, -2), main='') abline(a=0, b=1, lty=2) title(main='Jonsson et al. (2005) AER, Fig. 7 (p 47)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 8 (p 49)} \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.8\textwidth} <>= par(mfrow=c(2,2)) TL84.no.iso <- RemoveIsolatedNodes(TL84) TL86.no.iso <- RemoveIsolatedNodes(TL86) tl84.levels <- floor(TrophicHeight(TL84.no.iso)) tl86.levels <- floor(TrophicHeight(TL86.no.iso)) PlotNPyramid(TL84.no.iso, level=tl84.levels, main='', ylab='Trophic height') PlotNPyramid(TL86.no.iso, level=tl86.levels, main='') PlotBPyramid(TL84.no.iso, level=tl84.levels, main='', ylab='Trophic height') PlotBPyramid(TL86.no.iso, level=tl86.levels, main='') title(main='Jonsson et al. (2005) AER, Fig. 8 (p 49)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 10 (p 57)} \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=\textwidth} <>= par(mfrow=c(2,2)) PlotNvRankN(TL84, xlim=c(0,60), ylim=c(-2, 10), main='') PlotNvRankN(TL86, xlim=c(0,60), ylim=c(-2, 10), main='') PlotBvRankB(TL84, xlim=c(0,60), ylim=c(-8, -2), main='') PlotBvRankB(TL86, xlim=c(0,60), ylim=c(-8, -2), main='') title(main='Jonsson et al. (2005) AER, Fig. 10 (p 57)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 11 (p 60)} \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.8\textwidth} <>= par(mfrow=c(2,2)) PlotRankNPS(TL84, property='Log10N', rank.by='M', log10.rank=TRUE, xlim=c(0,2), ylim=c(-2, 10), ylab=Log10NLabel(TL84), main='') PlotRankNPS(TL86, property='Log10N', rank.by='M', log10.rank=TRUE, xlim=c(0,2), ylim=c(-2, 10), ylab=Log10NLabel(TL84), main='') PlotRankNPS(TL84, property='Log10Biomass', rank.by='M', log10.rank=TRUE, xlim=c(0,2), ylim=c(-8, -2), ylab=Log10BLabel(TL84), main='') PlotRankNPS(TL86, property='Log10Biomass', rank.by='M', log10.rank=TRUE, xlim=c(0,2), ylim=c(-8, -2), ylab=Log10BLabel(TL84), main='') title(main='Jonsson et al. (2005) AER, Fig. 11 (p 60)', outer=TRUE) @ \end{center} \newpage \subsubsection{\cite{JonssonEtAl2005AER}, Fig.\ 12 (p 61)} In this example, we write a small function that plots one community against another. \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} \begin{center} <>= PlotCommunityVCommunity <- function(a, b, property, xlim=NULL, ylim=NULL, ...) { a.nodes <- NP(a, 'node') b.nodes <- NP(b, 'node') all.nodes <- union(a.nodes, b.nodes) a.values <- NPS(a, property)[,property] names(a.values) <- a.nodes b.values <- NPS(b, property)[,property] names(b.values) <- b.nodes points <- PlaceMissingPoints(a.values[all.nodes], xlim, b.values[all.nodes], ylim) plot(points[,1], points[,2], xlim=xlim, ylim=ylim, ...) abline(a=0, b=1, lty=2) } par(mfrow=c(1,2)) PlotCommunityVCommunity(TL84, TL86, 'Log10N', xlim=c(-2,10), ylim=c(-2,10), xlab=~log[10]~(N~of~84), ylab=~log[10]~(N~of~86),pch=19) PlotCommunityVCommunity(TL84, TL86, 'Log10Biomass', xlim=c(-8,-2), ylim=c(-8,-2), xlab=~log[10]~(B~of~84), ylab=~log[10]~(B~of~86),pch=19) title(main='Jonsson et al. (2005) AER, Fig. 12 (p 61)', outer=TRUE) @ \end{center} \newpage \subsection{\cite{LayerEtAl2010AER}} \subsubsection{\cite{LayerEtAl2010AER}, Fig.\ 6 (p 282)} Four of the panels from the published figure. This example uses the `pHWebs' community collection dataset; see the `Collections' vignette for more information about community collections. \begin{center} \SweaveOpts{width=6,height=6} \setkeys{Gin}{width=0.7\textwidth} <>= data(pHWebs) par(mfrow=c(2,2)) for(community in pHWebs[1:2]) { PlotNvM(community, xlim=c(-15, 10), ylim=c(-5,15), main='', highlight.nodes=NULL) text(-15, 13, with(CPS(community), paste(title, ', pH ', pH, sep='')), adj=0, cex=1.5) tlps <- TLPS(community, node.properties='M') tlps <- tlps[!is.na(tlps$resource.M) & !is.na(tlps$consumer.M),] interaction.strength <- log10( (tlps$consumer.M / tlps$resource.M)^0.75 ) plot(density(interaction.strength), xlim=c(-4,14), ylim=c(0,0.6), main='', xlab=~log[10]((M[C]/M[R])^0.75)) rug(interaction.strength) } title(main='Layer et al. (2010) AER, Fig. 6 (p 282)', outer=TRUE) @ \end{center} \newpage \subsection{\cite{WoodwardEtAl2005AER}} \subsubsection{\cite{WoodwardEtAl2005AER}, Fig.\ 4. (p 108)} The left panel shows all trophic links. The right panel shows just links for which the resource has a larger body mass than the consumer. \begin{center} \SweaveOpts{width=8,height=4} \setkeys{Gin}{width=\textwidth} <>= data(BroadstoneStream) par(mfrow=c(1,2)) PlotMvN(BroadstoneStream, show.nodes.as='labels', label.cex=0.8, xlim=c(-2, 4.2), ylim=c(-6,2), main='', show.na=FALSE, highlight.links=NULL) abline(a=0, b=-1) tlps <- TLPS(BroadstoneStream, node.properties='M') lty <- rep(0, NumberOfTrophicLinks(BroadstoneStream)) lty[tlps$resource.M > tlps$consumer.M] <- 1 PlotMvN(BroadstoneStream, show.nodes.as='labels', label.cex=0.8, xlim=c(-2, 4.2), ylim=c(-6,2), main='', show.na=FALSE, highlight.links=NULL, link.lty=lty) abline(a=0, b=-1) title(main='Woodward et al. (2005) AER, Fig. 4 (p 108)', outer=TRUE) @ \end{center} \newpage \subsection{\cite{CohenEtAl2009PNAS}} \label{sec:tri_trophic_examples} \subsubsection{\cite{CohenEtAl2009PNAS}, Table 1 (p 22338)} The \code{NvMTriTrophicStatistics} function computes the tri-trophic statistics for a single community. The \code{NvMTriTrophicTable} function takes a collection of communities and assembles a summary table from the values returned by \code{NvMTriTrophicStatistics}. The function removes nodes lacking $M$ and/or $N$, cannibalistic links and isolated nodes from each community. The \code{data.frame} that it returns includes two sets of four network statistics - both with and without these removals. <<>>= collection <- CommunityCollection(list(TL84, TL86, YthanEstuary)) table <- NvMTriTrophicTable(collection) print(round(table,2)) @ \subsubsection{\cite{CohenEtAl2009PNAS}, Table S3 (p 8)} <<>>= res <- lapply(list(TL84, TL86, YthanEstuary), function(community) { community <- RemoveNodes(community, remove=with(NPS(community), node[is.na(M) | is.na(N)])) community <- RemoveCannibalisticLinks(community) community <- RemoveIsolatedNodes(community) chains <- ThreeNodeChains(community, node.properties='M') MR <- chains$bottom.M MI <- chains$intermediate.M MC <- chains$top.M lp <- TLPS(community, node.properties='M') return (c('MR<=MI<=MC'=sum(MR<=MI & MI<=MC), 'MR<=MCMC'=sum(lp$resource.M>lp$consumer.M), 'All links'=nrow(lp))) }) res <- do.call('cbind', res) colnames(res) <- c('TL84', 'TL86', 'Ythan Estuary') print(round(res,2)) @ \subsubsection{\cite{CohenEtAl2009PNAS}, Fig.\ 1, (p 22336)} \SweaveOpts{width=6,height=9} \setkeys{Gin}{width=6in} \begin{center} <>= par(mfrow=c(3,2)) for(community in list(TL84, TL86, YthanEstuary)) { community <- RemoveIsolatedNodes(community) pch <- rep(1, NumberOfNodes(community)) pch[IsIntermediateNode(community)] <- 20 pch[IsTopLevelNode(community)] <- 8 PlotNvM(community, col=1, highlight.nodes=NULL, show.web=FALSE, main='', pch=pch) PlotAuppervAlower(community, main='') } title(main='Cohen et al. (2009) PNAS, Fig. 1 (p 22336)', outer=TRUE) @ \end{center} \newpage \subsection{\cite{BersierEtAl2002Ecology}} \label{sec:quantitative_descriptors_examples} \subsubsection{\cite{BersierEtAl2002Ecology}, Table 1 (p 2401)} Only the first 6 rows of the table are shown here <<>>= data(ChesapeakeBay) res <- NodeQuantitativeDescriptors(ChesapeakeBay, 'biomass.flow') print(round(res[1:6,],2)) @ \subsubsection{\cite{BersierEtAl2002Ecology}, Table 2 (p 2402)} <<>>= res <- QuantitativeDescriptors(ChesapeakeBay, 'biomass.flow') print(round(res,3)) @ \bibliographystyle{plainnat} \bibliography{cheddar} \end{document}