---
title: "Usage with the paleoTS package"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Usage with the paleoTS package}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

```{r setup}
library(StratPal)
```

This vignette explain how to link the `StratPal` package with the `paleoTS` package ([Hunt 2006](#References)). For details on the underlying data structures, see `vignette("StratPal_docs")`.

## Quick summary

To combine the `StratPal` and the `paleoTS` packages,

1. simulate trait evolution using the functions with suffix `_sl` (specimen level), optionally specifying the additional parameters for intrapopulation variance and number of specimens
2. build your pipelines as before (e.g., as described in `vignette("phenotypic_evolution")`)
3. turn the result into `paleoTS` format using `reduce_to_paleoTS`
4. further analyze or plot the results using the `paleoTS` package

## Motivation

The `paleoTS` package allows to easily analyze paleontological time series, and `StratPal` can be connected to it to study how ecological, taphonomic, and stratigraphic effects change our inferences about the tempo and mode of evolution.

The `paleoTS` package defines its own `paleoTS` format to store time series. It is a summary format, as it combines trait values measured in multiple specimens into an intrapopulation variance. Because taphonomic and ecological effects act on individual specimens, they can not be modeled directly on the `paleoTS` format. To circumvent this problem, we introduce a `pre_paleoTS` format that models trait evolution on the specimen level. Stratigraphic, taphonomic, and ecological effects can act on this format. After all these effects have been applied, we can reduce this data format into the standard `paleoTS` format, and then analyze it with the standard toolbox of the `paloeTS` package.

## Modeling trait evolution on specimen level, complexity reduction, and plotting

Functions to model trait evolution on specimen level have the suffix `_sl`, standing for "specimen level". Internally, they are of S3 class `pre_paleoTS`. In addition to the options provided to simulate mean trait values, you can simulate strict stasis using `strict_stasis_sl`. All functions have the same parametrization as the as the other functions for simulating trait evolution, which simulate mean trait evolution. In addition, they take two additional parameters: `intrapop_var` for the variance of the population around the mean trait value, and `n_per_sample` for the number of specimens per sampling location.

`pre_paleoTS` results can be converted to the `paleoTS` format using the function `reduce_to_paleoTS`. After this you can plot them with the standard plotting procedure from paleoTS using `plot` (resp, `plot.paleoTS`). Note that you can not plot `pre_paloeTS` objects directly, only after conversion to `paleoTS`.

```{r fig.alt="plot of strict stasis"}
library(StratPal)
library(paleoTS)   # needed for plotting

strict_stasis_sl(t = 1:4) |>   # simulate strict stasis on specimen level in `pre_paleoTS` format
  reduce_to_paleoTS() |>       # convert pre_paleoTS to paleoTS
  plot()                       # plot
```

## Modeling of ecological, taphonomic, and stratigraphic effects

Modeling of niches, taphonomic effects, and stratigraphic biases works on `pre_paleoTS` objects identical to how it is described in the vignettes on event data and modeling phenotypic evolution. As example, we plot a random walk 2 km offshore with 5 specimens per sampling location:

```{r, fig.alt="plot of a random walk after stratigraphic transformation"}
library(admtools)                                  # load admtools for stratigraphic transformation
adm = tp_to_adm(t = scenarioA$t_myr,               # define age-depth model
                h = scenarioA$h_m[,"2km"],
                L_unit = "m",
                T_unit = "Myr")    
set.seed(42)                                      # set seed for reproducibility
seq(min_time(adm), max_time(adm), by = 0.01) |>   # sample every 0.01 Myr
  random_walk_sl(n_per_sample = 10) |>            # simulate random walk on specimen level
  time_to_strat(adm) |>                           # transform into stratigraphic domain
  reduce_to_paleoTS() |>                          # transform into paleoTS format
  plot()                                          # plot
```

Of course you can also immediately add functions from `paleoTS` to the pipeline, e.g. to fit models of phenotypic evolution from the simulated data:

```{r}
set.seed(42)                                      # set seed for reproducibility
seq(min_time(adm), max_time(adm), by = 0.01) |>   # sample every 0.01 Myr
  random_walk_sl(n_per_sample = 10) |>            # simulate random walk on specimen level
  time_to_strat(adm) |>                           # transform into stratigraphic domain
  reduce_to_paleoTS() |>                          # transform into paleoTS format
  fit3models()                                    # fit 3 models to time series
```

## References {#References}

* Hunt, Gene. 2006. “Fitting and Comparing Models of Phyletic Evolution: Random Walks and Beyond.” Paleobiology. https://doi.org/10.1666/05070.1.