bd.sim.constant: Constant rate Birth-Death simulation

Description Usage Arguments Value Author(s) Examples

View source: R/bd.sim.constant.R

Description

Simulates a species birth-death process with constant rates for any number of starting species. Allows for constraining results on the number of species at the end of the simulation, either total or extant. Returns an object containing vectors of speciation times, extinction times, parents (= species' mother species) and status at the end of the simulation (extant or not) for each species in the simulation. It may return true extinction times or simply information on whether species lived after the maximum simulation time. For time-varying and age-dependent simulations, see bd.sim.general. For a function that unites both scenarios, see bd.sim. Please note while time runs from 0 to tMax in the simulation, it returns speciation/extinction times as tMax (origin of the group) to 0 (the "present" and end of simulation), so as to conform to other packages in the literature.

Usage

1
2
3
4
5
6
7
8
9
bd.sim.constant(
  n0,
  pp,
  qq,
  tMax,
  nFinal = c(0, Inf),
  extOnly = FALSE,
  trueExt = FALSE
)

Arguments

n0

Initial number of species. Usually 1, in which case the simulation will describe the full diversification of a monophyletic lineage. Notice that when pp is less than or equal to qq, many simulations will go extinct before speciating even once. One way of generating simulations in this case is to increase n0. The user should notice this will simulate the diversification of a paraphyletic group.

pp

Speciation rate. Must be constant.

qq

Extinction rate, similar to above.

Note: pp and qq must always be greater than 0

tMax

Ending time of simulation, in million years after the clade origin. Any species still living after tMax is considered extant, and any species that would be generated after tMax is not born.

nFinal

A vector of length 2, indicating an interval of acceptable number of species at the end of the simulation. Default value is c(0, Inf), so that any number of species (including zero, the extinction of the whole clade) is accepted. If different from default value, the process will run until the number of total species reaches a number in the interval nFinal. Please note that extOnly can modify the meaning of this parameter. If extOnly = TRUE, the function will run repeatedly until the number of extant species at the end of the simulation lies within the nFinal interval. Note that using values other than the default for nFinal might condition simulation results.

extOnly

A logical indicating whether nFinal should be taken as an interval of the number of total or extant species during the simulation. If TRUE, the function will run repeatedly until the number of extant species at the end of the simulation lies within the nFinal interval. If FALSE (as default), it will run until the total number of species generated lies within that interval.

Note: The function returns NA if it runs for more than 100000 iterations without fulfilling the requirements of nFinal and extOnly.

trueExt

A logical used for bd.sim.general, indicating whether it should return true or truncated extinction times. When TRUE, time of extinction of extant species will be the true time, otherwise it will be NA if a species is alive at the end of the simulation.

Value

A list of vectors, as follows

TE

List of extinction times, with 0 as the time of extinction for extant species.

TS

List of speciation times, with NA as the time of speciation for species that started the simulation.

PAR

List of parents. Species that started the simulation have NA, while species that were generated during the simulation have their parent's number. Species are numbered as they are born.

EXTANT

List of logicals representing whether each species is extant.

Author(s)

Bruno do Rosario Petrucci.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
# we can show a couple scenarios

###
# first, extinction 0

# initial number of species
n0 <- 1

# maximum simulation time
tMax <- 40

# speciation
p <- 0.1

# extinction
q <- 0

# run the simulation
sim <- bd.sim.constant(n0, p, q, tMax, nFinal = c(2, Inf))

# we can plot the phylogeny to take a look
if (requireNamespace("ape", quietly = TRUE)) {
  phy <- make.phylo(sim)
  ape::plot.phylo(phy)
}

###
# now let us try to turn extinction up a bit

# initial number of species
n0 <- 1

# maximum simulation time
tMax <- 40

# speciation
p <- 0.1

# extinction
q <- 0.04

# run the simulation
sim <- bd.sim.constant(n0, p, q, tMax, nFinal = c(2, Inf))

# we can plot the phylogeny to take a look
if (requireNamespace("ape", quietly = TRUE)) {
  phy <- make.phylo(sim)
  ape::plot.phylo(phy)
}

###
# we can also try a pure-death process

# initial number of species - note the high number, so we get an appreciable
# sample size
n0 <- 100

# maximum simulation time
tMax <- 40

# speciation
p <- 0

# extinction
q <- 0.02

# run the simulation
sim <- bd.sim.constant(n0, p, q, tMax)

# of course in this case there are no phylogenies to plot

# note nFinal has to be sensible
## Not run: 
# this would return a warning, since it is virtually impossible to get 100
# species at a process with diversification rate -0.09 starting at n0 = 1
sim <- bd.sim.constant(1, pp = 0.01, qq = 1, tMax = 100, nFinal = c(100, Inf))

## End(Not run)

brpetrucci/paleobuddy documentation built on Aug. 8, 2020, 2:03 a.m.