simulate_deterministic_hbd {castor} | R Documentation |

Given a homogenous birth-death (HBD) model, i.e., with speciation rate *λ*, extinction rate *μ* and sampling fraction *ρ*, calculate various deterministic features of the model backwards in time, such as the total diversity over time. The speciation and extinction rates may be time-dependent. “Homogenous” refers to the assumption that, at any given moment in time, all lineages exhibit the same speciation/extinction rates (in the literature this is sometimes referred to simply as “birth-death model”; Morlon et al. 2011). “Deterministic” refers to the fact that all calculated properties are completely determined by the model's parameters (i.e. non-random), as if an infinitely large tree was generated (aka. “continuum limit”).

Alternatively to *λ*, one may provide the pulled diversification rate (PDR; Louca et al. 2018) and the speciation rate at some fixed age, *λ(τ_o)*. Similarly, alternatively to *μ*, one may provide the ratio of extinction over speciation rate, *μ/λ*. In either case, the time-profiles of *λ*, *μ*, *μ/λ* or the PDR are specified as piecewise polynomial functions (splines), defined on a discrete grid of ages.

simulate_deterministic_hbd(LTT0, oldest_age, age0 = 0, rho0 = 1, age_grid = NULL, lambda = NULL, mu = NULL, mu_over_lambda = NULL, PDR = NULL, lambda0 = NULL, splines_degree = 1, relative_dt = 1e-3, allow_unreal = FALSE)

`LTT0` |
The assumed number of sampled extant lineages at |

`oldest_age` |
Strictly positive numeric, specifying the oldest time before present (“age”) to include in the simulation. |

`age0` |
Non-negative numeric, specifying the age at which |

`rho0` |
Numeric between 0 (exclusive) and 1 (inclusive), specifying the sampling fraction of the tree at |

`age_grid` |
Numeric vector, listing discrete ages (time before present) on which either |

`lambda` |
Numeric vector, of the same size as |

`mu` |
Numeric vector, of the same size as |

`mu_over_lambda` |
Numeric vector, of the same size as |

`PDR` |
Numeric vector, of the same size as |

`lambda0` |
Non-negative numeric, specifying the speciation rate (in units 1/time) at |

`splines_degree` |
Integer, either 0,1,2 or 3, specifying the polynomial degree of the provided |

`relative_dt` |
Strictly positive numeric (unitless), specifying the maximum relative time step allowed for integration over time. Smaller values increase integration accuracy but increase computation time. Typical values are 0.0001-0.001. The default is usually sufficient. |

`allow_unreal` |
Logical, specifying whether HBD models with unrealistic parameters (e.g., negative |

This function supports the following alternative parameterizations of HBD models:

Using the speciation rate

*λ*and extinction rate*μ*.Using the speciation rate

*λ*and the ratio*μ/λ*.Using the pulled diversification rate (PDR), the extinction rate and the speciation rate given at some fixed

`age0`

(i.e.`lambda0`

).Using the PDR, the ratio

*μ/λ*and the speciation rate at some fixed`age0`

.

The PDR is defined as *PDR = λ-μ+λ^{-1}dλ/dτ*, where *τ* is age (time before present). To avoid ambiguities, only one of the above parameterizations is accepted at a time. The sampling fraction at `age0`

(i.e., `rho0`

) should always be provided; setting it to `NULL`

is equivalent to setting it to 1.

Note that in the literature the sampling fraction usually refers to the fraction of lineages extant at present-day that have been sampled (included in the tree); this present-day sampling fraction is then used to parameterize birth-death cladogenic models. The sampling fraction can however be generalized to past times, by defining it as the probability that a lineage extant at any given time is included in the tree. The simulation function presented here allows for specifying this generalized sampling fraction at any age of choice, not just present-day.

The simulated LTT refers to a hypothetical tree sampled at age `age_grid[1]`

, i.e. LTT(t) will be the number of lineages extant at age t that survived until age `age_grid[1]`

and have been sampled, given that the fraction of sampled extant lineages at `age0`

is `rho0`

. Similarly, the returned Pextinct(t) (see below) is the probability that a lineage extant at age t would not survive until `age_grid[1]`

. The same convention is used for the returned variables `Pmissing`

, `shadow_diversity`

, `PER`

, `PSR`

, `SER`

and `PND`

.

A named list with the following elements:

`success` |
Logical, indicating whether the calculation was successful. If |

`ages` |
Numerical vector of size NG, listing discrete ages (time before present) on which all returned time-curves are specified. Listed ages will be in ascending order, will cover exactly the range |

`total_diversity` |
Numerical vector of size NG, listing the predicted (deterministic) total diversity (number of extant species, denoted |

`shadow_diversity` |
Numerical vector of size NG, listing the predicted (deterministic) “shadow diversity” at the ages given in |

`Pmissing` |
Numeric vector of size NG, listing the probability that a lineage, extant at a given age, will be absent from the tree either due to extinction or due to incomplete sampling. |

`Pextinct` |
Numeric vector of size NG, listing the probability that a lineage, extant at a given age, will be fully extinct at present. Note that always |

`LTT` |
Numeric vector of size NG, listing the number of lineages represented in the tree at any given age, also known as “lineages-through-time” (LTT) curve. Note that |

`lambda` |
Numeric vector of size NG, listing the speciation rate (in units 1/time) at the ages given in |

`mu` |
Numeric vector of size NG, listing the extinction rate (in units 1/time) at the ages given in |

`diversification_rate` |
Numeric vector of size NG, listing the net diversification rate ( |

`PDR` |
Numeric vector of size NG, listing the pulled diversification rate (PDR, in units 1/time) at the ages given in |

`PND` |
Numeric vector of size NG, listing the pulled normalized diversity (PND, in units 1/time) at the ages given in |

`SER` |
Numeric vector of size NG, listing the “shadow extinction rate” (SER, in units 1/time) at the ages given in |

`PER` |
Numeric vector of size NG, listing the “pulled extinction rate” (PER, in units 1/time) at the ages given in |

`PSR` |
Numeric vector of size NG, listing the “pulled speciation rate” (PSR, in units 1/time) at the ages given in |

`rholambda0` |
Non-negative numeric, specifying the product of the sampling fraction and the speciation rate at |

Stilianos Louca

H. Morlon, T. L. Parsons, J. B. Plotkin (2011). Reconciling molecular phylogenies with the fossil record. Proceedings of the National Academy of Sciences. 108:16327-16332.

S. Louca et al. (2018). Bacterial diversification through geological time. Nature Ecology & Evolution. 2:1458-1467.

# define an HBD model with exponentially decreasing speciation/extinction rates Ntips = 1000 beta = 0.01 # exponential decay rate of lambda over time oldest_age= 10 age_grid = seq(from=0,to=oldest_age,by=0.1) # choose a sufficiently fine age grid lambda = 1*exp(beta*age_grid) # define lambda on the age grid mu = 0.2*lambda # assume similarly shaped but smaller mu # simulate deterministic HBD model simulation = simulate_deterministic_hbd(LTT0 = Ntips, oldest_age = oldest_age, rho0 = 0.5, age_grid = age_grid, lambda = lambda, mu = mu) # plot deterministic LTT plot( x = simulation$ages, y = simulation$LTT, type='l', main='dLTT', xlab='age', ylab='lineages')

[Package *castor* version 1.7.0 Index]