This RMarkdown document will walk you through the process of:
This approach is focused on a single record, but much of what is done here can be extended to multiple records using functions.
For this workshop element we only need four packages, neotoma2, dplyr, ggplot2 and Bchron. We’ll be loading a record from Neotoma, building a new chronology for the record, and then adding the chronology back to the record.
We’ll be using the R package pacman here (so really, we need five packages), to automatically load and install packages:
pacman::p_load(neotoma2, dplyr, ggplot2, Bchron)
We worked through the process for finding and downloading records using neotoma2 in the previous workshop. Assuming we found a record that we were interested in, we can go back and pull a single record using its datasetid. In this case, the dataset is for Stará Boleslav. Let’s start by pulling in the record and using the chronologies() helper function to look at the chronologies associated with the record:
stara <- get_downloads(24238)
## .
stara_chron <- chronologies(stara)
stara_chron %>% as.data.frame() %>%
DT::datatable(data = .,
options = list(scrollX = "100%"))
There are three chronologies here, but for whatever reason we’ve decided not to use any of them. We want to build a new one with the function Bchronology() from the Bchron package. First we want to see what chroncontrols we have for the prior chronologies. We’re going to select the chronologies used for chronology 14591 as our template.
chroncontrolscontrols <- chroncontrols(stara) %>%
dplyr::filter(chronologyid == 14591) %>%
arrange(depth)
controls %>% DT::datatable(data = .,
options = list(scrollX = "100%"))
We can look at other tools to decided how we want to manage the chroncontrols, for example, saving them and editing them using Excel or another spreadsheet program. We could add a new date by adding a new row. In this example we’re just going to modify the existing ages to provide better constraints at the core top. We are setting the core top to 0 calibrated years BP, and assuming an uncertainty of 2 years, and a thickness of 1cm.
To do these assignments we’re just directly modifying cells within the controls data.frame:
controls$chroncontrolage[1] <- 0
controls$agelimityounger[1] <- -2
controls$agelimitolder[1] <- 2
controls$thickness[1] <- 1
controls %>% DT::datatable(data = .,
options = list(scrollX = "100%"))
Once our chroncontrols table is updated, we extract the depths and analysisunitids from the dataset samples(). Pulling in both depths and analysisunitids is important because a single collection unit may have multiple datasets, which may have non-overlapping depth sequences. So, when adding sample ages back to a record we use the analysisunitid to make sure we are providing the correct assignment since depth may be specific to a single dataset.
# Get a two column data.frame with columns depth and analysisunitid.
# Sort the table by depth from top to bottom for "Bchronology"
predictDepths <- samples(stara) %>%
select(depth, analysisunitid) %>%
unique() %>%
arrange(depth)
# Pass the values from `controls`. We're assuming the difference between
# chroncontrolage and the agelimityounger is 1 SD.
newChron <- Bchron::Bchronology(ages = controls$chroncontrolage,
ageSds = abs(controls$agelimityounger -
controls$chroncontrolage),
calCurves = c("normal", rep("intcal20", 4)),
positionThicknesses = controls$thickness,
positions = controls$depth,
allowOutside = TRUE,
ids = controls$chroncontrolid)
# Predict ages at each depth for which we have samples. Returns a matrix.
newpredictions <- predict(newChron, predictDepths$depth)
plot(newChron) +
ggplot2::labs(
xlab = "Age (cal years BP)",
ylab = "Depth (cm)"
)
Age-depth model for Stará Boleslav, with probability distributions superimposed on the figure at each chronology control depth.
chronology and contact objectsGiven the new chronology, we want to add it to the sites object so that it becomes the default for any calls to samples(). To create the metadata for the new chronology, we use set_chronology() using the properties from the chronology table in Neotoma:
creators <- c(set_contact(givennames = "Simon James",
familyname = "Goring",
ORCID = "0000-0002-2700-4605"),
set_contact(givennames = "Socorro",
familyname = "Dominguez Vidaña",
ORCID = "0000-0002-7926-4935"))
newChronStara <- set_chronology(agemodel = "Bchron model",
contact = creators,
isdefault = 1,
ageboundolder = max(newpredictions),
ageboundyounger = min(newpredictions),
dateprepared = lubridate::today(),
modelagetype = "Calibrated radiocarbon years BP",
chronologyname = "Simon's example chronology",
chroncontrols = controls)
newChronStara$notes <- 'newChron <- Bchron::Bchronology(ages = controls$chroncontrolage,
ageSds = abs(controls$agelimityounger -
controls$chroncontrolage),
calCurves = c("normal", rep("intcal20", 4)),
positionThicknesses = controls$thickness,
positions = controls$depth,
allowOutside = TRUE,
ids = controls$chroncontrolid,
predictPositions = predictDepths)'
chronology to the collectionunitOnce we’ve created the chronology we need to apply it back to the collectionunit. We also need to add the predicted dates into the samples for each dataset associated with the collectionunit.
So:
stara that is accessible at stara[[1]]$collunits.add_chronology(), which takes the chronology object and a data.frame() of sample ages.samples object within the collectionunit.This is all bound up in the add_chronology() function, which takes the collectionunit, modifys it, and returns the newly updated collectionunit.
newSampleAges <- data.frame(predictDepths,
age = colMeans(newpredictions),
ageolder = colMeans(newpredictions) +
apply(newpredictions, 2, sd),
ageyounger = colMeans(newpredictions) -
apply(newpredictions, 2, sd),
agetype = "Calibrated radiocarbon years")
stara[[1]]$collunits[[1]] <- add_chronology(stara[[1]]$collunits[[1]],
newChronStara,
newSampleAges)
With this, we now have the updated collectionunit. Lets take a look at how this affects the age model overal. To pull the ages from the prior chronologies, we use the set_default() function to change the default chronology, and then extract ages, depths & analysisunits:
# The new chronology is currently the default chronology.
newages <- samples(stara) %>%
select(depth, analysisunitid, age) %>%
unique() %>%
arrange(depth) %>%
mutate(agecat = "new")
stara[[1]]$collunits[[1]]$chronologies <- set_default(stara[[1]]$collunits[[1]]$chronologies,
14591)
plotforages <- samples(stara) %>%
select(depth, analysisunitid, age) %>%
unique() %>%
arrange(depth) %>%
mutate(agecat = "old") %>%
bind_rows(newages)
And we can look at the difference visually:
ggplot(plotforages, aes(x = depth, y = age)) +
geom_path(aes(color = agecat)) +
theme_bw() +
xlab("Depth (cm)") +
ylab("Calibrated Years BP")
Differences in age representation between chronologies between existing chronologies and the new Bchron chronology.
So we can see the impact of the new chronology on the age model for the record, and we can make choices as to which model we want to use going forward. We can use this approach to create multiple new chronologies for a single record, tuning parameters within Bchronology(), or using Bacon and different parameters. Because the chronology is an R object we can save the objects for use in future sessions, and associate them with existin records, or we can re-run the models again.
From this notebook we have learned how to:
get_downloads())chronologies() and associated chronological controls (using chroncontrols())set_chronology())add_chronology())set_default())This approach is focused on a single record, but much of what is done here can be extended to multiple records using functions. We hope it’s been helpful!