sedproxy: Simulation of Sediment Archived Climate Proxy Records
Description
The sedproxy package provides functions to simulate sediment archived proxies
Author(s)
Maintainer: Andrew Dolman andrew.dolman@awi.de
Authors:
Thomas Laepple tlaepple@awi.de
See Also
Useful links:
Bioturbation weights
Description
For a given focal depth (or time), this function returns the probability that material collected from that depth was originally deposited at depth(s) z. In other words, that the material would have been found at depth z if there had been no bioturbation. It is the convolution of the depth solution from Berger and Heath (1968) with a uniform distribution to account for the width of the sediment layer from which samples were picked/extracted. It is a probability density function.
Usage
BioturbationWeights(
z,
focal.z,
layer.width = 1,
sed.acc.rate,
bio.depth,
scale = c("time", "depth")
)
Arguments
z
A vector of times or depths at which to evaluate the bioturbation weights
focal.z
The depth (or time) for which source dates are wanted
layer.width
the width of the sediment layer from which samples were taken, e.g. foraminifera were picked or alkenones were extracted, in cm. Defaults to 1 cm. If bio.depth and layer.width are both set to zero, each timepoint samples from a single year of the clim.signal, equivalent to sampling an annually laminated sediment core.
sed.acc.rate
Sediment accumulation rate in cm per 1000 years. Defaults to 50 cm per ka. Either a single value, or vector of same length as "timepoints"
bio.depth
Depth of the bioturbated layer in cm, defaults to 10 cm.
scale
Whether to scale depths by sediment accumulation rate to give positions in terms of time. Defaults to time.
Value
a numerical vector of weights.
References
Berger, W. H., & Heath, G. R. (1968). Vertical mixing in pelagic sediments. Journal of Marine Research, 26(2), 134–143.
Examples
z <- 0:10000
w <- BioturbationWeights(z, focal.z = 4000, layer.width = 1, sed.acc.rate = 5, bio.depth = 10)
plot(z, w, "l")
Simulate sediment archived proxy records from an input climate signal.
Description
ClimToProxyClim simulates the creation of a proxy climate
record from a climate signal that is assumed to be true.
The following aspects of proxy creation are currently modelled.
Seasonal bias in the encoding of a proxy due to the interaction between climate seasonality and any seasonality in the life cycle of the organism encoding the climate signal (e.g. Foraminifera for Mg/Ca ratios, or phytoplankton for Alkenone unsaturation indices).
Bioturbation of the sediment archived proxy. For each requested timepoint, the simulated proxy consists of a weighted mean of the climate signal over a time window that is determined by the sediment accumulation rate sed.acc.rate and the bioturbation depth bio.depth which defaults to 10 cm. The weights are given by the depth solution to an impulse response function (Berger and Heath, 1968).
Aliasing of seasonal and inter-annual climate variation onto to bioturbated (smoothed) signal. For proxies measured on a small number of discrete particles both seasonal and inter-annual climate variation is aliased into the proxy record. For example, Foraminifera have a life-cycle of approximately 1 month, so they record something like the mean temperature from a single month. If Mg/Ca is measured on e.g.
n.samples= 30 individuals, the measured proxy signal is a mean of 30 distinct monthly mean temperatures and will thus be a stochastic sample of the true mean climate.Measurement noise/error is added as a pure Gaussian white noise process with mean = 0, standard deviation =
sqrt(sigma.meas^2 + sigma.ind^2/n.samples).Additionally, a random bias can be added to each realisation of a proxy record. Bias is simulated as a Gaussian random variable with mean = 0, standard deviation =
meas.bias. The same randomly generated bias value is applied to all timepoints in a simulated proxy record, when multiple replicate proxies are generated (n.replicates > 1) each replicate has a different bias applied.
ClimToProxyClim returns one or more replicates of the final
simulated proxy as well as several intermediate stages (see section
Value below).
Usage
ClimToProxyClim(
clim.signal,
timepoints,
calibration.type = c("identity", "Uk37", "MgCa"),
calibration = switch(calibration.type, identity = NA, Uk37 = "Mueller global", MgCa =
"Ten planktonic species_350-500"),
slp.int.means = NULL,
slp.int.vcov = NULL,
noise.type = switch(calibration.type, identity = "additive", Uk37 = "additive", MgCa =
"multiplicative"),
plot.sig.res = 100,
habitat.weights = rep(1/ncol(clim.signal), ncol(clim.signal)),
habitat.wt.args = NULL,
bio.depth = 10,
sed.acc.rate = 50,
layer.width = 1,
sigma.meas = 0,
sigma.ind = 0,
meas.bias = 0,
scale.noise = switch(calibration.type, identity = FALSE, Uk37 = TRUE, MgCa = TRUE),
n.samples = Inf,
n.replicates = 1,
top.of.core = NULL,
n.bd = 3
)
Arguments
clim.signal
The "assumed true" climate signal, e.g. climate model
output or instrumental record. A ts object consisting of a
years x 12 (months) x n habitats (e.g. depths) matrix of temperatures. The
time series should be at annual resolution and in reverse, i.e. "most
recent timepoint first" order.
timepoints
The timepoints for which the proxy record is to be modelled
calibration.type
Type of proxy, e.g. Uk'37 or MgCa, to which the clim.signal is converted before the archiving and measurement of the proxy is simulated. Defaults to "identity" which means no conversion takes place.
calibration
The name of a specific calibration for which calibration parameters are provided by sedproxy. Currently applies only to calibration.type MgCa.
slp.int.means
Optional user supplied vector of values for the slope and intercept of the calibration function. Overides the defaults.
slp.int.vcov
Optional user supplied variance covariance matrix calibration parameters. Overides the defaults.
noise.type
Determines whether additive or multiplicative measurement noise is added. The appropriate type depends on the units of the proxy. Defaults to multiplicative for MgCa, additive for Uk'37 and identity (none) calibration types. Can be overidden with a string, "additive" or "multiplicative" in the case that pre-converted climate signal and measurement noise values are used in combination with an "identity" calibration type.
plot.sig.res
The resolution, in years, of the smoothed (block averaged) version of the input climate signal returned for plotting. This does not affect what the proxy model uses as input. If set to NA, no smoothed climate output is generated, this can speed up some simulations.
habitat.weights
Production weights for the proxy / proxy-carrier either as a vector of values with length = ncol(clim.signal), i.e. 1 weight for each month x habitat combination, a matrix of the same dimensions as the input climate signal matrix, or a function that produces an index of productivity as a function of temperature. Defaults to a vector of length = ncol(clim.signal) of equal weights.
habitat.wt.args
A named list of parameter values to be passed to a function named in habitat.weights.
bio.depth
Depth of the bioturbated layer in cm, defaults to 10 cm.
sed.acc.rate
Sediment accumulation rate in cm per 1000 years. Defaults to 50 cm per ka. Either a single value, or vector of same length as "timepoints"
layer.width
the width of the sediment layer from which samples were taken, e.g. foraminifera were picked or alkenones were extracted, in cm. Defaults to 1 cm. If bio.depth and layer.width are both set to zero, each timepoint samples from a single year of the clim.signal, equivalent to sampling an annually laminated sediment core.
sigma.meas
The standard deviation of the measurement error added to each simulated proxy value.
sigma.ind
The standard deviation of error between individuals (e.g. Forams) not otherwise modelled. This could included "vital effects" or aliasing of depth habitat variation not modelled via a depth resolved input climate signal and habitat weights. sigma.ind is scaled by n.samples before being combined with sigma.meas.
meas.bias
The amount of bias to add to each simulated proxy time-series. Each replicate proxy time-series has a constant bias added, drawn from a normal distribution with mean = 0, sd = meas.bias. Bias defaults to zero.
scale.noise
Scale noise to proxy units. Defaults to TRUE if calibration.type is not "identity"
n.samples
Number of e.g. Foraminifera sampled per timepoint, this can be either a single number, or a vector of length = timepoints. Can be set to Inf for non-discrete proxies, e.g. for Uk’37.
n.replicates
Number of replicate proxy time-series to simulate from the climate signal
top.of.core
The theoretical minimum age at the top of the core, ie. the year the core was sampled, defaults to the start of clim.in
n.bd
Number of multiples of the bioturbation width at which to truncate the bioturbation filter
Value
ClimToProxyClim returns an object of class "sedproxy.pfm", a list with three elements:
a dataframe
simulated.proxya dataframe
smoothed.signala dataframe
everything
The dataframe simulated.proxy contains a single realisation of the
final forward modelled proxy, as well as the intermediate stages and the
original climate signal at the requested timepoints.
The dataframe smoothed.signal contains a block averaged version the
input climate signal, defaults to 100 year means but this is set by the
parameter plot.sig.res. This is useful for plotting against the
resulting simulated proxy.
The dataframe everything contains all of the above but with multiple
replicates of the pseudo-proxy records if requested. The data are in
"long form", with the column "stage" inidcating the proxy stage or input
climate resolution and column "value" giving the values.
Named elements of the returned proxy record:
- timepoints
Requested timepoints
- clim.signal.ann
Input climate signal at requested timepoints at annual resolution
- clim.signal.smoothed
Input climate signal at regular time intervals and resolution = plot.sig.res
- clim.timepoints.ssr
Input climate signal at requested timepoints, smoothed to resolution = plot.sig.res
- proxy.bt
Climate signal after bioturbation
- proxy.bt.sb
Climate signal after bioturbation and habitat bias
- proxy.bt.sb.inf.b
Climate signal after bioturbation, habitat bias, and calibration bias
- proxy.bt.sb.inf.b.n
Climate signal after bioturbation, habitat bias, and measurement error
- proxy.bt.sb.sampY
Climate signal after bioturbation, habitat bias, and aliasing of inter-annual variation
- proxy.bt.sb.sampYM
Climate signal after bioturbation, habitat bias, and aliasing of inter-annual and intra-annual variation such as monthly temperatures or depth habitats
- proxy.bt.sb.sampYM.b
Climate signal after bioturbation, habitat bias, and aliasing of inter-annual and intra-annual variation such as monthly temperatures or depth habitats, and calibration bias
- proxy.bt.sb.sampYM.b.n
Climate signal after bioturbation, habitat bias, aliasing, and measurement error
- simulated.proxy
Final simulated pseudo-proxy, this will be same as proxy.bt.sb.inf.b.n when n.samples = Inf, and proxy.bt.sb.sampYM.b.n when n.samples is finite
- observed.proxy
True observed proxy (when supplied)
Examples
library(ggplot2)
set.seed(26052017)
clim.in <- ts(N41.t21k.climate[nrow(N41.t21k.climate):1,] - 273.15, start = -39)
PFM <- ClimToProxyClim(clim.signal = clim.in,
timepoints = round(N41.proxy$Published.age),
calibration.type = "identity",
habitat.weights = N41.G.ruber.seasonality,
sed.acc.rate = N41.proxy$Sed.acc.rate.cm.ka,
layer.width = 1,
sigma.meas = 0.46,
sigma.ind = 0, n.samples = Inf,
plot.sig.res = 10, meas.bias = 1,
n.replicates = 10)
PlotPFMs(PFM$everything, max.replicates = 1, stage.order = "seq") +
facet_wrap(~stage)
PlotPFMs(PFM$everything, max.replicates = 1, stage.order = "var")
PlotPFMs(PFM$everything, stage.order = "var", plot.stages = "all")
Foraminifer Growth Rate Function from Lombard et al. (2009)
Description
Implements the function for foraminifer growth rate as a function of temperature from Lombard et al. (2009) with parametrization from FAME 1.0 (Roche et al, 2018).
Usage
ForamGrowthfT(
foram = c("sacculifer", "bulloides", "pachy_d", "siphonifera", "universa", "pachy_s",
"dutertrei", "ruber"),
temperature_K,
norm = FALSE,
min.growth.thresh = 0
)
Arguments
foram
Name of foram.
temperature_K
Temperature in Kelvin
norm
Optional normalizing factor
min.growth.thresh
Sets a lower cutoff for growth as a proportion of the maximum growth rate for that taxon. For example in Roche et al (2018) a cutoff of 0.1 was used, meaning all growth rates less than 0.1*max were set to zero.
Value
A numerical vector or matrix with the same dimensions as the object passed to temperature_K. Units are daily growth rate, unless norm == TRUE.
References
Lombard, F., Labeyrie, L., Michel, E., Spero, H. J., and Lea, D. W.: Modelling the temperature dependent growth rates of planktic foraminifera, Marine Micropaleontology, 70, 1–7, https://doi.org/10.1016/j.marmicro.200809004, 2009.
Roche, D. M., Waelbroeck, C., Metcalfe, B. and Caley, T.: FAME (v1.0): a simple module to simulate the effect of planktonic foraminifer species-specific habitat on their oxygen isotopic content, Geosci. Model Dev. Discuss., 2017, 1–22, doi:10.5194/gmd-2017-251, 2017.
Examples
ForamGrowthfT(foram = 'ruber', temperature_K = (c(280, 290)), norm = 1)
Convert "everything" part of output from ClimToProxyClim to dataframe. Used internally.
Description
Convert "everything" part of output from ClimToProxyClim to dataframe. Used internally.
Usage
MakePFMDataframe(PFM)
Arguments
PFM
output from ClimToProxyClim
Value
a dataframe
Seasonality of Globigerinoides ruber at core MD97-2141
Description
Seasonality of Globigerinoides ruber at core MD97-2141 predicted by the PLAFOM model (Fraile et al. 2008).
Usage
N41.G.ruber.seasonality
Format
A vector of 12 values
Source
Fraile, I., Schulz, M., Mulitza, S., & Kucera, M. (2008). Predicting the global distribution of planktonic foraminifera using a dynamic ecosystem model. Biogeosciences, 5(3), 891–911.
Mg/Ca proxy based temperature reconstruction for core MD97-2141
Description
Mg/Ca proxy based temperature reconstruction for core MD97-2141, originally published in Rosenthal et al. (2003), extracted from Shakun et al. (2012).
Usage
N41.proxy
Format
A data frame with 216 rows and 4 variables:
Published.age Age in yr BP
Published.temperature Reconstructed temperature in degrees C
Sed.acc.rate.cm.ka Sediment accumulation rate in cm per kyr
Details
Published age and published temperature were extracted from Shakun, J. D., Clark, P. U., He, F., Marcott, S. A., Mix, A. C., Liu, Z., ... Bard, E. (2012). Global warming preceded by increasing carbon dioxide concentrations during the last deglaciation. Nature, 484(7392), 49–54. doi:10.1038/nature10915.
Sediment accumulation rates were estimated by fitting splines to the published age and depth estimates.
Source
Original reference: Rosenthal, Y., Oppo, D. W., & Linsley, B. K. (2003). The amplitude and phasing of climate change during the last deglaciation in the Sulu Sea, western equatorial Pacific. Geophysical Research Letters, 30(8), 1428. doi:10.1029/2002GL016612.
Metadata for datset N41.proxy
Description
Metadata for core MD97-2141 taken from Shakun et al. (2012)
Usage
N41.proxy.details
Format
A data frame with 1 rows and 17 variables:
Number Proxy ID number from Shakun et al. (2012)
ID.no Proxy ID number from Shakun et al. (2012) with prefix "N"
Core ID code of sediment core
Location Location of core
Proxy Proxy type
Lat Latutide of core
Lon Longitude of core
Elevation Depth of core-top below sea level in metres
Reference Original reference for the proxy record
Resolution Average time resolution in years
Calibration.ref Reference for Mg/Ca calibration
Calibration Mg/Ca calibration formula
Foram.sp Foram species analyses for proxy
Ref.14C Reference for carbon dating
Notes
Geo.cluster Coarse geographic location
Archive.type Type of proxy archive
Source
Shakun, J. D., Clark, P. U., He, F., Marcott, S. A., Mix, A. C., Liu, Z., ... Bard, E. (2012). Global warming preceded by increasing carbon dioxide concentrations during the last deglaciation. Nature, 484(7392), 49–54. doi:10.1038/nature10915
Climate (surface temperature) at core MD97-2141 from TraCE-21ka
Description
Modelled surface temperature at core MD97-2141. Model output from TraCE-21ka simulations.
Usage
N41.t21k.climate
Format
A matrix with 22040 rows and 12 columns
Source
Liu, Z., Otto-Bliesner, B. L., He, F., Brady, E. C., Tomas, R., Clark, P. U., ... Cheng, J. (2009). Transient Simulation of Last Deglaciation with a New Mechanism for Bølling-Allerød Warming. Science, 325(5938), 310–314. https://doi.org/10.1126/science.1171041
Plot forward modelled sedimentary proxies
Description
Plot forward modelled sedimentary proxies
Usage
PlotPFMs(
PFMs,
stage.order = c("var", "seq"),
plot.stages = c("default"),
max.replicates = 5,
colr.palette = "default",
alpha.palette = "default",
levl.labels = "default"
)
Arguments
PFMs
A dataframe of forward modelled proxies
stage.order
Controls the order in which proxy stages are plotted, either sequentially, "seq", or in order of variance, "var". Defaults to var.
plot.stages
Proxy stages to be plotted, "default", "all", or a custom character vector
max.replicates
Maximum number of replicates to plot at once
colr.palette
Colours for the proxy stages
alpha.palette
Alpha levels for the proxy stages
levl.labels
Labels for the proxy stages
Value
a ggplot object of class "gg" "ggplot"
Examples
library(ggplot2)
set.seed(26052017)
clim.in <- ts(N41.t21k.climate[nrow(N41.t21k.climate):1,] - 273.15)
PFM <- ClimToProxyClim(clim.signal = clim.in,
timepoints = round(N41.proxy$Published.age),
calibration.type = "identity",
habitat.weights = N41.G.ruber.seasonality,
sed.acc.rate = N41.proxy$Sed.acc.rate.cm.ka,
sigma.meas = 0.45,
sigma.ind = 0,
n.samples = Inf,
plot.sig.res = 10, meas.bias = 1,
n.replicates = 10)
PlotPFMs(PFM$everything, max.replicates = 1, stage.order = "seq") +
facet_wrap(~stage)
PlotPFMs(PFM$everything, max.replicates = 1, stage.order = "var")
PlotPFMs(PFM$everything, stage.order = "var", plot.stages = "all")
Convert between Temperature in Degrees C and Proxy Units
Description
A wrapper function for accessing proxy - temperature conversion functions
Usage
ProxyConversion(
temperature = NULL,
proxy.value = NULL,
calibration.type = "identity",
slp.int.means = NULL,
slp.int.vcov = NULL,
calibration = switch(calibration.type, identity = NA, Uk37 = "Mueller global", MgCa =
"Ten planktonic species_350-500"),
point.or.sample = c("point", "sample"),
n = 1
)
Arguments
temperature
Temperature in degrees C
proxy.value
Temperature in proxy units
calibration.type
Type of proxy, e.g. Uk37 or MgCa
slp.int.means
Optional user supplied vector of values for the slope and intercept of the calibration function. Overides the defaults.
slp.int.vcov
Optional user supplied variance covariance matrix calibration parameters. Overides the defaults.
calibration
The name of a specific calibration for which calibration parameters are provided by sedproxy. Currently applies only to calibration.type MgCa.
point.or.sample
Use the "best estimate" calibration parameters, or parameters sampled from the fitted calibration model
n
the number of replicate conversions to make in the case of sampled calibration parameters
Details
Valid entries for calibration are: "Ten planktonic species_350-500", "G. aequilateralis_350-500", "G. aequilateralis_500-1000", "G. conglobatus_350-500", "G. hirsuta_350-500", "G. inflata_350-500", "G. ruber pink_250-350", "G. ruber pink_350-500", "G. ruber white_250-350", "G. ruber white_350-500", "G. sacculifer with sac_350-500", "G. sacculifer without sac_350-500", "G. truncatulinoides_350-500", "G. truncatulinoides_500-1000", "N. dutertrei_350-500", "O. univesa_350-500", "P. obliquiloculata_350-500"
Value
a vector of temperatures or proxy values
Examples
# From temperature to UK'37
## With fixed calibration
ProxyConversion(temperature = c(10, 20), point.or.sample = "point",
calibration.type = "Uk37")
## With random calibration, 5 replicates
ProxyConversion(temperature = c(1, 2), n = 5, point.or.sample = "sample",
calibration.type = "Uk37")
## Back-transformation with same calibration
ProxyConversion(
proxy.value = as.vector(
ProxyConversion(
temperature = c(21, 22),
calibration.type = "Uk37",
point.or.sample = "point"
)
),
point.or.sample = "point",
calibration.type = "Uk37"
)
## Back-transformation with random calibration
ProxyConversion(
proxy.value = as.vector(
ProxyConversion(
temperature = c(21, 22),
calibration.type = "Uk37",
point.or.sample = "point"
)
)
,
n = 5,
point.or.sample = "sample",
calibration.type = "Uk37"
)
## Incompatible arguments
## Not run:
ProxyConversion(temperature = 1, proxy.value = 1)
## End(Not run)
Data for running the examples in package sedproxy
Description
A set of data objects for running the examples in sedproxy
Details
Example data
gisp2 ice core data at annual resolution
Description
gisp2 ice core data interpolated to annual resolution
Usage
gisp2.ann
Format
A data frame with 49885 rows and 3 variables:
age.yr.bpAge in years BP
temperaturetemperature in deg C
temperature.rescaledtemperature rescaled to resemble a d18O record as in Löwemark et al 2008
Details
GISP2 data from Alley, R.B.. 2004. GISP2 Ice Core Temperature and Accumulation Data. IGBP PAGES/World Data Center for Paleoclimatology Data Contribution Series 2004-013. NOAA/NGDC Paleoclimatology Program, Boulder CO, USA. Interpolated to annual resolution and additionally rescaled to resemble a d18O record. As in Fig.2 of Löwemark, L., Konstantinou, K. I. and Steinke, S.: Bias in foraminiferal multispecies reconstructions of paleohydrographic conditions caused by foraminiferal abundance variations and bioturbational mixing: A model approach, Marine Geology, 256(1–4), 101–106, doi:10.1016/j.margeo.200810005, 2008.
sedproxy parameters
Description
Parameters and variables required to generate a pseudo-proxy with ClimToProxyClim
Usage
param.tab
Format
A data frame with 12 rows and 4 variables:
Function argumentcharacter Argument name in
ClimToProxyClimDescriptioncharacter Description of argument and corresponding variable/parameter
Possible sourcescharacter Source or possible sources of values for arguments
Defaultcharacter Default values for arguments
Scussolini et al. (2013) Table 1
Description
Data from table 1 in Scussolini et al. (2013)
Usage
scussolini.tab1
Format
A dataframe with 22 rows and 6 columns
Source
Scussolini, P., van Sebille, E., & Durgadoo, J. V. (2013). Paleo Agulhas rings enter the subtropical gyre during the penultimate deglaciation. Climate of the Past, 9(6), 2631–2639. https://doi.org/10.5194/cp-9-2631-2013
Labels for proxy stages
Description
Labels for proxy stages. For plotting.
Usage
stage.labels
Format
A named character vector
Description of proxy stages
Description
A description of the proxy stages in the output of ClimToProxyClim
and default labels, colours and order for plotting
Usage
stages.key
Format
A data frame with 13 rows and 6 variables:
stage proxy stage
label label for proxy stage
description description of proxy stage
plot.order default plotting order of stages
plotting.colour default colour for plotting
plotting.alpha default alpha level for plotting