fxTWAPLS: An Improved Version of WA-PLS

Overview

The goal of this package is to provide an improved version of WA-PLS by including the tolerances of taxa and the frequency of the sampled climate variable. This package also provides a way of leave-out cross-validation that removes both the test site and sites that are both geographically close and climatically close for each cycle, to avoid the risk of pseudo-replication.

Installation

You can install the released version of fxTWAPLS from CRAN with:

install.packages("fxTWAPLS")

And the development version from GitHub with:

install.packages("remotes")
remotes::install_github("special-uor/fxTWAPLS", "dev")

Publications

• Latest: Liu, M., Shen, Y., González-Sampériz, P., Gil-Romera, G., ter Braak, C. J. F., Prentice, I. C., and Harrison, S. P.: Holocene climates of the Iberian Peninsula: pollen-based reconstructions of changes in the west-east gradient of temperature and moisture, Clim. Past Discuss. [preprint], https://doi.org/10.5194/cp-2021-174, in review, 2021.- fxTWAPLS v0.1.0

  install.packages("remotes")
  remotes::install_github("special-uor/fxTWAPLS@v0.1.0")

• Liu Mengmeng, Prentice Iain Colin, ter Braak Cajo J. F., Harrison Sandy P.. 2020 An improved statistical approach for reconstructing past climates from biotic assemblages. Proc. R. Soc. A. 476: 20200346. https://doi.org/10.1098/rspa.2020.0346 - fxTWAPLS v0.0.2

  install.packages("remotes")
  remotes::install_github("special-uor/fxTWAPLS@v0.0.2")

Notes

The following functions can be executed in parallel:

• cv.pr.w
• cv.w
• get_distance
• get_pseudo
• sse.sample

To do so, include the cpus parameter. For example:

cv_pr_tf_Tmin2 <- fxTWAPLS::cv.pr.w(
  taxa,
  modern_pollen$Tmin,
  nPLS = 5,
  fxTWAPLS::TWAPLS.w2,
  fxTWAPLS::TWAPLS.predict.w,
  pseudo_Tmin,
  usefx = TRUE,
  fx_method = "pspline",
  bin = 0.02,
  cpus = 2
)

Optionally, a progress bar can be displayed for long computations. Just “pipe” the function call to fxTWAPLS::pb().

`%>%` <- magrittr::`%>%`
cv_pr_tf_Tmin2 <- fxTWAPLS::cv.pr.w(
  taxa,
  modern_pollen$Tmin,
  nPLS = 5,
  fxTWAPLS::TWAPLS.w2,
  fxTWAPLS::TWAPLS.predict.w,
  pseudo_Tmin,
  usefx = TRUE,
  fx_method = "pspline",
  bin = 0.02,
  cpus = 2
) %>%
  fxTWAPLS::pb()

Alternatively, if you are not familiar with the “pipe” operator, you can run the following code:

cv_pr_tf_Tmin2 <- fxTWAPLS::pb(
  fxTWAPLS::cv.pr.w(
    taxa,
    modern_pollen$Tmin,
    nPLS = 5,
    fxTWAPLS::TWAPLS.w2,
    fxTWAPLS::TWAPLS.predict.w,
    pseudo_Tmin,
    usefx = TRUE,
    fx_method = "pspline",
    bin = 0.02,
    cpus = 2
  )
)

Example

Training

# Load modern data
modern_pollen <- read.csv("/path/to/modern_pollen.csv")

# Extract modern pollen taxa
taxaColMin <- which(colnames(modern_pollen) == "taxa0")
taxaColMax <- which(colnames(modern_pollen) == "taxaN")
taxa <- modern_pollen[, taxaColMin:taxaColMax]

# Set the binwidth to get the sampling frequency of the climate (fx),
# the fit is almost insenitive to binwidth when choosing pspline method.
bin <- 0.02

# Use fxTWAPLSv2 to train
fit_tf_Tmin2 <- fxTWAPLS::TWAPLS.w2(
  taxa,
  modern_pollen$Tmin,
  nPLS = 5,
  usefx = TRUE,
  fx_method = "pspline",
  bin = bin
)

Cross validation

# Set CPUS to run in parallel
CPUS <- 6

# Import pipe operator to use with the progress bar
`%>%` <- magrittr::`%>%`

# Get the location information of each sample
point <- modern_pollen[, c("Long", "Lat")]

# Get the distance between each point
dist <- fxTWAPLS::get_distance(point, cpus = CPUS)

# Get the pseudo sites (which are both geographically close and climatically
# close to the test site) which should be removed in cross validation
pseudo_Tmin <- fxTWAPLS::get_pseudo(
  dist,
  modern_pollen$Tmin,
  cpus = CPUS
)

# Leave-out cross validation
cv_pr_tf_Tmin2 <- fxTWAPLS::cv.pr.w(
  taxa,
  modern_pollen$Tmin,
  nPLS = 5,
  fxTWAPLS::TWAPLS.w2,
  fxTWAPLS::TWAPLS.predict.w,
  pseudo_Tmin,
  usefx = TRUE,
  fx_method = "pspline",
  bin = bin,
  cpus = CPUS,
  test_mode = FALSE
) %>%
  fxTWAPLS::pb()

# Random t test to the cross validation result
rand_pr_tf_Tmin2 <-
  fxTWAPLS::rand.t.test.w(cv_pr_tf_Tmin2, n.perm = 999)

Reconstruction

# Load fossil data
Holocene <- read.csv("/path/to/Holocene.csv")

# Extract fossil pollen taxa
taxaColMin <- which(colnames(Holocene) == "taxa0")
taxaColMax <- which(colnames(Holocene) == "taxaN")
core <- Holocene[, taxaColMin:taxaColMax]

# Choose nsig (the last significant number of components) based on the p-value
nsig <- 3

# Predict
fossil_tf_Tmin2 <- fxTWAPLS::TWAPLS.predict.w(fit_tf_Tmin2, core)

# Get the sample specific errors
sse_tf_Tmin2 <- fxTWAPLS::sse.sample(
  modern_taxa = taxa,
  modern_climate = modern_pollen$Tmin,
  fossil_taxa = core,
  trainfun = fxTWAPLS::TWAPLS.w2,
  predictfun = fxTWAPLS::TWAPLS.predict.w,
  nboot = nboot,
  nPLS = 5,
  nsig = nsig,
  usefx = TRUE,
  fx_method = "pspline",
  bin = bin,
  cpus = CPUS
) %>%
  fxTWAPLS::pb()
# Output
recon_result <-
  cbind.data.frame(
    recon_Tmin = fossil_tf_Tmin2[["fit"]][, nsig],
    sse_recon_Tmin = sse_tf_Tmin2
  )