| Type: | Package |
| Title: | A Stable Isotope Mixing Model |
| Version: | 0.5.1.217 |
| Date: | 2024-10-15 |
| URL: | https://github.com/andrewcparnell/simmr,https://andrewcparnell.github.io/simmr/ |
| BugReports: | https://github.com/andrewcparnell/simmr/issues |
| Language: | en-US |
| Description: | Fits Stable Isotope Mixing Models (SIMMs) and is meant as a longer term replacement to the previous widely-used package SIAR. SIMMs are used to infer dietary proportions of organisms consuming various food sources from observations on the stable isotope values taken from the organisms' tissue samples. However SIMMs can also be used in other scenarios, such as in sediment mixing or the composition of fatty acids. The main functions are simmr_load() and simmr_mcmc(). The two vignettes contain a quick start and a full listing of all the features. The methods used are detailed in the papers Parnell et al 2010 <doi:10.1371/journal.pone.0009672>, and Parnell et al 2013 <doi:10.1002/env.2221>. |
| Depends: | R (≥ 3.5.0), R2jags, ggplot2 |
| Imports: | compositions, boot, reshape2, graphics, stats, viridis,bayesplot, checkmate, Rcpp, GGally |
| Encoding: | UTF-8 |
| License: | GPL-2 |GPL-3 [expanded from: GPL (≥ 2)] |
| LazyData: | TRUE |
| Suggests: | knitr, rmarkdown, readxl, testthat, covr, vdiffr, tibble,ggnewscale |
| VignetteBuilder: | knitr |
| NeedsCompilation: | yes |
| Packaged: | 2024-10-16 14:45:08 UTC; emmagovan |
| RoxygenNote: | 7.3.2 |
| Repository: | CRAN |
| Date/Publication: | 2024-10-16 15:10:02 UTC |
| LinkingTo: | Rcpp, RcppArmadillo, RcppDist, |
| Author: | Emma Govan [cre, aut], Andrew Parnell [aut] |
| Maintainer: | Emma Govan <emma.govan.2021@mumail.ie> |
simmr: A package for fitting stable isotope mixing models via JAGS andFFVB in R
Description
This package runs asimple Stable Isotope Mixing Model (SIMM) and is meant as a longer termreplacement to the previous function SIAR.. These are used to infer dietaryproportions of organisms consuming various food sources from observations onthe stable isotope values taken from the organisms' tissue samples. HoweverSIMMs can also be used in other scenarios, such as in sediment mixing or thecomposition of fatty acids. The main functions aresimmr_load,simmr_mcmc, andsimmr_ffvb. The help files contain examples of the use of this package. See also the vignette for a longer walkthrough.
Details
An even longer term replacement for properly running SIMMs is MixSIAR, whichallows for more detailed random effects and the inclusion of covariates.
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
References
Andrew C. Parnell, Donald L. Phillips, Stuart Bearhop, Brice X.Semmens, Eric J. Ward, Jonathan W. Moore, Andrew L. Jackson, Jonathan Grey,David J. Kelly, and Richard Inger. Bayesian stable isotope mixing models.Environmetrics, 24(6):387–399, 2013.
Andrew C Parnell, Richard Inger, Stuart Bearhop, and Andrew L Jackson.Source partitioning using stable isotopes: coping with too much variation.PLoS ONE, 5(3):5, 2010.
Examples
# A first example with 2 tracers (isotopes), 10 observations, and 4 food sourcesdata(geese_data_day1)simmr_in <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_in)# MCMC runsimmr_out <- simmr_mcmc(simmr_in)# Check convergence - values should all be close to 1summary(simmr_out, type = "diagnostics")# Look at outputsummary(simmr_out, type = "statistics")# Look at influence of priorsprior_viz(simmr_out)# Plot outputplot(simmr_out, type = "histogram")Combine the dietary proportions from two food sources after running simmr
Description
This function takes in an object of classsimmr_output and combinestwo of the food sources. It works for single and multiple group data.
Usage
combine_sources( simmr_out, to_combine = NULL, new_source_name = "combined_source")Arguments
simmr_out | An object of class |
to_combine | The names of exactly two sources. These should match thenames given to |
new_source_name | A name to give to the new combined source. |
Details
Often two sources either (1) lie in similar location on the iso-space plot,or (2) are very similar in phylogenetic terms. In case (1) it is common toexperience high (negative) posterior correlations between the sources.Combining them can reduce this correlation and improve precision of theestimates. In case (2) we might wish to determine the joint amount eaten ofthe two sources when combined. This function thus combines two sources aftera run ofsimmr_mcmc orsimmr_ffvb (known asa posteriori combination). The new object can then be called withplot.simmr_input orplot.simmr_output toproduce iso-space plots of summaries of the output after combination.
Value
A newsimmr_output object
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>, Emma Govan
See Also
Seesimmr_mcmc andsimmr_ffvb andthe associated vignette for examples.
Examples
# The datadata(geese_data)# Load into simmrsimmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# Printsimmr_1# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Print itprint(simmr_1_out)# Summarysummary(simmr_1_out)summary(simmr_1_out, type = "diagnostics")summary(simmr_1_out, type = "correlations")summary(simmr_1_out, type = "statistics")ans <- summary(simmr_1_out, type = c("quantiles", "statistics"))# Plotplot(simmr_1_out)plot(simmr_1_out, type = "boxplot")plot(simmr_1_out, type = "histogram")plot(simmr_1_out, type = "density")plot(simmr_1_out, type = "matrix")simmr_out_combine <- combine_sources(simmr_1_out, to_combine = c("U.lactuca", "Enteromorpha"), new_source_name = "U.lac+Ent")plot(simmr_out_combine$input)plot(simmr_out_combine, type = "boxplot", title = "simmr output: combined sources")Compare dietary proportions for a single source across different groups
Description
This function takes in an object of classsimmr_output and createsprobabilistic comparisons for a given source and a set of at least twogroups.
Usage
compare_groups( simmr_out, source_name = simmr_out$input$source_names[1], groups = 1:2, plot = TRUE)Arguments
simmr_out | An object of class |
source_name | The name of a source. This should match the names exactlygiven to |
groups | The integer values of the group numbers to be compared. Atleast two groups must be specified. |
plot | A logical value specifying whether plots should be produced ornot. |
Details
When two groups are specified, the function produces a direct calculation ofthe probability that one group is bigger than the other. When more than twogroups are given, the function produces a set of most likely probabilisticorderings for each combination of groups. The function produces boxplots bydefault and also allows for the storage of the output for further analysis ifrequired.
Value
If there are two groups, a vector containing the differences betweenthe two groups proportions for that source. If there are multiple groups, alist containing the following fields:
Ordering | The differentpossible orderings of the dietary proportions across groups |
out_all | The dietary proportions for this source across the groups specified ascolumns in a matrix |
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
See Also
Seesimmr_mcmc for complete examples.
Examples
## Not run: data(geese_data)simmr_in <- with( geese_data, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means, group = groups ))# Printsimmr_in# Plotplot(simmr_in, group = 1:8, xlab = expression(paste(delta^13, "C (per mille)", sep = "")), ylab = expression(paste(delta^15, "N (per mille)", sep = "")), title = "Isospace plot of Inger et al Geese data")# Run MCMC for each groupsimmr_out <- simmr_ffvb(simmr_in)# Print outputsimmr_out# Summarise outputsummary(simmr_out, type = "quantiles", group = 1)summary(simmr_out, type = "quantiles", group = c(1, 3))summary(simmr_out, type = c("quantiles", "statistics"), group = c(1, 3))# Plot - only a single group allowedplot(simmr_out, type = "boxplot", group = 2, title = "simmr output group 2")plot(simmr_out, type = c("density", "matrix"), grp = 6, title = "simmr output group 6")# Compare groupscompare_groups(simmr_out, source = "Zostera", groups = 1:2)compare_groups(simmr_out, source = "Zostera", groups = 1:3)compare_groups(simmr_out, source = "U.lactuca", groups = c(4:5, 7, 2))## End(Not run)Compare dietary proportions between multiple sources
Description
This function takes in an object of classsimmr_output and createsprobabilistic comparisons between the supplied sources. The group number canalso be specified.
Usage
compare_sources( simmr_out, source_names = simmr_out$input$source_names, group = 1, plot = TRUE)Arguments
simmr_out | An object of class |
source_names | The names of at least two sources. These should matchthe names exactly given to |
group | The integer values of the group numbers to be compared. If notspecified assumes the first or only group |
plot | A logical value specifying whether plots should be produced ornot. |
Details
When two sources are specified, the function produces a direct calculationof the probability that the dietary proportion for one source is bigger thanthe other. When more than two sources are given, the function produces a setof most likely probabilistic orderings for each combination of sources. Thefunction produces boxplots by default and also allows for the storage of theoutput for further analysis if required.
Value
If there are two sources, a vector containing the differencesbetween the two dietary proportion proportions for these two sources. Ifthere are multiple sources, a list containing the following fields:
Ordering | The different possible orderings of the dietary proportionsacross sources |
out_all | The dietary proportions for these sourcesspecified as columns in a matrix |
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
See Also
Seesimmr_mcmc for complete examples.
Examples
data(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# Printsimmr_1# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Print itprint(simmr_1_out)# Summarysummary(simmr_1_out)summary(simmr_1_out, type = "diagnostics")summary(simmr_1_out, type = "correlations")summary(simmr_1_out, type = "statistics")ans <- summary(simmr_1_out, type = c("quantiles", "statistics"))# Plotplot(simmr_1_out, type = "boxplot")plot(simmr_1_out, type = "histogram")plot(simmr_1_out, type = "density")plot(simmr_1_out, type = "matrix")# Compare two sourcescompare_sources(simmr_1_out, source_names = c("Zostera", "Grass"))# Compare multiple sourcescompare_sources(simmr_1_out)Geese stable isotope mixing data set
Description
A real Geese data set with 251 observations on 2 isotopes, with 4 sources, and with corrections/trophic enrichment factors (TEFs or TDFs), and concentration dependence means. Taken from Inger et al (2016). See link for paper
Usage
geese_dataFormat
A list with the following elements
- mixtures
A two column matrix containing delta 13C and delta 15N values respectively
- source_names
A character vector of the food source names
- tracer_names
A character vector of the tracer names (d13C, d15N, d34S)
- source_means
A matrix of source mean values for the tracers in the same order as
mixturesabove- source_sds
A matrix of source sd values for the tracers in the same order as
mixturesabove- correction_means
A matrix of TEFs mean values for the tracers in the same order as
mixturesabove- correction_sds
A matrix of TEFs sd values for the tracers in the same order as
mixturesabove- concentration_means
A matrix of concentration dependence mean values for the tracers in the same order as
mixturesabove
...
@seealsosimmr_mcmc for an example where it is used
Source
<doi:10.1111/j.1365-2656.2006.01142.x>
A smaller version of the Geese stable isotope mixing data set
Description
A real Geese data set with 9 observations on 2 isotopes, with 4 sources, and with corrections/trophic enrichment factors (TEFs or TDFs), and concentration dependence means. Taken from Inger et al (2016). See link for paper
Usage
geese_data_day1Format
A list with the following elements
- mixtures
A two column matrix containing delta 13C and delta 15N values respectively
- source_names
A character vector of the food source names
- tracer_names
A character vector of the tracer names (d13C, d15N, d34S)
- source_means
A matrix of source mean values for the tracers in the same order as
mixturesabove- source_sds
A matrix of source sd values for the tracers in the same order as
mixturesabove- correction_means
A matrix of TEFs mean values for the tracers in the same order as
mixturesabove- correction_sds
A matrix of TEFs sd values for the tracers in the same order as
mixturesabove- concentration_means
A matrix of concentration dependence mean values for the tracers in the same order as
mixturesabove
...
@seealsosimmr_mcmc for an example where it is used
Source
<doi:10.1111/j.1365-2656.2006.01142.x>
Plot thesimmr_input data created fromsimmr_load
Description
This function creates iso-space (AKA tracer-space or delta-space) plots.They are vital in determining whether the data are suitable for running in aSIMM.
Usage
## S3 method for class 'simmr_input'plot( x, tracers = c(1, 2), title = "Tracers plot", xlab = colnames(x$mixtures)[tracers[1]], ylab = colnames(x$mixtures)[tracers[2]], sigmas = 1, group = 1:x$n_groups, mix_name = "Mixtures", ggargs = NULL, colour = TRUE, ...)Arguments
x | An object of class |
tracers | The choice of tracers to plot. If there are more than twotracers, it is recommended to plot every pair of tracers to determinewhether the mixtures lie in the mixing polygon defined by the sources |
title | A title for the graph |
xlab | The x-axis label. By default this is assumed to be delta-13C butcan be made richer if required. See examples below. |
ylab | The y-axis label. By default this is assumed to be delta-15N inper mil but can be changed as with the x-axis label |
sigmas | The number of standard deviations to plot on the sourcevalues. Defaults to 1. |
group | Which groups to plot. Can be a single group or multiple groups |
mix_name | A optional string containing the name of the mixtureobjects, e.g. Geese. |
ggargs | Extra arguments to be included in the ggplot (e.g. axis limits) |
colour | If TRUE (default) creates a plot. If not, puts the plot inblack and white |
... | Not used |
Details
It is desirable to have the vast majority of the mixture observations to beinside the convex hull defined by the food sources. When there are more thantwo tracers (as in one of the examples below) it is recommended to plot allthe different pairs of the food sources. See the vignette for furtherdetails of richer plots.
Value
isospace plot
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
See Also
Seeplot.simmr_output for plotting the output of asimmr run. Seesimmr_mcmc for running a simmr object once theiso-space is deemed acceptable.
Examples
# A simple example with 10 observations, 4 food sources and 2 tracersdata(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)### A more complicated example with 30 obs, 3 tracers and 4 sourcesdata(simmr_data_2)simmr_3 <- with( simmr_data_2, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plot 3 times - first default d13C vs d15Nplot(simmr_3)# Now plot d15N vs d34Splot(simmr_3, tracers = c(2, 3))# and finally d13C vs d34Splot(simmr_3, tracers = c(1, 3))# See vignette('simmr') for fancier x-axis labels# An example with multiple groups - the Geese data from Inger et al 2006data(geese_data)simmr_4 <- with( geese_data, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means, group = groups ))# Printsimmr_4# Plotplot(simmr_4, xlab = expression(paste(delta^13, "C (per mille)", sep = "")), ylab = expression(paste(delta^15, "N (per mille)", sep = "")), title = "Isospace plot of Inger et al Geese data") #'Plot different features of an object created fromsimmr_mcmcorsimmr_ffvb.
Description
This function allows for 4 different types of plots of the simmr outputcreated fromsimmr_mcmc orsimmr_ffvb. Thetypes are: histogram, kernel density plot, matrix plot (most useful) andboxplot. There are some minor customisation options.
Usage
## S3 method for class 'simmr_output'plot( x, type = c("isospace", "histogram", "density", "matrix", "boxplot"), group = 1, binwidth = 0.05, alpha = 0.5, title = if (length(group) == 1) { "simmr output plot" } else { paste("simmr output plot: group", group) }, ggargs = NULL, ...)Arguments
x | An object of class |
type | The type of plot required. Can be one or more of 'histogram','density', 'matrix', or 'boxplot' |
group | Which group(s) to plot. |
binwidth | The width of the bins for the histogram. Defaults to 0.05 |
alpha | The degree of transparency of the plots. Not relevant formatrix plots |
title | The title of the plot. |
ggargs | Extra arguments to be included in the ggplot (e.g. axis limits) |
... | Currently not used |
Details
The matrix plot should form a necessary part of any SIMM analysis since itallows the user to judge which sources are identifiable by the model.Further detail about these plots is provided in the vignette.Some code fromhttps://stackoverflow.com/questions/14711550/is-there-a-way-to-change-the-color-palette-for-ggallyggpairs-using-ggplotaccessed March 2023
Value
one or more of 'histogram', 'density', 'matrix', or 'boxplot'
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>, Emma Govan
See Also
Seesimmr_mcmc andsimmr_ffvb forcreating objects suitable for this function, and many more examples. Seealsosimmr_load for creating simmr objects,plot.simmr_input for creating isospace plots,summary.simmr_output for summarising output.
Examples
# A simple example with 10 observations, 2 tracers and 4 sources# The datadata(geese_data)# Load into simmrsimmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Plotplot(simmr_1_out) # Creates all 4 plotsplot(simmr_1_out, type = "boxplot")plot(simmr_1_out, type = "histogram")plot(simmr_1_out, type = "density")plot(simmr_1_out, type = "matrix")Plot the posterior predictive distribution for a simmr run
Description
This function takes the output fromsimmr_mcmcand plots the posterior predictive distribution to enable visualisation of model fit. The simulated posterior predicted values are returned as part of the object and can be saved for external use
Usage
posterior_predictive(simmr_out, group = 1, prob = 0.5, plot_ppc = TRUE)Arguments
simmr_out | A run of the simmr model from |
group | Which group to run it for (currently only numeric rather than group names) |
prob | The probability interval for the posterior predictives. The default is 0.5 (i.e. 50pc intervals) |
plot_ppc | Whether to create a bayesplot of the posterior predictive or not. |
Value
plot of posterior predictives and simulated values
Examples
data(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# Printsimmr_1# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Prior predictivepost_pred <- posterior_predictive(simmr_1_out)Print simmr input object
Description
Print simmr input object
Usage
## S3 method for class 'simmr_input'print(x, ...)Arguments
x | An object of class |
... | Other arguments (not supported) |
Value
A neat presentation of your simmr object.
Print simmr output object
Description
Print simmr output object
Usage
## S3 method for class 'simmr_output'print(x, ...)Arguments
x | An object of class |
... | Other arguments (not supported) |
Value
Returns a neat summary of the object
See Also
simmr_mcmc andsimmr_ffvb for creatingsimmr_output objects
Plot the prior distribution for a simmr run
Description
This function takes the output fromsimmr_mcmc orsimmr_ffvb and plots the prior distribution to enable visualinspection. This can be used by itself or together withposterior_predictive to visually evaluate the influence ofthe prior on the posterior distribution.
Usage
prior_viz( simmr_out, group = 1, plot = TRUE, include_posterior = TRUE, n_sims = 10000, scales = "free")Arguments
simmr_out | A run of the simmr model from |
group | Which group to run it for (currently only numeric rather than group names) |
plot | Whether to create a density plot of the prior or not. The simulated prior values are returned as part of the object |
include_posterior | Whether to include the posterior distribution on top of the priors. Defaults to TRUE |
n_sims | The number of simulations from the prior distribution |
scales | The type of scale from |
Value
A list containingplot: the ggplot object (useful if requires customisation), andsim: the simulated prior values which can be compared with the posterior densities
Examples
data(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# Printsimmr_1# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Prior predictiveprior <- prior_viz(simmr_1_out)head(prior$p_prior_sim)summary(prior$p_prior_sim)A simple fake stable isotope mixing data set
Description
A simple fake data set with 10 observations on 2 isotopes, with 4 sources, and with corrections/trophic enrichment factors (TEFs or TDFs), and concentration dependence means
Usage
simmr_data_1Format
A list with the following elements
- mixtures
A two column matrix containing delta 13C and delta 15N values respectively
- source_names
A character vector of the food source names
- tracer_names
A character vector of the tracer names (d13C and d15N)
- source_means
A matrix of source mean values for the tracers in the same order as
mixturesabove- source_sds
A matrix of source sd values for the tracers in the same order as
mixturesabove- correction_means
A matrix of TEFs mean values for the tracers in the same order as
mixturesabove- correction_sds
A matrix of TEFs sd values for the tracers in the same order as
mixturesabove- concentration_means
A matrix of concentration dependence mean values for the tracers in the same order as
mixturesabove
...
@seealsosimmr_mcmc for an example where it is used
A 3-isotope fake stable isotope mixing data set
Description
A fake data set with 30 observations on 3 isotopes, with 4 sources, and with corrections/trophic enrichment factors (TEFs or TDFs), and concentration dependence means
Usage
simmr_data_2Format
A list with the following elements
- mixtures
A three column matrix containing delta 13C, delta 15N, and delta 34S values respectively
- source_names
A character vector of the food source names
- tracer_names
A character vector of the tracer names (d13C, d15N, d34S)
- source_means
A matrix of source mean values for the tracers in the same order as
mixturesabove- source_sds
A matrix of source sd values for the tracers in the same order as
mixturesabove- correction_means
A matrix of TEFs mean values for the tracers in the same order as
mixturesabove- correction_sds
A matrix of TEFs sd values for the tracers in the same order as
mixturesabove- concentration_means
A matrix of concentration dependence mean values for the tracers in the same order as
mixturesabove
...
@seealsosimmr_mcmc for an example where it is used
Function to allow informative prior distribution to be included in simmr
Description
The mainsimmr_mcmc function allows for a prior distributionto be set for the dietary proportions. The prior distribution is specifiedby transforming the dietary proportions using the centralised log ratio(CLR). Thesimmr_elicit andsimmr_elicitfunctions allows the user to specifyprior means and standard deviations for each of the dietary proportions, andthen finds CLR-transformed values suitable for input intosimmr_mcmc.
Usage
simmr_elicit( n_sources, proportion_means = rep(1/n_sources, n_sources), proportion_sds = rep(0.1, n_sources), n_sims = 1000)Arguments
n_sources | The number of sources required |
proportion_means | The desired prior proportion means. These should sumto 1. Should be a vector of length |
proportion_sds | The desired prior proportions standard deviations.These have no restricted sum but should be reasonable estimates for aproportion. |
n_sims | The number of simulations for which to run the optimisationroutine. |
Details
The function takes the desired proportion means and standard deviations,and fits an optimised least squares to the means and standard deviations inturn to produced CLR-transformed estimates for use insimmr_mcmc. Using prior information in SIMMs is highlydesirable given the restricted nature of the inference. The priorinformation might come from previous studies, other experiments, or otherobservations of e.g. animal behaviour.
Due to the nature of the restricted space over which the dietary proportionscan span, and the fact that this function uses numerical optimisation, theprocedure will not match the target dietary proportion means and standarddeviations exactly. If this problem is severe, try increasing then_sims value.
Value
A list object with two components
mean | The best estimates ofthe mean to use in |
sd | The best estimates of the standard deviations to use in |
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
Examples
# Data set: 10 observations, 2 tracers, 4 sourcesdata(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Look at the prior influenceprior_viz(simmr_1_out)# Summarysummary(simmr_1_out, "quantiles")# A bit vague:# 2.5% 25% 50% 75% 97.5%# Source A 0.029 0.115 0.203 0.312 0.498# Source B 0.146 0.232 0.284 0.338 0.453# Source C 0.216 0.255 0.275 0.296 0.342# Source D 0.032 0.123 0.205 0.299 0.465# Now suppose I had prior information that:# proportion means = 0.5,0.2,0.2,0.1# proportion sds = 0.08,0.02,0.01,0.02prior <- simmr_elicit(4, c(0.5, 0.2, 0.2, 0.1), c(0.08, 0.02, 0.01, 0.02))simmr_1a_out <- simmr_mcmc(simmr_1, prior_control = list(means = prior$mean, sd = prior$sd, sigma_shape = c(3,3), sigma_rate = c(3/50, 3/50)))#' # Look at the prior influence nowprior_viz(simmr_1a_out)summary(simmr_1a_out, "quantiles")# Much more precise:# 2.5% 25% 50% 75% 97.5%# Source A 0.441 0.494 0.523 0.553 0.610# Source B 0.144 0.173 0.188 0.204 0.236# Source C 0.160 0.183 0.196 0.207 0.228# Source D 0.060 0.079 0.091 0.105 0.135Run asimmr_input object through the Fixed Form VariationalBayes(FFVB) function
Description
This is the main function of simmr. It takes asimmr_input objectcreated viasimmr_load, runs it in fixed formVariational Bayes to determine the dietary proportions, and thenoutputs asimmr_output object for further analysis and plottingviasummary.simmr_output andplot.simmr_output.
Usage
simmr_ffvb( simmr_in, prior_control = list(mu_0 = rep(0, simmr_in$n_sources), sigma_0 = 1), ffvb_control = list(n_output = 3600, S = 100, P = 10, beta_1 = 0.9, beta_2 = 0.9, tau = 100, eps_0 = 0.0225, t_W = 50))Arguments
simmr_in | An object created via the function |
prior_control | A list of values including arguments named |
ffvb_control | A list of values including arguments named |
Value
An object of classsimmr_output with two named top-levelcomponents:
input | The |
output | A set of outputs produced bythe FFVB function. These can be analysed using the |
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>, Emma Govan
References
Andrew C. Parnell, Donald L. Phillips, Stuart Bearhop, Brice X.Semmens, Eric J. Ward, Jonathan W. Moore, Andrew L. Jackson, Jonathan Grey,David J. Kelly, and Richard Inger. Bayesian stable isotope mixing models.Environmetrics, 24(6):387–399, 2013.
Andrew C Parnell, Richard Inger, Stuart Bearhop, and Andrew L Jackson.Source partitioning using stable isotopes: coping with too much variation.PLoS ONE, 5(3):5, 2010.
See Also
simmr_load for creating objects suitable for thisfunction,plot.simmr_input for creating isospace plots,summary.simmr_output for summarising output, andplot.simmr_output for plotting output.
Examples
## Not run: ## See the package vignette for a detailed run through of these 4 examples# Data set 1: 10 obs on 2 isos, 4 sources, with tefs and concdepdata(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# Printsimmr_1# FFVB runsimmr_1_out <- simmr_ffvb(simmr_1)# Print itprint(simmr_1_out)# Summarysummary(simmr_1_out, type = "correlations")summary(simmr_1_out, type = "statistics")ans <- summary(simmr_1_out, type = c("quantiles", "statistics"))# Plotplot(simmr_1_out, type = "boxplot")plot(simmr_1_out, type = "histogram")plot(simmr_1_out, type = "density")plot(simmr_1_out, type = "matrix")# Compare two sourcescompare_sources(simmr_1_out, source_names = c("Zostera", "Enteromorpha"))# Compare multiple sourcescompare_sources(simmr_1_out)###################################################################################### A version with just one observationdata(geese_data_day1)simmr_2 <- with( geese_data_day1, simmr_load( mixtures = mixtures[1, , drop = FALSE], source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_2)# FFVB run - automatically detects the single observationsimmr_2_out <- simmr_ffvb(simmr_2)# Print itprint(simmr_2_out)# Summarysummary(simmr_2_out)ans <- summary(simmr_2_out, type = c("quantiles"))# Plotplot(simmr_2_out)plot(simmr_2_out, type = "boxplot")plot(simmr_2_out, type = "histogram")plot(simmr_2_out, type = "density")plot(simmr_2_out, type = "matrix")###################################################################################### Data set 2: 3 isotopes (d13C, d15N and d34S), 30 observations, 4 sourcesdata(simmr_data_2)simmr_3 <- with( simmr_data_2, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Get summaryprint(simmr_3)# Plot 3 timesplot(simmr_3)plot(simmr_3, tracers = c(2, 3))plot(simmr_3, tracers = c(1, 3))# See vignette('simmr') for fancier axis labels# FFVB runsimmr_3_out <- simmr_ffvb(simmr_3)# Print itprint(simmr_3_out)# Summarysummary(simmr_3_out)summary(simmr_3_out, type = "quantiles")summary(simmr_3_out, type = "correlations")# Plotplot(simmr_3_out)plot(simmr_3_out, type = "boxplot")plot(simmr_3_out, type = "histogram")plot(simmr_3_out, type = "density")plot(simmr_3_out, type = "matrix")################################################################# Data set 5 - Multiple groups Geese data from Inger et al 2006# Do this in raw data format - Note that there's quite a few mixtures!data(geese_data)simmr_5 <- with( geese_data, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means, group = groups ))# Plotplot(simmr_5, xlab = expression(paste(delta^13, "C (per mille)", sep = "")), ylab = expression(paste(delta^15, "N (per mille)", sep = "")), title = "Isospace plot of Inger et al Geese data")# Run MCMC for each groupsimmr_5_out <- simmr_ffvb(simmr_5)# Summarise outputsummary(simmr_5_out, type = "quantiles", group = 1)summary(simmr_5_out, type = "quantiles", group = c(1, 3))summary(simmr_5_out, type = c("quantiles", "statistics"), group = c(1, 3))# Plot - only a single group allowedplot(simmr_5_out, type = "boxplot", group = 2, title = "simmr output group 2")plot(simmr_5_out, type = c("density", "matrix"), grp = 6, title = "simmr output group 6")# Compare sources within a groupcompare_sources(simmr_5_out, source_names = c("Zostera", "U.lactuca"), group = 2)compare_sources(simmr_5_out, group = 2)# Compare between groupscompare_groups(simmr_5_out, source = "Zostera", groups = 1:2)compare_groups(simmr_5_out, source = "Zostera", groups = 1:3)compare_groups(simmr_5_out, source = "U.lactuca", groups = c(4:5, 7, 2))## End(Not run)Function to load in simmr data and check for errors
Description
This function takes in the mixture data, food source means and standarddeviations, and (optionally) correction factor means and standarddeviations, and concentration proportions. It performs some (non-exhaustive)checking of the data to make sure it will run through simmr. It outputs anobject of classsimmr_input.
Usage
simmr_load( mixtures, source_names, source_means, source_sds, correction_means = NULL, correction_sds = NULL, concentration_means = NULL, group = NULL)Arguments
mixtures | The mixture data given as a matrix where the number of rowsis the number of observations and the number of columns is the number oftracers (usually isotopes) |
source_names | The names of the sources given as a character string |
source_means | The means of the source values, given as a matrix wherethe number of rows is the number of sources and the number of columns is thenumber of tracers |
source_sds | The standard deviations of the source values, given as amatrix where the number of rows is the number of sources and the number ofcolumns is the number of tracers |
correction_means | The means of the correction values, given as amatrix where the number of rows is the number of sources and the number ofcolumns is the number of tracers. If not provided these are set to 0. |
correction_sds | The standard deviations of the correction values,given as a matrix where the number of rows is the number of sources and thenumber of columns is the number of tracers. If not provided these are set to0. |
concentration_means | The means of the concentration values, given as amatrix where the number of rows is the number of sources and the number ofcolumns is the number of tracers. These should be between 0 and 1. If notprovided these are all set to 1. |
group | A grouping variable. These can be a character or factor variable |
Details
For standard stable isotope mixture modelling, the mixture matrix willcontain a row for each individual and a column for each isotopic value.simmr will allow for any number of isotopes and any number ofobservations, within computational limits. The source means/sds should beprovided for each food source on each isotope. The correction means (usuallytrophic enrichment factors) can be set as zero if required, and should be ofthe same shape as the source values. The concentration dependence meansshould be estimated values of the proportion of each element in the foodsource in question and should be given in proportion format between 0 and 1.At present there is no means to include concentration standard deviations.
Value
An object of classsimmr_input with the following elements:
mixtures | The mixture data |
source_neams | Source means |
sources_sds | Source standard deviations |
correction_means | Correction means |
correction_sds | Correction standard deviations |
concentration_means | Concentration dependence means |
n_obs | The number of observations |
n_tracers | The number oftracers/isotopes |
n_sources | The number of sources |
n_groups | The number of groups |
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
See Also
Seesimmr_mcmc for complete examples.
Examples
# A simple example with 10 observations, 2 tracers and 4 sourcesdata(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))print(simmr_1)Run asimmr_input object through the main simmr Markov chain MonteCarlo (MCMC) function
Description
This is the main function of simmr. It takes asimmr_input objectcreated viasimmr_load, runs an MCMC to determine the dietaryproportions, and then outputs asimmr_output object for furtheranalysis and plotting viasummary.simmr_output andplot.simmr_output.
Usage
simmr_mcmc( simmr_in, prior_control = list(means = rep(0, simmr_in$n_sources), sd = rep(1, simmr_in$n_sources), sigma_shape = rep(3, simmr_in$n_tracers), sigma_rate = rep(3/50, simmr_in$n_tracers)), mcmc_control = list(iter = 10000, burn = 1000, thin = 10, n.chain = 4))Arguments
simmr_in | An object created via the function |
prior_control | A list of values including arguments named: |
mcmc_control | A list of values including arguments named |
Details
If, after runningsimmr_mcmc the convergence diagnostics insummary.simmr_output are not satisfactory, the values ofiter,burn andthin inmcmc_control should beincreased by a factor of 10.
Value
An object of classsimmr_output with two named top-levelcomponents:
input | The |
output | A set of MCMC chains of class |
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>
References
Andrew C. Parnell, Donald L. Phillips, Stuart Bearhop, Brice X.Semmens, Eric J. Ward, Jonathan W. Moore, Andrew L. Jackson, Jonathan Grey,David J. Kelly, and Richard Inger. Bayesian stable isotope mixing models.Environmetrics, 24(6):387–399, 2013.
Andrew C Parnell, Richard Inger, Stuart Bearhop, and Andrew L Jackson.Source partitioning using stable isotopes: coping with too much variation.PLoS ONE, 5(3):5, 2010.
See Also
simmr_load for creating objects suitable for thisfunction,plot.simmr_input for creating isospace plots,summary.simmr_output for summarising output, andplot.simmr_output for plotting output.
Examples
## Not run: ## See the package vignette for a detailed run through of these 4 examples# Data set 1: 10 obs on 2 isos, 4 sources, with tefs and concdepdata(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# Printsimmr_1# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Print itprint(simmr_1_out)# Summarysummary(simmr_1_out, type = "diagnostics")summary(simmr_1_out, type = "correlations")summary(simmr_1_out, type = "statistics")ans <- summary(simmr_1_out, type = c("quantiles", "statistics"))# Plotplot(simmr_1_out, type = "boxplot")plot(simmr_1_out, type = "histogram")plot(simmr_1_out, type = "density")plot(simmr_1_out, type = "matrix")# Compare two sourcescompare_sources(simmr_1_out, source_names = c("Zostera", "Enteromorpha"))# Compare multiple sourcescompare_sources(simmr_1_out)###################################################################################### A version with just one observationdata(geese_data_day1)simmr_2 <- with( geese_data_day1, simmr_load( mixtures = mixtures[1, , drop = FALSE], source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_2)# MCMC run - automatically detects the single observationsimmr_2_out <- simmr_mcmc(simmr_2)# Print itprint(simmr_2_out)# Summarysummary(simmr_2_out)summary(simmr_2_out, type = "diagnostics")ans <- summary(simmr_2_out, type = c("quantiles"))# Plotplot(simmr_2_out)plot(simmr_2_out, type = "boxplot")plot(simmr_2_out, type = "histogram")plot(simmr_2_out, type = "density")plot(simmr_2_out, type = "matrix")###################################################################################### Data set 2: 3 isotopes (d13C, d15N and d34S), 30 observations, 4 sourcesdata(simmr_data_2)simmr_3 <- with( simmr_data_2, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Get summaryprint(simmr_3)# Plot 3 timesplot(simmr_3)plot(simmr_3, tracers = c(2, 3))plot(simmr_3, tracers = c(1, 3))# See vignette('simmr') for fancier axis labels# MCMC runsimmr_3_out <- simmr_mcmc(simmr_3)# Print itprint(simmr_3_out)# Summarysummary(simmr_3_out)summary(simmr_3_out, type = "diagnostics")summary(simmr_3_out, type = "quantiles")summary(simmr_3_out, type = "correlations")# Plotplot(simmr_3_out)plot(simmr_3_out, type = "boxplot")plot(simmr_3_out, type = "histogram")plot(simmr_3_out, type = "density")plot(simmr_3_out, type = "matrix")###################################################################################### Data set 5 - Multiple groups Geese data from Inger et al 2006# Do this in raw data format - Note that there's quite a few mixtures!data(geese_data)simmr_5 <- with( geese_data, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means, group = groups ))# Plotplot(simmr_5, xlab = expression(paste(delta^13, "C (per mille)", sep = "")), ylab = expression(paste(delta^15, "N (per mille)", sep = "")), title = "Isospace plot of Inger et al Geese data")# Run MCMC for each groupsimmr_5_out <- simmr_mcmc(simmr_5)# Summarise outputsummary(simmr_5_out, type = "quantiles", group = 1)summary(simmr_5_out, type = "quantiles", group = c(1, 3))summary(simmr_5_out, type = c("quantiles", "statistics"), group = c(1, 3))# Plot - only a single group allowedplot(simmr_5_out, type = "boxplot", group = 2, title = "simmr output group 2")plot(simmr_5_out, type = c("density", "matrix"), grp = 6, title = "simmr output group 6")# Compare sources within a groupcompare_sources(simmr_5_out, source_names = c("Zostera", "U.lactuca"), group = 2)compare_sources(simmr_5_out, group = 2)# Compare between groupscompare_groups(simmr_5_out, source = "Zostera", groups = 1:2)compare_groups(simmr_5_out, source = "Zostera", groups = 1:3)compare_groups(simmr_5_out, source = "U.lactuca", groups = c(4:5, 7, 2))## End(Not run)An artificial data set used to indicate effect of priors
Description
A fake box data set identified by Fry (2014) as a failing of SIMMsSee the link for more interpretation of these data and the output
Usage
square_dataFormat
A list with the following elements
- mixtures
A two column matrix containing delta 13C and delta 15N values respectively
- source_names
A character vector of the food source names
- tracer_names
A character vector of the tracer names (d13C, d15N)
- source_means
A matrix of source mean values for the tracers in the same order as
mixturesabove- source_sds
A matrix of source sd values for the tracers in the same order as
mixturesabove- correction_means
A matrix of TEFs mean values for the tracers in the same order as
mixturesabove- correction_sds
A matrix of TEFs sd values for the tracers in the same order as
mixturesabove- concentration_means
A matrix of concentration dependence mean values for the tracers in the same order as
mixturesabove
...
@seealsosimmr_mcmc for an example where it is used
Source
<doi:10.3354/meps10535>
Summarises the output created withsimmr_mcmc orsimmr_ffvb
Description
Produces textual summaries and convergence diagnostics for an object createdwithsimmr_mcmc orsimmr_ffvb. The differentoptions are: 'diagnostics' which produces Brooks-Gelman-Rubin diagnosticsto assess MCMC convergence, 'quantiles' which produces credible intervalsfor the parameters, 'statistics' which produces means and standarddeviations, and 'correlations' which produces correlations between theparameters.
Usage
## S3 method for class 'simmr_output'summary( object, type = c("diagnostics", "quantiles", "statistics", "correlations"), group = 1, ...)Arguments
object | An object of class |
type | The type of output required. At least none of 'diagnostics','quantiles', 'statistics', or 'correlations'. |
group | Which group or groups the output is required for. |
... | Not used |
Details
The quantile output allows easy calculation of 95 per cent credibleintervals of the posterior dietary proportions. The correlations, along withthe matrix plot inplot.simmr_output allow the user to judgewhich sources are non-identifiable. The Gelman diagnostic values should beclose to 1 to ensure satisfactory convergence.
When multiple groups are included, the output automatically includes theresults for all groups.
Value
A list containing the following components:
gelman | Theconvergence diagnostics |
quantiles | The quantiles of each parameterfrom the posterior distribution |
statistics | The means and standarddeviations of each parameter |
correlations | The posteriorcorrelations between the parameters |
Note that this object is reportedsilently so will be discarded unless the function is called with an objectas in the example below.
Author(s)
Andrew Parnell <andrew.parnell@mu.ie>, Emma Govan
See Also
Seesimmr_mcmc andsimmr_ffvbforcreating objects suitable for this function, and many more examples.See alsosimmr_load for creating simmr objects,plot.simmr_input for creating isospace plots,plot.simmr_output for plotting output.
Examples
# A simple example with 10 observations, 2 tracers and 4 sources# The datadata(geese_data_day1)simmr_1 <- with( geese_data_day1, simmr_load( mixtures = mixtures, source_names = source_names, source_means = source_means, source_sds = source_sds, correction_means = correction_means, correction_sds = correction_sds, concentration_means = concentration_means ))# Plotplot(simmr_1)# MCMC runsimmr_1_out <- simmr_mcmc(simmr_1)# Summarisesummary(simmr_1_out) # This outputs all the summariessummary(simmr_1_out, type = "diagnostics") # Just the diagnostics# Store the output in anans <- summary(simmr_1_out, type = c("quantiles", "statistics"))