Title: | Simulation and Plots for Fossil and Taxonomy Data |
---|---|
Description: | Simulating and plotting taxonomy and fossil data on phylogenetic trees under mechanistic models of speciation, preservation and sampling. |
Authors: | Rachel Warnock [aut, cph], Joelle Barido-Sottani [aut, cre, cph], Walker Pett [aut, cph], O'Reilly Joseph [aut, cph], Ugnė Stolz [aut, cph] |
Maintainer: | Joelle Barido-Sottani <[email protected]> |
License: | GPL-3 |
Version: | 2.4.1 |
Built: | 2024-12-18 15:26:32 UTC |
Source: | https://github.com/fossilsim/fossilsim |
Transforms a tree and fossils into a sampled tree in beast-usable format and writes it in Newick format. Designed to work with FBD.
beast.fbd.format(tree, fossils, rho = 1, sampled_tips = NULL, ...)
beast.fbd.format(tree, fossils, rho = 1, sampled_tips = NULL, ...)
tree |
Complete tree. |
fossils |
fossils dataframe. |
rho |
Sampling probability of extant tips. Default 1, will be disregarded if sampled_tips is not null. |
sampled_tips |
List of tip labels corresponding to sampled extant tips. |
... |
Additional parameters will be passed to ape::write.tree |
Output of write.tree.
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # output for BEAST beast.fbd.format(t, f) # output on the console ## Not run: beast.fbd.format(t, f, file="example.tre") # output in file ## End(Not run)
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # output for BEAST beast.fbd.format(t, f) # output on the console ## Not run: beast.fbd.format(t, f, file="example.tre") # output in file ## End(Not run)
Count the total number of fossils
count.fossils(fossils)
count.fossils(fossils)
fossils |
Fossils object. |
Number of extinct samples.
Count the total number of fossils per interval
count.fossils.binned(fossils, interval.ages)
count.fossils.binned(fossils, interval.ages)
fossils |
Fossils object. |
interval.ages |
Vector of stratigraphic interval ages, starting with the minimum age of the youngest interval and ending with the maximum age of the oldest interval. |
Vector of extinct samples corresponding to each interval. Note the last value corresponds to the number of samples > the maximum age of the oldest interval.
Creates an fbdrange object from provided tree and data
fbdrange(tree, data)
fbdrange(tree, data)
tree |
|
data |
associated data to be included in the fbdrange object. Needs to contain |
the new fbdrange object
Create a fossil record object. The input is taken to be a dataframe or list.
fossils(data = NULL, from.taxonomy = FALSE) as.fossils(data, from.taxonomy = FALSE) is.fossils(data)
fossils(data = NULL, from.taxonomy = FALSE) as.fossils(data, from.taxonomy = FALSE) is.fossils(data)
data |
Dataframe or list of sampled fossils. See Details for the list of required fields. If NULL, the function creates an empty fossils object. |
from.taxonomy |
Boolean indicating whether the fossils were sampled using a taxonomy object, as opposed to a tree object. Default = FALSE. |
The fossil record object contains 4 fields for each fossil with the following information:
sp
the label of the corresponding species. This label matches the edge labels in the corresponding phylo object or the species labels in the corresponding taxonomy object if additional taxonomic
information was provided
edge
the label of the sampled node or tip in the phylogeny, i.e the node at the end of the edge along which the fossil was sampled
hmin
the age of the fossil or the youngest bound of the time interval in which the fossil was sampled
hmax
the oldest bound of the time interval in which the fossil was sampled.
This is equal to hmin
if exact sampling times are known
If complete = FALSE, only the extant taxa are used to construct the taxon constraints, resulting in a DPPDIV style analysis in which the extant topology is fixed and fossils can float in the tree. The resulting output uses the stronglyMonophyletic taxon constraint on the root, this means that all fossil taxa will be sampled in the crown group, and never in a position below the root.
fossils.to.BEAST.constraints( fossils, tree, file = "BEASTconstraints.xml", complete = FALSE, tree.name = "beastTree" )
fossils.to.BEAST.constraints( fossils, tree, file = "BEASTconstraints.xml", complete = FALSE, tree.name = "beastTree" )
fossils |
an object of class "fossils" that corresponds to fossil occurrences for the "tree" argument. |
tree |
an object of class "phylo", representing the tree upon which the fossil occurrences were simulated. |
file |
the name of the file to which the constraints will be written, defaults to "BEASTconstraints.xml". |
complete |
logical, if TRUE then taxon constraints are built for the complete tree, if FALSE then constraints are built for the crown clades only. Default value is FALSE. |
tree.name |
the name of the tree as used in the BEAST2 xml format. |
NULL.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) fossils.to.BEAST.constraints(f, t, file = tempfile(), complete = TRUE)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) fossils.to.BEAST.constraints(f, t, file = tempfile(), complete = TRUE)
If complete = FALSE, only the extant taxa are used to construct the tree, resulting in a DPPDIV style analysis in which the extant topology is fixed and fossils can float in the tree.
fossils.to.BEAST.start.tree(tree, fossils, complete = FALSE)
fossils.to.BEAST.start.tree(tree, fossils, complete = FALSE)
tree |
an object of class "phylo", representing the tree upon which the fossil occurrences were simulated. |
fossils |
an object of class "fossils" that corresponds to fossil occurrences for the "tree" argument. |
complete |
logical, if TRUE then the tree are built for the complete tree, if FALSE then the tree is built for the crown clades only. |
a string representing the starting tree in newick format.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) fossils.to.BEAST.start.tree(t,f, complete = FALSE)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) fossils.to.BEAST.start.tree(t,f, complete = FALSE)
Transforms a fossils dataframe and either taxonomy or tree into a fossilRecordSimulation object from package paleotree.
fossils.to.paleotree.record(fossils, tree = NULL, taxonomy = NULL)
fossils.to.paleotree.record(fossils, tree = NULL, taxonomy = NULL)
fossils |
fossils object |
tree |
phylo object containing the tree. If provided and taxonomy = NULL, all speciation is assumed symmetric |
taxonomy |
taxonomy object. If both tree and taxonomy are provided, only taxonomy will be used. |
The converted paleotree record
taxonomy
, fossils
, paleotree.record.to.fossils
# simulate tree t = ape::rtree(6) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) f = sim.fossils.poisson(2, taxonomy = s) # transform format record = fossils.to.paleotree.record(f, taxonomy = s)
# simulate tree t = ape::rtree(6) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) f = sim.fossils.poisson(2, taxonomy = s) # transform format record = fossils.to.paleotree.record(f, taxonomy = s)
Generate output in the format used by the program PyRate
fossils.to.pyrate( fossils, python = TRUE, traits = NULL, cutoff = NULL, random = FALSE, min = NULL, exclude.extant.singletons = TRUE, file = "", use.sp.names = FALSE )
fossils.to.pyrate( fossils, python = TRUE, traits = NULL, cutoff = NULL, random = FALSE, min = NULL, exclude.extant.singletons = TRUE, file = "", use.sp.names = FALSE )
fossils |
Fossils object. |
python |
If TRUE the function outputs the data in the python format used by PyRate (default). If FALSE the function outputs a tab-delimited table used by tools associated with PyRate. |
traits |
Vector of trait values equal to the number of unique species in the fossils dataframe.
The order should correspond to the order in which they appear in |
cutoff |
Exclude occurrences with age uncertainty greater than this value i.e. |
random |
If TRUE use a random number from within the interval U(hmin, hmax) for specimen ages,
otherwise use the midpoint of this interval (default). Applicable only when |
min |
Value used to represent the minimum possible interval age of extinct specimens with |
exclude.extant.singletons |
If TRUE exclude species that have extant samples only (default = TRUE). |
file |
Output file name. |
use.sp.names |
If TRUE use the value in fossils$sp as the complete taxon name, otherwise the function adds the prefix "taxa" (default = FALSE). |
set.seed(123) # simulate tree t = ape::rtree(6) # assign a max age based on tree height max.age = tree.max(t) # define a set of non-uniform length intervals times = c(0, sort(runif(3, min = 0, max = max.age)), max.age) rates = c(1,2,3,4) # simulate fossils reflect age uncertainty f = sim.fossils.intervals(tree = t, interval.ages = times, rates = rates, use.exact.times = FALSE) # simulate extant samples rho = 1 f = sim.extant.samples(f, t, rho = 1) plot(f, t) # generate input files for pyrate fossils.to.pyrate(f) fossils.to.pyrate(f, python = FALSE) # add trait values traits = runif(length(unique(f$sp))) fossils.to.pyrate(f, traits = traits)
set.seed(123) # simulate tree t = ape::rtree(6) # assign a max age based on tree height max.age = tree.max(t) # define a set of non-uniform length intervals times = c(0, sort(runif(3, min = 0, max = max.age)), max.age) rates = c(1,2,3,4) # simulate fossils reflect age uncertainty f = sim.fossils.intervals(tree = t, interval.ages = times, rates = rates, use.exact.times = FALSE) # simulate extant samples rho = 1 f = sim.extant.samples(f, t, rho = 1) plot(f, t) # generate input files for pyrate fossils.to.pyrate(f) fossils.to.pyrate(f, python = FALSE) # add trait values traits = runif(length(unique(f$sp))) fossils.to.pyrate(f, traits = traits)
This package provides functions for simulating both taxonomy and fossil data from an existing phylogeny.
Taxonomy can be simulated in FossilSim under a mixed model of speciation that can incorporate three modes of speciation – budding (or asymmetric), bifurcating (or symmetric) and anagenetic – in addition to cryptic speciation. A description of the resulting taxonomy objects and simulation functions can be found in the "Simulating taxonomy" vignette.
Fossils can be simulated from a phylogeny or a taxonomy under a model of constant fossil recovery or time-dependent, environment-dependent or species-dependent fossil recovery. A description of the resulting fossil objects and simulation functions can be found in the "Simulating fossils" vignette.
Both taxonomy and fossil objects are provided with custom plotting functions that highlight important features of the simulated objects
along the original phylogeny. More details about these functions can be found in the vignettes or by calling
?plot.taxonomy
and ?plot.fossils
.
FossilSim is designed to use phylogenies in the ape format. It provides functions to convert to and from the fossilRecordSimulation format used by the package paleotree (see the vignette "Converting from and to paleotree format"), as well as functions to convert to the zero-edge format used by BEAST2 and RevBayes (see the vignette "Exporting sampled ancestor trees").
Maintainer: Joelle Barido-Sottani [email protected] [copyright holder]
Authors:
Rachel Warnock [email protected] [copyright holder]
Walker Pett [copyright holder]
O'Reilly Joseph [copyright holder]
Ugnė Stolz [copyright holder]
# simulate a tree using TreeSim conditioned on tip number t = TreeSim::sim.bd.taxa(n = 10, numbsim = 1, lambda = 1, mu = 0.2)[[1]] # simulate taxonomy under mixed speciation s = sim.taxonomy(tree = t, beta = 0.5, lambda.a = 1, kappa = 0.1) # plot the result plot(s, tree = t, legend.position = "topleft") # simulate fossils using the phylogeny and a constant fossil recovery rate f = sim.fossils.poisson(rate = 3, tree = t) # plot the result plot(f, tree = t) # simulate fossils using the taxonomy and a constant fossil recovery rate f = sim.fossils.poisson(rate = 3, taxonomy = s) # plot the result plot(f, tree = t, taxonomy = s, show.taxonomy = TRUE) # simulate fossils using time-dependent fossil recovery rates f = sim.fossils.intervals(tree = t, rates = c(1, 0.1, 1, 0.1), max.age = tree.max(t), strata = 4) # plot the result, with the time intervals plot(f, tree = t, show.strata = TRUE, max.age = tree.max(t), strata = 4)
# simulate a tree using TreeSim conditioned on tip number t = TreeSim::sim.bd.taxa(n = 10, numbsim = 1, lambda = 1, mu = 0.2)[[1]] # simulate taxonomy under mixed speciation s = sim.taxonomy(tree = t, beta = 0.5, lambda.a = 1, kappa = 0.1) # plot the result plot(s, tree = t, legend.position = "topleft") # simulate fossils using the phylogeny and a constant fossil recovery rate f = sim.fossils.poisson(rate = 3, tree = t) # plot the result plot(f, tree = t) # simulate fossils using the taxonomy and a constant fossil recovery rate f = sim.fossils.poisson(rate = 3, taxonomy = s) # plot the result plot(f, tree = t, taxonomy = s, show.taxonomy = TRUE) # simulate fossils using time-dependent fossil recovery rates f = sim.fossils.intervals(tree = t, rates = c(1, 0.1, 1, 0.1), max.age = tree.max(t), strata = 4) # plot the result, with the time intervals plot(f, tree = t, show.strata = TRUE, max.age = tree.max(t), strata = 4)
Import fbdrange object from file
get_fbdrange_from_file(input_file)
get_fbdrange_from_file(input_file)
input_file |
path to a tree file in Nexus format containing range and orientation metadata for the trees (for instance, a file generated by the BEAST2 package sRanges) |
an fbdrange object
tree_file <- system.file("extdata", "fbdrange.trees", package = "FossilSim") fbdr <- get_fbdrange_from_file(tree_file)
tree_file <- system.file("extdata", "fbdrange.trees", package = "FossilSim") fbdr <- get_fbdrange_from_file(tree_file)
Obtain the tips that define each node in a tree
get.tip.descs(tree)
get.tip.descs(tree)
tree |
an object of class "Phylo". |
A list of vectors, with one entry for each node consisting of the tip labels that define that node.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] get.tip.descs(t)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] get.tip.descs(t)
The returned tree is in paleotree format, with zero-length edges leading to tips at bifurcation and anagenetic events. Fossils and taxonomy are only specified on non-zero-length edges. The label assigned to the parent of the origin or root will be zero.
paleotree.record.to.fossils(record)
paleotree.record.to.fossils(record)
record |
fossilRecordSimulation object. |
A list containing the converted tree, taxonomy and fossils
taxonomy
, fossils
, fossils.to.paleotree.record
if (requireNamespace("paleotree", quietly = TRUE)) { # simulate record record = paleotree::simFossilRecord(p=0.1, q=0.1,r=0.1, nruns=1, nTotalTaxa=c(30,40), nExtant=0, nSamp = c(5,25)) # transform format l_tf = paleotree.record.to.fossils(record) l_tf$tree l_tf$taxonomy l_tf$fossils }
if (requireNamespace("paleotree", quietly = TRUE)) { # simulate record record = paleotree::simFossilRecord(p=0.1, q=0.1,r=0.1, nruns=1, nTotalTaxa=c(30,40), nExtant=0, nSamp = c(5,25)) # transform format l_tf = paleotree.record.to.fossils(record) l_tf$tree l_tf$taxonomy l_tf$fossils }
If "ext.tree" is not supplied, this function will find the direct ancestral node for each of the supplied fossil samples. If "ext.tree" is supplied, this function will find the direct ancestral node for each fossil in "ext.tree". This second behaviour is used for placing fossils simulated on a complete Birth-Death tree in the extant-only counterpart tree. This results in fossil samples being placed in the crown clades of the tree upon which they were simulated. When "ext.tree" is supplied, any fossil samples appearing before the MRCA of the crown group are discarded.
place.fossils(tree, fossils, ext.tree)
place.fossils(tree, fossils, ext.tree)
tree |
an object of class "Phylo". |
fossils |
an object of class "fossils" that corresponds to fossil occurrences for the "tree" argument. |
ext.tree |
an object of class "Phylo" representing the extant counterpart to "tree", this can be obtained with prune.fossil.tips(tree). |
a vector of node numbers corresponding to the direct ancestor of each fossil sample in "fossils".
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) place.fossils(t,f)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) place.fossils(t,f)
Plot oriented tree with stratigraphic ranges
## S3 method for class 'fbdrange' plot(x, smart.labels = FALSE, ...)
## S3 method for class 'fbdrange' plot(x, smart.labels = FALSE, ...)
x |
object of type |
smart.labels |
whether to label the ranges (default |
... |
other arguments to be passed to the plot labels (does nothing if |
a ggtree plot which can be combined with any other commands from [ggplot2] or [ggtree]
tree_file <- system.file("extdata", "fbdrange.trees", package = "FossilSim") fbdr <- get_fbdrange_from_file(tree_file) p <- plot(fbdr, smart.labels = TRUE)
tree_file <- system.file("extdata", "fbdrange.trees", package = "FossilSim") fbdr <- get_fbdrange_from_file(tree_file) p <- plot(fbdr, smart.labels = TRUE)
This function is adapted from the ape function plot.phylo
used to plot phylogenetic trees.
The function can be used to plot simulated fossils (show.fossils = TRUE
), with or without the corresponding tree (show.tree = TRUE
),
stratigraphic intervals (show.strata = TRUE
), stratigraphic ranges (show.ranges = TRUE
) and sampling proxy data (show.proxy = TRUE
).
Interval ages can be specified as a vector (interval.ages
) or a uniform set of interval ages can be specified using the
number of intervals (strata
) and maximum interval age (max
), where interval length .
If no maximum age is specified, the function calculates a maximum interval age slightly older than the root edge (or root age if
root.edge = FALSE
),
using the function tree.max()
.
## S3 method for class 'fossils' plot( x, tree, show.fossils = TRUE, show.tree = TRUE, show.ranges = FALSE, show.strata = FALSE, strata = 1, max.age = NULL, interval.ages = NULL, binned = FALSE, show.axis = TRUE, show.proxy = FALSE, proxy.data = NULL, show.preferred.environ = FALSE, preferred.environ = NULL, show.taxonomy = FALSE, taxonomy = NULL, show.unknown = FALSE, rho = 1, root.edge = TRUE, hide.edge = FALSE, edge.width = 1, show.tip.label = FALSE, align.tip.label = FALSE, reconstructed = FALSE, fossil.col = 1, range.col = rgb(0, 0, 1), extant.col = 1, taxa.palette = "viridis", col.axis = "gray35", cex = 1.2, pch = 18, ... )
## S3 method for class 'fossils' plot( x, tree, show.fossils = TRUE, show.tree = TRUE, show.ranges = FALSE, show.strata = FALSE, strata = 1, max.age = NULL, interval.ages = NULL, binned = FALSE, show.axis = TRUE, show.proxy = FALSE, proxy.data = NULL, show.preferred.environ = FALSE, preferred.environ = NULL, show.taxonomy = FALSE, taxonomy = NULL, show.unknown = FALSE, rho = 1, root.edge = TRUE, hide.edge = FALSE, edge.width = 1, show.tip.label = FALSE, align.tip.label = FALSE, reconstructed = FALSE, fossil.col = 1, range.col = rgb(0, 0, 1), extant.col = 1, taxa.palette = "viridis", col.axis = "gray35", cex = 1.2, pch = 18, ... )
x |
Fossils object. |
tree |
Phylo object. |
show.fossils |
If TRUE plot fossils (default = TRUE). |
show.tree |
If TRUE plot the tree (default = TRUE). |
show.ranges |
If TRUE plot stratigraphic ranges (default = FALSE). If |
show.strata |
If TRUE plot strata (default = FALSE). |
strata |
Number of stratigraphic intervals (default = 1). |
max.age |
Maximum age of a set of equal length intervals. If no value is specified ( |
interval.ages |
Vector of stratigraphic interval ages, starting with the minimum age of the youngest interval and ending with the maximum age of the oldest interval. |
binned |
If TRUE fossils are plotted at the mid point of each interval. |
show.axis |
If TRUE plot x-axis (default = TRUE). |
show.proxy |
If TRUE add profile of sampling data to plot (e.g. rates in time-dependent rates model) (default = FALSE). |
proxy.data |
Vector of sampling proxy data (default = NULL). Should be as long as the number of stratigraphic intervals. |
show.preferred.environ |
If TRUE add species preferred environmental value (e.g. water depth) (default = FALSE). Only works if combined with |
preferred.environ |
Preferred environmental value (e.g. water depth). Currently only one value can be shown. |
show.taxonomy |
If TRUE highlight species taxonomy. |
taxonomy |
Taxonomy object. |
show.unknown |
If TRUE plot fossils with unknown taxonomic affiliation (i.e. sp = NA) (default = FALSE). |
rho |
Extant species sampling probability (default = 1). Will be disregarded if fossils object already contains extant samples. |
root.edge |
If TRUE include the root edge (default = TRUE). |
hide.edge |
If TRUE hide the root edge but still incorporate it into the automatic timescale (default = FALSE). |
edge.width |
A numeric vector giving the width of the branches of the plotted phylogeny. These are taken to be in the same order as the component edge of |
show.tip.label |
Whether to show the tip labels on the phylogeny (defaults to FALSE). |
align.tip.label |
A logical value or an integer. If TRUE, the tips are aligned and dotted lines are drawn between the tips of the tree and the labels. If an integer, the tips are aligned and this gives the type of the lines (following |
reconstructed |
If TRUE plot the reconstructed tree. If fossils object contains no extant samples, the function assumes rho = 1 and includes all species at the present. |
fossil.col |
Colour of fossil occurrences. A vector equal to the length of the fossils object can be used to assign different colours. |
range.col |
Colour of stratigraphic ranges. |
extant.col |
Colour of extant samples. If |
taxa.palette |
Colour palette used if |
col.axis |
Colour of the time scale axis (default = "gray35"). |
cex |
Numeric value giving the factor used to scale the points representing the fossils when |
pch |
Numeric value giving the symbol used for the points representing the fossils when |
... |
Additional parameters to be passed to |
set.seed(123) ## simulate tree t = TreeSim::sim.bd.taxa(8, 1, 1, 0.3)[[1]] ## simulate fossils under a Poisson sampling process f = sim.fossils.poisson(rate = 3, tree = t) plot(f, t) # add a set of equal length strata plot(f, t, show.strata = TRUE, strata = 4) # show stratigraphic ranges plot(f, t, show.strata = TRUE, strata = 4, show.ranges = TRUE) ## simulate fossils and highlight taxonomy s = sim.taxonomy(t, 0.5, 1) f = sim.fossils.poisson(rate = 3, taxonomy = s) plot(f, t, taxonomy = s, show.taxonomy = TRUE, show.ranges = TRUE) ## simulate fossils under a non-uniform model of preservation # assign a max interval based on tree height max.age = tree.max(t) times = c(0, 0.3, 1, max.age) rates = c(4, 1, 0.1) f = sim.fossils.intervals(t, interval.ages = times, rates = rates) plot(f, t, show.strata = TRUE, interval.ages = times) # add proxy data plot(f, t, show.strata = TRUE, interval.ages = times, show.proxy = TRUE, proxy.data = rates)
set.seed(123) ## simulate tree t = TreeSim::sim.bd.taxa(8, 1, 1, 0.3)[[1]] ## simulate fossils under a Poisson sampling process f = sim.fossils.poisson(rate = 3, tree = t) plot(f, t) # add a set of equal length strata plot(f, t, show.strata = TRUE, strata = 4) # show stratigraphic ranges plot(f, t, show.strata = TRUE, strata = 4, show.ranges = TRUE) ## simulate fossils and highlight taxonomy s = sim.taxonomy(t, 0.5, 1) f = sim.fossils.poisson(rate = 3, taxonomy = s) plot(f, t, taxonomy = s, show.taxonomy = TRUE, show.ranges = TRUE) ## simulate fossils under a non-uniform model of preservation # assign a max interval based on tree height max.age = tree.max(t) times = c(0, 0.3, 1, max.age) rates = c(4, 1, 0.1) f = sim.fossils.intervals(t, interval.ages = times, rates = rates) plot(f, t, show.strata = TRUE, interval.ages = times) # add proxy data plot(f, t, show.strata = TRUE, interval.ages = times, show.proxy = TRUE, proxy.data = rates)
This function is adapted from the ape function plot.phylo
used to plot phylogenetic trees.
The function can be used to plot simulated taxonomy along with the corresponding tree.
## S3 method for class 'taxonomy' plot( x, tree, show.mode = TRUE, show.legend = TRUE, legend.position = "bottomleft", root.edge = TRUE, hide.edge = FALSE, edge.width = 1, show.tip.label = FALSE, align.tip.label = FALSE, taxa.palette = "viridis", cex = 1.2, ... )
## S3 method for class 'taxonomy' plot( x, tree, show.mode = TRUE, show.legend = TRUE, legend.position = "bottomleft", root.edge = TRUE, hide.edge = FALSE, edge.width = 1, show.tip.label = FALSE, align.tip.label = FALSE, taxa.palette = "viridis", cex = 1.2, ... )
x |
Taxonomy object. |
tree |
Phylo object. |
show.mode |
Indicate speciation mode. |
show.legend |
Add a legend for the symbols indicating different speciation modes. |
legend.position |
Position of the legend. Options include |
root.edge |
If TRUE include the root edge (default = TRUE). |
hide.edge |
If TRUE hide the root edge but still incorporate it into the automatic timescale (default = FALSE). |
edge.width |
A numeric vector giving the width of the branches of the plotted phylogeny. These are taken to be in the same order as the component edge of |
show.tip.label |
Whether to show the tip labels on the phylogeny (defaults to FALSE). |
align.tip.label |
A logical value or an integer. If TRUE, the tips are aligned and dotted lines are drawn between the tips of the tree and the labels. If an integer, the tips are aligned and this gives the type of the lines (lty). |
taxa.palette |
Colour palette used for taxa. Colours are assigned to taxa using the function |
cex |
Numeric value giving the factor used to scale the points representing the fossils when |
... |
Additional parameters to be passed to |
set.seed(123) ## simulate tree t = TreeSim::sim.bd.taxa(8, 1, 1, 0.3)[[1]] ## simulate taxonomy s = sim.taxonomy(t, 0.5, 1) ## plot the output plot(s, t)
set.seed(123) ## simulate tree t = TreeSim::sim.bd.taxa(8, 1, 1, 0.3)[[1]] ## simulate taxonomy s = sim.taxonomy(t, 0.5, 1) ## plot the output plot(s, t)
Remove fossil lineages from a tree
prune.fossil.tips(tree)
prune.fossil.tips(tree)
tree |
an object of class "Phylo". |
an object of class "Phylo". If fossil lineages were found in the tree these will be pruned, if not then the original tree is returned.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] prune.fossil.tips(t)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] prune.fossil.tips(t)
Removes all intermediate fossils from a combined tree and labels the first and last fossils of each lineage. Can be used with sampled or complete trees. If only one fossil is present for a particular species it is labelled as first.
prune.fossils(tree)
prune.fossils(tree)
tree |
Combined tree with fossils. |
Tree with pruned fossils.
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # transform format t2 = SAtree.from.fossils(t,f)$tree # prune fossils t4 = prune.fossils(t2) # or transform to sampled tree first t3 = sampled.tree.from.combined(t2) t4 = prune.fossils(t3) plot(t4)
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # transform format t2 = SAtree.from.fossils(t,f)$tree # prune fossils t4 = prune.fossils(t2) # or transform to sampled tree first t3 = sampled.tree.from.combined(t2) t4 = prune.fossils(t3) plot(t4)
Make an asymmetric stratigraphic range plot from a tree object of class phylo
rangeplot.asymmetric(x, complete = FALSE, ...)
rangeplot.asymmetric(x, complete = FALSE, ...)
x |
phylo object to plot. |
complete |
Plot unsampled species. |
... |
Additional parameters to be passed to |
This function assumes all speciation events are asymmetric.
tree = sim.fbd.taxa(n = 10, numbsim = 1, lambda = 3, mu = 2, psi = 1, complete = TRUE)[[1]] rangeplot.asymmetric(tree, complete=TRUE)
tree = sim.fbd.taxa(n = 10, numbsim = 1, lambda = 3, mu = 2, psi = 1, complete = TRUE)[[1]] rangeplot.asymmetric(tree, complete=TRUE)
This function uses edge identifiers (edge
) and fossil sampling times (hmin
) to reassign fossil species identifiers (sp, origin
) using an existing taxonomy object.
It can only be used if exact fossil sampling times are known (i.e. hmin = hmax
), otherwise edges containing multiple species may be indistinguishable.
reconcile.fossils.taxonomy(fossils, taxonomy)
reconcile.fossils.taxonomy(fossils, taxonomy)
fossils |
Fossils object. |
taxonomy |
Taxonomy object. |
An object of class fossils.
# simulate tree t = ape::rtree(6) # simulate fossils using the tree rate = 2 f = sim.fossils.poisson(rate, tree = t) plot(f, t) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) f = reconcile.fossils.taxonomy(f, s) plot(f, t)
# simulate tree t = ape::rtree(6) # simulate fossils using the tree rate = 2 f = sim.fossils.poisson(rate, tree = t) plot(f, t) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) f = reconcile.fossils.taxonomy(f, s) plot(f, t)
Note that for datasets containing extinct only samples (& rho = 0) the ages output are scaled so that the youngest sample = 0.
reconstructed.tree.fossils.objects(fossils, tree, rho = 1)
reconstructed.tree.fossils.objects(fossils, tree, rho = 1)
fossils |
Fossils object. |
tree |
Tree object. |
rho |
Extant species sampling probability. Default = 1, will be disregarded if fossils object already contains extant samples. |
A list containing the tree and fossil objects.
# simulate tree birth = 0.1 death = 0.05 tips = 10 t = TreeSim::sim.bd.taxa(tips, 1, birth, death)[[1]] # simulate fossils f = sim.fossils.poisson(rate = 0.3, tree = t) # simulate extant samples f = sim.extant.samples(f, tree = t, rho = 0.5) # plot the complete tree plot(f,t) # generate tree & fossil objects corresponding to the reconstructed tree out = reconstructed.tree.fossils.objects(f, t) f.reconst = out$fossils t.reconst = out$tree # plot the reconstructed tree plot(f.reconst, t.reconst)
# simulate tree birth = 0.1 death = 0.05 tips = 10 t = TreeSim::sim.bd.taxa(tips, 1, birth, death)[[1]] # simulate fossils f = sim.fossils.poisson(rate = 0.3, tree = t) # simulate extant samples f = sim.extant.samples(f, tree = t, rho = 0.5) # plot the complete tree plot(f,t) # generate tree & fossil objects corresponding to the reconstructed tree out = reconstructed.tree.fossils.objects(f, t) f.reconst = out$fossils t.reconst = out$tree # plot the reconstructed tree plot(f.reconst, t.reconst)
Remove fossil samples that occur in the stem
remove.stem.fossils(fossils, tree)
remove.stem.fossils(fossils, tree)
fossils |
an object of class "fossils" that corresponds to fossil occurrences for "tree". |
tree |
an object of class "Phylo". |
an object of class "fossils", containing only the fossil samples that occur in the crown.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t) remove.stem.fossils(f, t)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t) remove.stem.fossils(f, t)
Remove stem lineages from a tree.
remove.stem.lineages(tree)
remove.stem.lineages(tree)
tree |
an object of class "Phylo". |
an object of class "Phylo", if stem lineages were found in the tree these will be pruned; if not then the original tree is returned.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] remove.stem.lineages(t)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] remove.stem.lineages(t)
Removes all unsampled lineages from a combined tree. Extinct tips are only sampled if they are fossils. With default settings all extant tips are sampled.
sampled.tree.from.combined(tree, rho = 1, sampled_tips = NULL)
sampled.tree.from.combined(tree, rho = 1, sampled_tips = NULL)
tree |
Combined tree with fossils. |
rho |
Sampling probability of extant tips. Default 1, will be disregarded if sampled_tips is not null. |
sampled_tips |
List of tip labels corresponding to sampled extant tips. |
Sampled tree with fossils.
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # transform format t2 = SAtree.from.fossils(t,f)$tree # transform to sampled tree t3 = sampled.tree.from.combined(t2) plot(t3)
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # transform format t2 = SAtree.from.fossils(t,f)$tree # transform to sampled tree t3 = sampled.tree.from.combined(t2) plot(t3)
Converts a phylo object to SAtree, without modification of tip labels.
SAtree(tree, complete = TRUE)
SAtree(tree, complete = TRUE)
tree |
Phylo object. |
complete |
Whether the tree is complete. Default TRUE. If the tree is not complete, then all fossil tips correspond to fossil samples, otherwise only sampled ancestors are considered samples. |
Transforms a tree and fossils dataframe to a combined SA tree. Sampled ancestors are represented as tips on zero-length edges to maintain compatibility with the ape format. Tip labels are set to "species id"_"index", where the most recent tip of a given species receives index 1 and indices increase towards the past.
SAtree.from.fossils(tree, fossils)
SAtree.from.fossils(tree, fossils)
tree |
Phylo object. |
fossils |
Fossils object. |
A list of 'tree', the SA tree integrating the fossils, and 'fossils', the fossils object updated with the tip label of each sample.
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # transform format t2 = SAtree.from.fossils(t,f) plot(t2$tree)
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(rate = 2, tree = t) # transform format t2 = SAtree.from.fossils(t,f) plot(t2$tree)
Simulate anagenetic species on a taxonomy object
sim.anagenetic.species(tree, species, lambda.a)
sim.anagenetic.species(tree, species, lambda.a)
tree |
Phylo object. |
species |
Taxonomy object. |
lambda.a |
Rate of anagenetic speciation. Default = 0. |
Object of class taxonomy.
t = ape::rtree(10) sp = sim.taxonomy(t, 1) sim.anagenetic.species(t, sp, 0.1)
t = ape::rtree(10) sp = sim.taxonomy(t, 1) sim.anagenetic.species(t, sp, 0.1)
Simulate cryptic species on a taxonomy object
sim.cryptic.species(species, kappa)
sim.cryptic.species(species, kappa)
species |
Taxonomy object. |
kappa |
Probability that speciation event is cryptic. |
An object of class taxonomy. Note the origin or root can not be cryptic.
t = ape::rtree(10) sp = sim.taxonomy(t, 1) sim.cryptic.species(sp, 0.5)
t = ape::rtree(10) sp = sim.taxonomy(t, 1) sim.cryptic.species(sp, 0.5)
Include extant samples in the fossil object, with optional rho sampling.
sim.extant.samples(fossils, tree = NULL, taxonomy = NULL, rho = 1, tol = NULL)
sim.extant.samples(fossils, tree = NULL, taxonomy = NULL, rho = 1, tol = NULL)
fossils |
Fossils object. |
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
rho |
Extant species sampling probability. Can be a single value or a vector. Vector entries will be applied to extant tips in the order in which they appear in the taxonomy object. |
tol |
Rounding error tolerance for tip ages. |
An object of class fossils containing extant tip samples equal to the age of the tips (i.e. 0.0).
# simulate tree lambda = 0.1 mu = 0.05 tips = 8 t = TreeSim::sim.bd.taxa(tips, 1, lambda, mu)[[1]] # simulate fossils f = sim.fossils.poisson(0.5, t) # simulate extant samples f = sim.extant.samples(f, t, rho = 0.5) plot(f, t)
# simulate tree lambda = 0.1 mu = 0.05 tips = 8 t = TreeSim::sim.bd.taxa(tips, 1, lambda, mu)[[1]] # simulate fossils f = sim.fossils.poisson(0.5, t) # simulate extant samples f = sim.extant.samples(f, t, rho = 0.5) plot(f, t)
sim.fbd.age: Simulating fossilized birth-death trees of a fixed age.
sim.fbd.age( age, numbsim, lambda, mu, psi, frac = 1, mrca = FALSE, complete = FALSE, K = 0 )
sim.fbd.age( age, numbsim, lambda, mu, psi, frac = 1, mrca = FALSE, complete = FALSE, K = 0 )
age |
Time since origin / most recent common ancestor. |
numbsim |
Number of trees to simulate. |
lambda |
Speciation rate. |
mu |
Extinction rate. |
psi |
Fossil sampling rate. |
frac |
Extant sampling fraction: The actual (simulated) number of tips is n, but only n*frac tips are included in the sampled tree (incomplete sampling). |
mrca |
If mrca=FALSE: age is the time since origin. If mrca=TRUE: age is the time since the most recent common ancestor of all sampled tips. |
complete |
whether to return the complete tree (with non-sampled lineages) or the reconstructed tree (with unsampled lineages removed). |
K |
If K = 0 (default), then lambda is constant. If K>0, density-dependent speciation is assumed, with speciation rate = lambda(1-m/K) when there are m living species. |
Array of 'numbsim' SAtrees with the time since origin / most recent common ancestor being 'age'. If the tree goes extinct or no tips are sampled (only possible when mrca = FALSE), return value is '0'. If only one extant and no extinct tips are sampled, return value is '1'.
age = 1 lambda = 2.0 mu = 0.5 psi = 0.6 numbsim = 2 if (requireNamespace("TreeSim", quietly = TRUE)) { sim.fbd.age(age, numbsim, lambda, mu, psi) }
age = 1 lambda = 2.0 mu = 0.5 psi = 0.6 numbsim = 2 if (requireNamespace("TreeSim", quietly = TRUE)) { sim.fbd.age(age, numbsim, lambda, mu, psi) }
sim.fbd.rateshift.taxa: Simulating fossilized birth death trees incorporating rate shifts.
sim.fbd.rateshift.taxa(n, numbsim, lambda, mu, psi, times, complete = FALSE)
sim.fbd.rateshift.taxa(n, numbsim, lambda, mu, psi, times, complete = FALSE)
n |
Number of extant sampled tips. |
numbsim |
Number of trees to simulate. |
lambda |
Vector of speciation rates, the rate in entry i is the speciation rate prior (ancestral) to time times[i]. |
mu |
Vector of extinction rates, the rate in entry i is the extinction rate prior (ancestral) to time times[i]. |
psi |
Vector of fossil sampling rates, the rate in entry i is the fossil sampling rate prior (ancestral) to time times[i]. |
times |
Vector of mass extinction and rate shift times. Time is 0 today and increasing going backwards in time. Specify the vector as times[i]<times[i+1]. times[1]=0 (today). |
complete |
whether to return the complete tree (with non-sampled lineages) or the reconstructed tree (with unsampled lineages removed). |
List of numbsim simulated SAtrees with n extant sampled tips.
n = 10 numbsim = 1 if (requireNamespace("TreeSim", quietly = TRUE)) { sim.fbd.rateshift.taxa(n, numbsim, lambda = c(2,1), mu = c(0,0.3), psi = c(1,0.1), times = c(0,0.3)) }
n = 10 numbsim = 1 if (requireNamespace("TreeSim", quietly = TRUE)) { sim.fbd.rateshift.taxa(n, numbsim, lambda = c(2,1), mu = c(0,0.3), psi = c(1,0.1), times = c(0,0.3)) }
sim.fbd.taxa: Simulating fossilized birth-death trees on a fixed number of extant taxa.
sim.fbd.taxa(n, numbsim, lambda, mu, psi, frac = 1, complete = FALSE)
sim.fbd.taxa(n, numbsim, lambda, mu, psi, frac = 1, complete = FALSE)
n |
Number of extant sampled tips. |
numbsim |
Number of trees to simulate. |
lambda |
Speciation rate. |
mu |
Extinction rate. |
psi |
Fossil sampling rate. |
frac |
Extant sampling fraction. When complete = FALSE, the actual (simulated) number of extant tips is n/frac, but only n tips are included in the result (incomplete sampling). When complete = TRUE: all unsampled lineages are included, i.e. the final tree has n/frac extant tips. |
complete |
whether to return the complete tree (with non-sampled lineages) or the reconstructed tree (with unsampled lineages removed). |
List of numbsim simulated SAtrees with n extant sampled tips.
n = 10 lambda = 2.0 mu = 0.5 psi = 0.6 numbsim = 2 if (requireNamespace("TreeSim", quietly = TRUE)) { sim.fbd.taxa(n, numbsim, lambda, mu, psi) }
n = 10 lambda = 2.0 mu = 0.5 psi = 0.6 numbsim = 2 if (requireNamespace("TreeSim", quietly = TRUE)) { sim.fbd.taxa(n, numbsim, lambda, mu, psi) }
This function uses a three parameter Gaussian model to simulate non-uniform fossil recovery along a specified phylogeny. Preservation varies with respect to water depth, which is a useful for proxy for changes in the depositional environment. The per interval probability of sampling is
where PA is species peak abundance, PD is preferred depth, DT is depth tolerance and d is current water depth.
PD is the depth at which the species is most likely to be found and is equivalent to the mean of the distribution.
PA is the probability of sampling an occurrence at this depth.
DT is the potential of a species to be found at a range of depths and is equivalent to the standard deviation.
Although here fossil recovery is described with respect to water depth, the model could be applied in the context of any environmental gradient.
The model uses a probability of collecting a fossil within a given interval, rather than a rate.
To simulate discrete fossil sampling events and times within each interval we need to convert the probability into a rate
(use.rates = TRUE
). This is done using the formula
where t is the interval length.
One caveat of this approach is that the model cannot use a probability of 1, as it would correspond to rate = infinity.
In this instance we use an approximation for probabilities = 1 (e.g. pr.1.approx = 0.999
).
Non-uniform interval ages can be specified as a vector (interval.ages
) or a uniform set of interval ages can be specified using
maximum interval age (max.age
) and the number of intervals (strata
), where interval length .
A vector of values can be specified for the model parameters PA, PD and DT to allow for variation across lineages.
If a vector is provided, each entry will apply to each unique species in the order in which they appear in the taxonomy object (if taxonomy is provided),
or to each unique edge in the order in which they appear in the tree object.
If the tree object has a root edge (root.edge
), the first entry in the vector will apply to this edge.
Fossils can be simulated for a phylo (tree
) or taxonomy (taxonomy
) object.
If both are specified, the function uses taxonomy.
If no taxonomic information is provided, the function assumes all speciation is symmetric (i.e. bifurcating, beta = 1
).
sim.fossils.environment( tree = NULL, taxonomy = NULL, interval.ages = NULL, max.age = NULL, strata = NULL, proxy.data = NULL, PD = 0.5, DT = 0.5, PA = 0.5, root.edge = TRUE, use.rates = FALSE, pr.1.approx = 0.999, use.exact.times = TRUE )
sim.fossils.environment( tree = NULL, taxonomy = NULL, interval.ages = NULL, max.age = NULL, strata = NULL, proxy.data = NULL, PD = 0.5, DT = 0.5, PA = 0.5, root.edge = TRUE, use.rates = FALSE, pr.1.approx = 0.999, use.exact.times = TRUE )
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
interval.ages |
Vector of stratigraphic interval ages, starting with the minimum age of the youngest interval and ending with the maximum age of the oldest interval. |
max.age |
Maximum age of the oldest stratigraphic interval or age at the base of the basin. |
strata |
Number of stratigraphic intervals. |
proxy.data |
Vector of relative water depth or other proxy data. The first number corresponds to the youngest interval. The length of the vector should be 1 less than the length of interval.ages. |
PD |
Preferred depth parameter value or a vector of values. |
DT |
Depth tolerance parameter value or a vector of values. |
PA |
Peak abundance parameter value or a vector of values. |
root.edge |
If TRUE include the root edge. Default = TRUE. |
use.rates |
If TRUE convert per interval sampling probability into a per interval Poisson rate. Default = FALSE. |
pr.1.approx |
Value used to approximate sampling probabilities = 1 when use.rates = TRUE. |
use.exact.times |
If TRUE use exact sampling times. If FALSE |
An object of class fossils.
Holland, S.M. 1995. The stratigraphic distribution of fossils. Paleobiology 21: 92-109.
sim.fossils.poisson
, sim.fossils.intervals
, sim.trait.values
# simulate tree t = ape::rtree(6) # assign a max age based on tree height max.age = tree.max(t) # generate water depth profile strata = 7 wd = sim.gradient(strata) # simulate fossils using tree & max.age and strata f = sim.fossils.environment(t, max.age = max.age, strata = strata, proxy.data = wd, PD = 0.5, DT = 1, PA = 1) plot(f, t, show.proxy = TRUE, proxy.data = wd, strata = strata, show.strata = TRUE) # simulate fossils using taxonomy & interval.ages s = sim.taxonomy(t, 0.1, 0.1, 1) times = seq(0, max.age, length.out = strata + 1) f = sim.fossils.environment(taxonomy = s, interval.ages = times, proxy.data = wd, PD = 0.5, DT = 1, PA = 1) plot(f, t, strata = strata, binned = TRUE) # simulate fossils with variable preservation across lineages dist = function() {runif(1)} PD = sim.trait.values(1, taxonomy = s, model = "independent", dist = dist, change.pr = 0.1) f = sim.fossils.environment(taxonomy = s, interval.ages = times, proxy.data = wd, PD = PD, DT = 1, PA = 1) plot(f, t, strata = strata, binned = TRUE)
# simulate tree t = ape::rtree(6) # assign a max age based on tree height max.age = tree.max(t) # generate water depth profile strata = 7 wd = sim.gradient(strata) # simulate fossils using tree & max.age and strata f = sim.fossils.environment(t, max.age = max.age, strata = strata, proxy.data = wd, PD = 0.5, DT = 1, PA = 1) plot(f, t, show.proxy = TRUE, proxy.data = wd, strata = strata, show.strata = TRUE) # simulate fossils using taxonomy & interval.ages s = sim.taxonomy(t, 0.1, 0.1, 1) times = seq(0, max.age, length.out = strata + 1) f = sim.fossils.environment(taxonomy = s, interval.ages = times, proxy.data = wd, PD = 0.5, DT = 1, PA = 1) plot(f, t, strata = strata, binned = TRUE) # simulate fossils with variable preservation across lineages dist = function() {runif(1)} PD = sim.trait.values(1, taxonomy = s, model = "independent", dist = dist, change.pr = 0.1) f = sim.fossils.environment(taxonomy = s, interval.ages = times, proxy.data = wd, PD = PD, DT = 1, PA = 1) plot(f, t, strata = strata, binned = TRUE)
Intervals can be specified by specifying the interval boundaries using interval.ages
or specifying both max.age
and strata
.
In the second scenario all intervals will be of equal length.
Preservation can be specified using rates
, which represent the rates of a Poisson process in each interval,
or probabilities
, which represent the probabilities of sampling per interval.
When using probabilities
, at most one fossil per species will be sampled per interval.
Fossils can be simulated for a phylo (tree
) or taxonomy (taxonomy
) object.
If both are specified, the function uses taxonomy.
If no taxonomic information is provided, the function assumes all speciation is symmetric (i.e. bifurcating, beta = 1
).
sim.fossils.intervals( tree = NULL, taxonomy = NULL, fossils = NULL, interval.ages = NULL, max.age = NULL, strata = NULL, probabilities = NULL, rates = NULL, ignore.taxonomy = FALSE, root.edge = TRUE, use.exact.times = TRUE )
sim.fossils.intervals( tree = NULL, taxonomy = NULL, fossils = NULL, interval.ages = NULL, max.age = NULL, strata = NULL, probabilities = NULL, rates = NULL, ignore.taxonomy = FALSE, root.edge = TRUE, use.exact.times = TRUE )
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
fossils |
Append fossils to to an existing fossils object. |
interval.ages |
Vector of stratigraphic interval ages, starting with the minimum age of the youngest interval and ending with the maximum age of the oldest interval. |
max.age |
Maximum age of the oldest stratigraphic interval. |
strata |
Number of stratigraphic intervals. |
probabilities |
Probability of sampling/preservation in each interval. The number of probabilities should match the number of intervals and the first entry should correspond to youngest interval. |
rates |
Poisson sampling rate for each interval. The number of rates should match the number of intervals and the first entry should correspond to youngest interval. |
ignore.taxonomy |
Ignore species taxonomy (returns sp = NA). Default = FALSE. |
root.edge |
If TRUE include the root edge. Default = TRUE. |
use.exact.times |
If TRUE use exact sampling times. If FALSE |
An object of class fossils.
sim.fossils.poisson
, sim.fossils.environment
# simulate tree t = ape::rtree(6) # assign a max age based on tree height max.age = tree.max(t) # simulate fossils using max.age and strata & probabilities strata = 4 probability = rep(0.7, 4) f = sim.fossils.intervals(t, max.age = max.age, strata = strata, probabilities = probability) plot(f, t, strata = strata, show.strata = TRUE) # simulate fossils using interval.ages & rates times = c(0, sort(runif(3, min = 0, max = max.age)), max.age) rates = c(5, 0, 5, 0) f = sim.fossils.intervals(t, interval.ages = times, rates = rates) plot(f, t, interval.ages = times, show.strata = TRUE) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.1, 0.1, 1) f = sim.fossils.intervals(taxonomy = s, interval.ages = times, rates = rates) plot(f, t, interval.ages = times, show.strata = TRUE) # append fossils to an existing fossils object new.rates = rates * 2 f2 = sim.fossils.intervals(taxonomy = s, fossils = f, interval.ages = times, rates = new.rates)
# simulate tree t = ape::rtree(6) # assign a max age based on tree height max.age = tree.max(t) # simulate fossils using max.age and strata & probabilities strata = 4 probability = rep(0.7, 4) f = sim.fossils.intervals(t, max.age = max.age, strata = strata, probabilities = probability) plot(f, t, strata = strata, show.strata = TRUE) # simulate fossils using interval.ages & rates times = c(0, sort(runif(3, min = 0, max = max.age)), max.age) rates = c(5, 0, 5, 0) f = sim.fossils.intervals(t, interval.ages = times, rates = rates) plot(f, t, interval.ages = times, show.strata = TRUE) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.1, 0.1, 1) f = sim.fossils.intervals(taxonomy = s, interval.ages = times, rates = rates) plot(f, t, interval.ages = times, show.strata = TRUE) # append fossils to an existing fossils object new.rates = rates * 2 f2 = sim.fossils.intervals(taxonomy = s, fossils = f, interval.ages = times, rates = new.rates)
Fossils can be simulated for a phylo (tree
) or taxonomy (taxonomy
) object.
If both are specified, the function uses taxonomy.
If no taxonomic information is provided, the function assumes all speciation is symmetric (i.e. bifurcating, beta = 1
).
A vector of rates can be specified to allow for rate variation across lineages.
If a vector is provided, each entry will apply to each unique species in the order in which they appear in the taxonomy object (if taxonomy is provided),
or to each unique edge in the order in which they appear in the tree object.
If the tree object has a root edge (root.edge
), the first entry in the rates vector should correspond to this edge.
sim.fossils.poisson( rate, tree = NULL, taxonomy = NULL, fossils = NULL, ignore.taxonomy = FALSE, root.edge = TRUE )
sim.fossils.poisson( rate, tree = NULL, taxonomy = NULL, fossils = NULL, ignore.taxonomy = FALSE, root.edge = TRUE )
rate |
A single Poisson sampling rate or a vector of rates. |
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
fossils |
Append fossils to to an existing fossils object. |
ignore.taxonomy |
Ignore species taxonomy (returns sp = NA). Default = FALSE. |
root.edge |
If TRUE include the root edge. Default = TRUE. |
An object of class fossils.
sim.fossils.intervals
, sim.fossils.environment
, sim.trait.values
# simulate tree t = ape::rtree(6) # simulate fossils using the tree rate = 2 f = sim.fossils.poisson(rate, tree = t) plot(f, t) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) f = sim.fossils.poisson(rate, taxonomy = s) plot(f, t) # simulate fossils with autocorrelated rate variation across lineages rates = sim.trait.values(init = rate, taxonomy = s, v = 1) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # append fossils to an existing fossils object rate = 1 f1 = sim.fossils.poisson(rate, tree = t) plot(f1, t) rate = 2 f2 = sim.fossils.poisson(rate, tree = t, fossils = f1) plot(f2, t) f3 = sim.fossils.poisson(rate, tree = t, fossils = f2, ignore.taxonomy = TRUE) plot(f3, t, show.unknown = TRUE)
# simulate tree t = ape::rtree(6) # simulate fossils using the tree rate = 2 f = sim.fossils.poisson(rate, tree = t) plot(f, t) # simulate fossils using taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) f = sim.fossils.poisson(rate, taxonomy = s) plot(f, t) # simulate fossils with autocorrelated rate variation across lineages rates = sim.trait.values(init = rate, taxonomy = s, v = 1) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # append fossils to an existing fossils object rate = 1 f1 = sim.fossils.poisson(rate, tree = t) plot(f1, t) rate = 2 f2 = sim.fossils.poisson(rate, tree = t, fossils = f1) plot(f2, t) f3 = sim.fossils.poisson(rate, tree = t, fossils = f2, ignore.taxonomy = TRUE) plot(f3, t, show.unknown = TRUE)
Function returns a vector using the sine wave function
for a given set of intervals.
This vector can be used as a gradient to simulate fossils under an environment-dependent model of fossil recovery using the
function
sim.fossils.environment
.
sim.gradient(strata, depth = 2, cycles = 2)
sim.gradient(strata, depth = 2, cycles = 2)
strata |
Number of stratigraphic intervals. |
depth |
Maximum water depth. |
cycles |
Number of cycles (e.g. transgressions and regressions). |
vector of sampled water depths.
strata = 100 wd = sim.gradient(strata) plot(wd, type="l")
strata = 100 wd = sim.gradient(strata) plot(wd, type="l")
Reassign exact fossil ages using the minimum and maximum ages of a set of stratigraphic intervals.
If use.species.ages = TRUE
the function will respect species durations and will not
return minimum and maximum bounds that may be younger or older than the species durations.
This requires supplying a phylo or taxonomy object.
sim.interval.ages( fossils, tree = NULL, taxonomy = NULL, interval.ages = NULL, max.age = NULL, strata = NULL, use.species.ages = FALSE, root.edge = TRUE, sim.extant = FALSE )
sim.interval.ages( fossils, tree = NULL, taxonomy = NULL, interval.ages = NULL, max.age = NULL, strata = NULL, use.species.ages = FALSE, root.edge = TRUE, sim.extant = FALSE )
fossils |
Fossil object. |
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
interval.ages |
Vector of stratigraphic interval ages, starting with the minimum age of the youngest interval and ending with the maximum age of the oldest interval. |
max.age |
Maximum age of the oldest stratigraphic interval. |
strata |
Number of stratigraphic intervals. |
use.species.ages |
If TRUE reassigned fossil ages will respect the speciation times. Default = FALSE. |
root.edge |
If TRUE include root edge. |
sim.extant |
If TRUE simulate age uncertainty for extant samples as well, default FALSE. |
An object of class fossils.
# simulate tree t = ape::rtree(8) # simulate fossils rate = 2 f = sim.fossils.poisson(rate, t) plot(f, t) # assign a max age based on tree height max.age = tree.max(t) # define intervals times = seq(0, max.age, length.out = 5) # reassign ages f = sim.interval.ages(f, t, interval.ages = times) # plot output plot(f, t, interval.ages = times)
# simulate tree t = ape::rtree(8) # simulate fossils rate = 2 f = sim.fossils.poisson(rate, t) plot(f, t) # assign a max age based on tree height max.age = tree.max(t) # define intervals times = seq(0, max.age, length.out = 5) # reassign ages f = sim.interval.ages(f, t, interval.ages = times) # plot output plot(f, t, interval.ages = times)
Simulate a taxonomy object relating species identity to a phylo object under a mixed model of speciation.
Anagenetic and cryptic species can also be added later using the sim.anagenetic.species
and sim.cryptic.species
functions.
sim.taxonomy(tree, beta = 0, lambda.a = 0, kappa = 0, root.edge = TRUE)
sim.taxonomy(tree, beta = 0, lambda.a = 0, kappa = 0, root.edge = TRUE)
tree |
Phylo object. |
beta |
Probability of bifurcating speciation. Default = 0. |
lambda.a |
Rate of anagenetic speciation. Default = 0. |
kappa |
Probability that speciation event is cryptic. Default = 0. |
root.edge |
If TRUE include root edge. Default = TRUE. |
An object of class taxonomy.
t = ape::rtree(10) sim.taxonomy(t, 0.5, 0.1, 0.5)
t = ape::rtree(10) sim.taxonomy(t, 0.5, 0.1, 0.5)
Include extant and extinct tip samples in the fossil object, with optional rho sampling.
sim.tip.samples(fossils, tree, taxonomy = NULL, rho = 1)
sim.tip.samples(fossils, tree, taxonomy = NULL, rho = 1)
fossils |
Fossils object. |
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
rho |
Tip sampling probability. Can be a single value or a vector. Vector entries will be applied to extant tips in the order in which they appear in the taxonomy object. |
An object of class fossils containing extant or extinct tip samples equal to the age of the tips.
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(2, t) # simulate tip samples f = sim.tip.samples(f, t, rho = 0.5) plot(f, t)
# simulate tree t = ape::rtree(6) # simulate fossils f = sim.fossils.poisson(2, t) # simulate tip samples f = sim.tip.samples(f, t, rho = 0.5) plot(f, t)
Fossil recovery rates or other parameter values can be simulated for a phylo (tree
) or taxonomy (taxonomy
) object.
Under the autocorrelated
model, trait values evolve along lineages according to a Brownian motion process, where the strength of the relationship between ancestor and descendant values is determined by the parameter (
v
).
If is small values will be more similar between ancestor and descendants, and if
is zero all trait values will be equal.
For a given species
with ancestor
, a new trait value
is drawn from a lognormal distribution with
where and
is the lineage duration of the species.
This fossil recovery model is described in Heath et al. (2014) and is equivalent to the autocorrelated relaxed clock model described in Kishino et al. (2001).
Under the
BM
and OU
models traits are simulated under a standard Brownian motion or Ornstein-Uhlenbeck process with rate parameter (
v
).
The OU model has the additional parameter alpha
, which determines the strength with which trait values are attracted to the mean. Note the init
argument will specify both the value at the root and the mean of the process under the OU model.
Under the independent
model a new trait value is drawn for each species from any valid user-specified distribution (dist
).
change.pr
is the probability that a trait value will change at each speciation event.
If change.pr = 1
trait values will be updated at every speciation events.
Finally, traits can be simulated under the standard Lewis Mk model (Mk
), with symmetric rates of change. The rate is specified using v
and number of states using k
.
sim.trait.values( init = 1, tree = NULL, taxonomy = NULL, root.edge = TRUE, model = "autocorrelated", v = 0.01, alpha = 0.1, min.value = -Inf, max.value = Inf, dist = function() { runif(1, 0, 2) }, change.pr = 1, k = 2 )
sim.trait.values( init = 1, tree = NULL, taxonomy = NULL, root.edge = TRUE, model = "autocorrelated", v = 0.01, alpha = 0.1, min.value = -Inf, max.value = Inf, dist = function() { runif(1, 0, 2) }, change.pr = 1, k = 2 )
init |
Initial value at the origin or root of the phylo or taxonomy object. Default = 1. |
tree |
Phylo object. |
taxonomy |
Taxonomy object. |
root.edge |
If TRUE include the root edge. Default = TRUE. |
model |
Model used to simulate rate variation across lineages. Options include "autocorrelated" (default), "BM" (Brownian motion), "OU" (Ornstein-Uhlenbeck), "independent" or the Lewis "Mk" model. |
v |
Brownian motion parameter |
alpha |
Ornstein-Uhlenbeck parameter |
min.value |
Min trait value allowed under the BM and OU models. Default = -Inf. |
max.value |
Max trait value allowed under the BM and OU models. Default = Inf. |
dist |
Distribution of trait values used to draw new values under the "independent" model. This parameter is ignored if |
change.pr |
Probability that trait values change at speciation events. Default = 1. |
k |
Number of states used for the Mk model. Default = 2. |
A vector of parameter values.
Values are output for each species in the order in which they appear in the taxonomy object (if taxonomy was provided) or for each edge in the order in which they appear in the tree object.
If the tree object has a root edge (root.edge
), the first entry in the vector will correspond to this edge.
Heath et al. 2014. The fossilized birth-death process for coherent calibration of divergence-time estimates. PNAS 111:E2957-E2966.
Kishino et al. 2001. Performance of a divergence time estimation method under a probabilistic model of rate evolution MBE 18:352-361.
# simulate tree t = ape::rtree(6) # simulate taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) # simulate rates under the autocorrelated trait values model rate = 2 rates = sim.trait.values(rate, taxonomy = s, v = 1) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # simulate rates under the independent trait values model dist = function() { rlnorm(1, log(rate), 1) } rates = sim.trait.values(rate, taxonomy = s, model = "independent", dist = dist) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # simulate rates under the independent trait values model with infrequent changes rates = sim.trait.values(rate, taxonomy = s, model = "independent", dist = dist, change.pr = 0.1) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # simulate traits under Brownian motion and convert into rates traits = sim.trait.values(0, taxonomy = s, model = "BM", v = 2) # function for translating states into rates translate.states = function(traits, low, high) sapply(traits, function(t) if(t < 0) low else high) # sampling rates low = 0.1 high = 2 rates = translate.states(traits, low, high) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, tree = t)
# simulate tree t = ape::rtree(6) # simulate taxonomy s = sim.taxonomy(t, 0.5, 1, 0.5) # simulate rates under the autocorrelated trait values model rate = 2 rates = sim.trait.values(rate, taxonomy = s, v = 1) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # simulate rates under the independent trait values model dist = function() { rlnorm(1, log(rate), 1) } rates = sim.trait.values(rate, taxonomy = s, model = "independent", dist = dist) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # simulate rates under the independent trait values model with infrequent changes rates = sim.trait.values(rate, taxonomy = s, model = "independent", dist = dist, change.pr = 0.1) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, t) # simulate traits under Brownian motion and convert into rates traits = sim.trait.values(0, taxonomy = s, model = "BM", v = 2) # function for translating states into rates translate.states = function(traits, low, high) sapply(traits, function(t) if(t < 0) low else high) # sampling rates low = 0.1 high = 2 rates = translate.states(traits, low, high) f = sim.fossils.poisson(rates, taxonomy = s) plot(f, tree = t)
Find a species' end (i.e extinction) time from a taxonomy object
species.end(species, taxonomy)
species.end(species, taxonomy)
species |
Species id (as written in |
taxonomy |
Taxonomy object. |
End time.
Find a species' start (i.e speciation) time from a taxonomy object
species.start(species, taxonomy)
species.start(species, taxonomy)
species |
Species id (as written in |
taxonomy |
Taxonomy object. |
Start time.
Obtain a subsample of fossil occurrences containing the oldest fossil sample in each node of the tree.
subsample.fossils.oldest(fossils, tree, complete = TRUE)
subsample.fossils.oldest(fossils, tree, complete = TRUE)
fossils |
an object of class "fossils" that corresponds to fossil occurrences. |
tree |
an object of class "Phylo", representing the tree upon which the fossil occurrences were simulated. |
complete |
logical, if TRUE the oldest sample from each clade in the complete tree is returned, if FALSE the oldest sample from each clade in the extant only counterpart tree is returned. |
an object of class "fossils" containing the subsampled fossil occurrences.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.oldest(f, t, complete = FALSE)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.oldest(f, t, complete = FALSE)
Obtain a subsample of fossil occurrences containing the oldest and youngest fossil sample found at each node of the tree.
subsample.fossils.oldest.and.youngest(fossils, tree, complete = TRUE)
subsample.fossils.oldest.and.youngest(fossils, tree, complete = TRUE)
fossils |
an object of class "fossils" that corresponds to fossil occurrences. |
tree |
an object of class "Phylo", representing the tree upon which the fossil occurrences were simulated. |
complete |
logical, if TRUE the oldest and youngest sample from each clade in the complete tree is returned, if FALSE the oldest and youngest sample from each clade in the extant only counterpart tree is returned. |
an object of class "fossils" containing the subsampled fossil occurrences.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.oldest.and.youngest(f, t, complete = FALSE)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.oldest.and.youngest(f, t, complete = FALSE)
Obtain a uniform random sample of fossil occurrences.
subsample.fossils.uniform(fossils, proportion)
subsample.fossils.uniform(fossils, proportion)
fossils |
an object of class "fossils" that corresponds to fossil occurrences. |
proportion |
the proportion of all fossil samples to return in the subsample. |
an object of class "fossils" containing the subsampled fossil occurrences.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.uniform(f, 0.5)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.uniform(f, 0.5)
Obtain a subsample of fossil occurrences containing the youngest fossil sample in each node of the tree.
subsample.fossils.youngest(fossils, tree, complete = TRUE)
subsample.fossils.youngest(fossils, tree, complete = TRUE)
fossils |
an object of class "fossils" that corresponds to fossil occurrences. |
tree |
an object of class "Phylo", representing the tree upon which the fossil occurrences were simulated. |
complete |
logical, if TRUE the youngest sample from each clade in the complete tree is returned, if FALSE the youngest sample from each clade in the extant only counterpart tree is returned. |
an object of class "fossils" containing the subsampled fossil occurrences.
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.youngest(f, t, complete = FALSE)
t = TreeSim::sim.bd.taxa(10, 1, 0.1, 0.05)[[1]] f = sim.fossils.poisson(0.1, t, root.edge = FALSE) subsample.fossils.youngest(f, t, complete = FALSE)
Display taxonomy object
## S3 method for class 'taxonomy' summary(object, max.length = 50, round.x = 12, details = TRUE, ...)
## S3 method for class 'taxonomy' summary(object, max.length = 50, round.x = 12, details = TRUE, ...)
object |
Taxonomy object. |
max.length |
Max number of rows to print out. |
round.x |
Number of decimal places to be used for species and edge ages. |
details |
If TRUE include summary statistics. |
... |
Additional parameters to be passed to |
Create a taxonomy object relating species identity to a phylo object.
taxonomy(data) as.taxonomy(data) is.taxonomy(data)
taxonomy(data) as.taxonomy(data) is.taxonomy(data)
data |
Dataframe of species taxonomy. See Details for the list of required fields. |
The taxonomy object includes the following 6 fields for each edge in the corresponding phylo object:
sp
true species identity label.
If all species originated via budding or bifurcation this will always correspond to the terminal-most edge label (i.e. the youngest node) associated with each species.
This is not the case if the data set also contains anagenetic species, when multiple species may be associated with a single edge
edge
edge label of the branch in the corresponding phylo object.
Note some species may be associated with multiple edges
parent
= ancestor of species sp
. Parent labels follow the same convention as species.
The label assigned to the parent of the origin or root will be zero
start
= start time of the corresponding edge
and/or origin time of the species.
If the corresponding edge is also the oldest edge associated with the species this value will equal the species origination time.
If speciation mode is asymmetric or symmetric the speciation time will match the start time of the corresponding edge.
If speciation mode is anagenetic the speciation time will be younger than the start time of the corresponding edge
end
= end time of the corresponding edge
and/or end time of the species.
If the corresponding edge is also the youngest edge associated with the species this value will equal the species end time.
Unless the species end time coincides with an anagenetic speciation event, the species end time will match the end time of the corresponding edge.
If the species end time coincides with an anagenetic speciation event, the speciation time will be older than the end time of the corresponding edge
mode
= speciation mode. "o" = origin or "r" = root (the edge/species that began the process).
"b" = asymmetric or budding speciation. "s" = symmetric or bifurcating speciation. "a" = anagenetic speciation
Optional fields:
cryptic
TRUE if the speciation event was cryptic. If missing the function assumes cryptic = FALSE
cryptic.id
= cryptic species identity. If cryptic = TRUE
cryptic.id
will differ from the true species identity sp
Function returns the the root age or the origin time (if root.edge = TRUE
).
tree.max(tree, root.edge = TRUE)
tree.max(tree, root.edge = TRUE)
tree |
Phylo object. |
root.edge |
If TRUE include the root edge (default = TRUE). |
max age
t = ape::rtree(6) tree.max(t, root.edge = FALSE)
t = ape::rtree(6) tree.max(t, root.edge = FALSE)