A Not so Simple Workflow
+Simon Goring, Socorro Dominguez Vidaña
+2023-06-27
+ +Building New Chronologies
+This RMarkdown document will walk you through the process of:
+-
+
- Downloading a single record +
- Examining the chronologies for that record and associated +chronological controls +
- Creating a new chronology for the record +
- Adding the chronology to the record +
- 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 Lago Grande di Monticchio.
+Let’s start by pulling in the record and using the
+chronologies() helper function to look at the chronologies
+associated with the record:
# We could also search for Lago Monticchio:
+# monticchio <- get_sites(sitename = "Lago Grande di Monticchio", datasettype = "pollen") %>%
+# get_downloads()
+# But we know the datasetid so we can directly call get_downloads with the datasetid:
+monticchio <- get_downloads(4275)
+monticchio_chron <- chronologies(monticchio)
+monticchio_chron %>% as.data.frame()## .Lago Monticchio has three chronologies and we’ve decided not to use
+them because we want to standardize our methods. We want to build a new
+one with the function Bchronology() from the Bchron
+package.
It’s worth pointing out the isdefault column here.
+Neotoma provides the opportunity to link multiple chronologies to a
+single record. This lets researchers add their chronologies when they
+publish new studies. For example, the Wang et
+al. chronology comes from a set of Bayesian chronologies
+published by Yue Wang (Wang et al.,
+2019). For each age type (Radiocarbon years BP, Calendar years BP,
+&cetera) there is a default chronology that defines the model for
+the date interpolation. There is also a hierarchy for the default
+chronologies. By default Neotoma assigns the age model using
+calendar years the highest priority, then calibrated
+radiocarbon years, then radiocarbon years. You can see the
+order in practice if we look at the content of
+get_table("age types").
Extract chroncontrols
+We’re going to select chronology 25746 as our template.
+This is a Power interpolation model that _____ generated. To generate a
+new chronology for this record we want to see which chronological
+control points were used for the record. We will extract all the
+chroncontrols, filter by the chronologyid, and then arrange them by
+depth:
# Extract the chronological controls used in the original chronology:
+controls <- chroncontrols(monticchio) %>%
+ dplyr::filter(chronologyid == 25746) %>%
+ arrange(depth)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 -55 calibrated years BP, and assuming an uncertainty of 2 +years, and a thickness of 2cm.
+This generally won’t change too much, and I have no real basis for +doing this explicitly, but this is simply for illustration.
+To do these assignments we’re just directly modifying cells within
+the controls data.frame:
# Directly assign the values
+controls$chroncontrolage[1] <- 20000
+controls$agelimityounger[1] <- 15000
+controls$agelimitolder[1] <- 29000
+controls$thickness[1] <- 16Extract 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(monticchio) %>%
+ select(depth, analysisunitid) %>%
+ unique() %>%
+ arrange(depth)
+
+# Pass the values from `controls`. We're assuming the difference between
+# chroncontrolage and the agelimityounger is 1 SD.
+# Note that for the parameter 'calCurves' we are using a "normal"
+# distribution for the modern sample (core top) and choosing the
+# IntCal20 curve for the other two radiocarbon dates.
+
+newChron <- Bchron::Bchronology(ages = controls$chroncontrolage,
+ ageSds = abs(controls$agelimityounger -
+ controls$chroncontrolage),
+ calCurves = c("normal", rep("intcal20", 5)),
+ positionThicknesses = controls$thickness,
+ positions = controls$depth,
+ predictPositions = predictDepths$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)"
+ )
+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:
# Add information about the people who generated the new chronology:
+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"))
+
+# Add information about the chronology:
+newChronmonticchio <- set_chronology(agemodel = "Bchron model",
+ contact = creators,
+ isdefault = 1,
+ ageboundolder = max(newpredictions),
+ ageboundyounger = min(newpredictions),
+ dateprepared = lubridate::today(),
+ modelagetype = "Calendar years BP",
+ chronologyname = "Simon's example chronology",
+ chroncontrols = controls)
+
+
+newChronmonticchio$notes <- 'newChron <- Bchron::Bchronology(ages = controls$chroncontrolage,
+ ageSds = abs(controls$agelimityounger -
+ controls$chroncontrolage),
+ calCurves = c("normal", rep("intcal20", 2)),
+ 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:
+-
+
- we have a collectionunit in
monticchiothat is +accessible atmonticchio[[1]]$collunits.
+ - We can use the function
add_chronology(), which takes +the chronology object and adata.frame()of sample +ages.
+ - The predicted dates associated with the new chronology need to be
+transferred to each
samplesobject 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 = "Calendar years")
+
+monticchio[[1]]$collunits[[1]] <- add_chronology(monticchio[[1]]$collunits[[1]],
+ newChronmonticchio,
+ 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(monticchio) %>%
+ select(depth, analysisunitid, age) %>%
+ unique() %>%
+ arrange(depth) %>%
+ mutate(agecat = "new")
+
+monticchio[[1]]$collunits[[1]]$chronologies <- set_default(monticchio[[1]]$collunits[[1]]$chronologies,
+ 25746)
+plotforages <- samples(monticchio) %>%
+ 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")
+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:
+-
+
- Download a single record (the monticchio record using
+
get_downloads())
+ - Examining the chronologies for the record (using
+
chronologies()and associated chronological controls (using +chroncontrols())
+ - Creating a new chronology for the record (using
+
set_chronology())
+ - Adding the chronology to the record (using
+
add_chronology())
+ - 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!
+
+
+
+
+