Building New Chronologies

This RMarkdown document will walk you through the process of:

  1. Downloading a single record
  2. Examining the chronologies for that record and associated chronological controls
  3. Creating a new chronology for the record
  4. Adding the chronology to the record
  5. Switching between default chronologies

This approach is focused on a single record, but much of what is done here can be extended to multiple records using functions.

Load Libraries

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)

Loading Datasets

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.

Extract chroncontrols

controls <- 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%"))

Extract Depth & Analysis Unit IDs

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.

Age-depth model for Stará Boleslav, with probability distributions superimposed on the figure at each chronology control depth.

Creating the New chronology and contact objects

Given 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)'

Adding the chronology to the collectionunit

Once 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:

  1. we have a collectionunit in stara that is accessible at stara[[1]]$collunits.
  2. We can use the function add_chronology(), which takes the chronology object and a data.frame() of sample ages.
  3. The predicted dates associated with the new chronology need to be transferred to each 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.

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.

Summary

From this notebook we have learned how to:

  1. Download a single record (the Stara record using get_downloads())
  2. Examining the chronologies for the record (using chronologies() and associated chronological controls (using chroncontrols())
  3. Creating a new chronology for the record (using set_chronology())
  4. Adding the chronology to the record (using add_chronology())
  5. Switching between default chronologies (using 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!