knitr::opts_chunk$set( collapse = TRUE, comment = "#>" )
This is a tutorial covering the basic features of bayfoxr. This tutorial assumes that you have a recent version of R installed and are familiar with the R language.
bayfoxr is a suite of linear Bayesian calibration models for planktic core top foraminiferal δ18O (δ18Oc) and sea surface temperature (SST). These calibrations are especially useful because they capture the uncertainty in the relationship between modern SSTs and coretop δ18Oc. This package is a companion to a paper currently under preparation for the journal "Paleoceanography and Paleoclimatology".
Please cite our work if you use bayfoxr in your research. We have a paper currently in preparation and I'll be sure to update this section with the citation as soon as the paper is out.
To cite the code repository directly use:
Malevich, Steven B., 2019. bayfoxr. \<https://github.com/brews/bayfoxr >.
Alternatively, you can cite the package in R's CRAN repository. You can see this information by running
citation("bayfoxr") in an R session.
You can install bayfoxr from the CRAN repository with:
Development for bayfoxr is hosted at https://github.com/brews/bayfoxr. You can download and install directly from this github repository with the
This will give you the development versions of the package. It may be unstable. I only recommend this for advanced users.
Once bayfoxr has been installed, you can load it into your R session like any other R package:
bayfoxr revolves around two main functions,
predict_seatemp(). Unsurprisingly, these functions help us to predict δ18Oc or SST. Let's walk through a very simple case with some example sediment core data from Bass River. This data is bundled with bayfoxr, but feel free to use your own data -- just be sure it doesn't contain missing values.
data("bassriver") # Load the "bassriver" dataframe. head(bassriver) # Take a quick look at what the data looks like.
The example data are marine core samples from John et al. (2008). The dataframe has two columns: "depth", giving down-core depth in meters, and "d18o", foraminiferal (Morozovella spp.) calcite δ18O samples (‰ VPDB). The core samples cover the Paleocene-Eocene thermal maximum (PETM).
Now run this information through the
sst <- predict_seatemp(bassriver$d18o, d18osw = 0.0, prior_mean = 30.0, prior_std = 20.0)
The predict function spits out a
prediction object. Note that we need to specify δ18O for seawater in units ‰ VSMOW (
d18osw), and a prior mean and standard deviation for our SST inference.
sst variable contains an ensemble, or empirical distribution, rather than single prediction points because the calibration is a Bayesian regression model. This ensemble is in
sst[['ensemble']]. Here we get median and 90% interval for the prediction:
quants <- quantile(sst, probs = c(0.05, 0.50, 0.95)) head(quants) # Just see the top of the data...
We can also make a quick and dirty plot to visualize
predictplot(x = bassriver$depth, y = sst, ylim = c(20, 40), ylab = "SST (°C)", xlab = "Depth (m)")
You can see more options for
help(predictplot). This applies to any of the functions I've mentioned in this tutorial.
Of course, predictions with the
predict_d18oc() function are very similar to what we've already seen. Here we get a δ18Oc prediction from some made-up SST values. Note that we don't need to specify a prior mean or standard deviation with
d18oc <- predict_d18oc(c(24, 25, 23), d18osw = 0.0)
We can use with this prediction object just like before:
print(quantile(d18oc)) predictplot(y=d18oc, ylab="δ18Oc (‰ VPDB)", xlab='Sample', ylim = c(-3.5, 0))
There are four calibration models available for the
predict_seatemp() prediction functions.
The "pooled annual" model is the simplest case, as it has all foraminiferal species pooled together and calibrated against annual SST. The "hierarchical annual" model is similar, but accounts for some species-specific differences in calibration parameters. There are also the "pooled seasonal" and "hierarchical seasonal" models. These are similar to those described above, but use SSTs that are averaged to reflect general temperature preferences of foraminifera. See our 2019 paper for further details on these calibration models and their implementation.
predict_seatemp() use the pooled annual calibration model. To use a seasonal model set the
seasonal_seatemp argument to
TRUE. To use a hierarchical model, pass a foraminiferal species name string to the
foram argument. So, for example, our earlier
d18oc <- predict_d18oc(c(24, 25, 23), d18osw = 0.0)
uses the pooled annual model. To get other variations:
#d18oc <- predict_d18oc(c(24, 25, 23), d18osw = 0.0, seasonal_seatemp = TRUE) # Pooled seasonal d18oc <- predict_d18oc(c(24, 25, 23), d18osw = 0.0, foram = "G. bulloides") # Hierarchical annual for bulloides # Hierarchical seasonal for bulloides: d18oc <- predict_d18oc(c(24, 25, 23), d18osw = 0.0, seasonal_seatemp = TRUE, foram = "G. bulloides")
help(predict_seatemp) for more details including a full list of supported foram species.
Generally, the hierarchical annual, hierarchical seasonal, and pooled annual models are recommended. See our 2019 paper for further details
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.