simulate_deterministic_hbds {castor} | R Documentation |

Given a homogenous birth-death-sampling (HBDS) model, i.e., with speciation rate *λ*, extinction rate *μ*, continuous (Poissonian) sampling rate *ψ* and retention probability *κ*, calculate various deterministic features of the model backwards in time, such as the total population size and the LTT over time.
Continuously sampled lineages are kept in the pool of extant lineages at probability *κ*. The variables *λ*, *μ*, *ψ* and *κ* may depend on time.
In addition, the model can include concentrated sampling attempts at a finite set of discrete time points *t_1,..,t_m*. “Homogenous” refers to the assumption that, at any given moment in time, all lineages exhibit the same speciation/extinction/sampling rates and retention proabbility. Every HBDS model is thus defined based on the values that *λ*, *μ*, *ψ* and *κ* take over time, as well as the sampling probabilities *ψ_1,..,ψ_m* and retention probabilities *κ_1,..,κ_m* during the concentrated sampling attempts. Special cases of this model are sometimes known as “birth-death-skyline plots” in the literature (Stadler 2013). In epidemiology, these models are often used to describe the phylogenies of viral strains sampled over the course of the epidemic. A “concentrated sampling attempt” is a brief but intensified sampling period that lasted much less than the typical timescales of speciation/extinction. “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”).
The time-profiles of *λ*, *μ*, *ψ* and *κ* are specified as piecewise polynomial curves (splines), defined on a discrete grid of ages.

simulate_deterministic_hbds(age_grid = NULL, lambda = NULL, mu = NULL, psi = NULL, kappa = NULL, splines_degree = 1, CSA_ages = NULL, CSA_probs = NULL, CSA_kappas = NULL, requested_ages = NULL, age0 = 0, N0 = NULL, LTT0 = NULL, ODE_relative_dt = 0.001, ODE_relative_dy = 1e-4)

`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 |

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

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

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

`CSA_ages` |
Optional numeric vector, listing the ages of concentrated sampling attempts, in ascending order. Concentrated sampling is performed in addition to any continuous (Poissonian) sampling specified by |

`CSA_probs` |
Optional numeric vector of the same size as |

`CSA_kappas` |
Optional numeric vector of the same size as |

`requested_ages` |
Optional numeric vector, listing ages (in ascending order) at which the various model variables are requested. If |

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

`N0` |
Positive numeric, specifying the number of extant species (sampled or not) at |

`LTT0` |
Positive numeric, specifying the number of lineages present in the tree at |

`ODE_relative_dt` |
Positive unitless number, specifying the default relative time step for internally used ordinary differential equation solvers. Typical values are 0.01-0.001. |

`ODE_relative_dy` |
Positive unitless number, specifying the relative difference between subsequent simulated and interpolated values, in internally used ODE solvers. Typical values are 1e-2 to 1e-5. A smaller |

The simulated LTT refers to a hypothetical tree sampled at age 0, i.e. LTT(t) will be the number of lineages extant at age t that survived and were sampled until by the present day. Note that if a concentrated sampling attempt occurs at age *τ*, then LTT(*τ*) is the number of lineages in the tree right before the occurrence of the sampling attempt, i.e., in the limit where *τ* is approached from above.

Note that in this function age always refers to time before present, i.e., present day age is 0, and age increases towards the root.

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. Will be equal to |

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

`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 |

`nLTT` |
Numeric vector of size NG, listing values of the normalized LTT at ages |

`Pmissing` |
Numeric vector of size NG, listing the probability that a lineage, extant at a given age, will not be represented in the tree. |

`lambda` |
Numeric vector of size NG, listing the speciation rates at the ages given in |

`mu` |
Numeric vector of size NG, listing the extinctions rates at the ages given in |

`psi` |
Numeric vector of size NG, listing the Poissonian sampling rates at the ages given in |

`kappa` |
Numeric vector of size NG, listing the retention probabilities (for continuously sampled lineages) at the ages given in |

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

`IPRP` |
Numeric vector of size NG, listing the age-integrated pulled diversification rate 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 |

`PRP` |
Numeric vector of size NG, listing the “pulled retention probability” (PRP) at the ages given in |

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

`branching_density` |
Numeric vector of size NG, listing the deterministic branching density (PSR * nLTT, in units nodes/time) at the ages given in |

`sampling_density` |
Numeric vector of size NG, listing the deterministic sampling density ( |

`lambda_psi` |
Numeric vector of size NG, listing the product of the speciation rate and Poissonian sampling rate (in units 1/time^2) at the ages given in |

`kappa_psi` |
Numeric vector of size NG, listing the product of the continuous sampling rate and the continuous retention probability (in units 1/time) at the ages given in |

`Reff` |
Numeric vector of size NG, listing the effective reproduction ratio ( |

`removal_rate` |
Numeric vector of size NG, listing the total removal rate ( |

`sampling_proportion` |
Numeric vector of size NG, listing the instantaneous sampling proportion ( |

`CSA_pulled_probs` |
Numeric vector of size NG, listing the pulled concentrated sampling probabilities, |

`CSA_psis` |
Numeric vector of size NG, listing the continuous (Poissonian) sampling rates during the concentrated sampling attempts, |

`CSA_PSRs` |
Numeric vector of size NG, listing the pulled speciation rates during the concentrated sampling attempts. |

Stilianos Louca

T. Stadler, D. Kuehnert, S. Bonhoeffer, A. J. Drummond (2013). Birth-death skyline plot reveals temporal changes of epidemic spread in HIV and hepatitis C virus (HCV). PNAS. 110:228-233.

`generate_tree_hbds`

,
`fit_hbds_model_parametric`

,
`simulate_deterministic_hbd`

# define an HBDS model with exponentially decreasing speciation/extinction rates # and constant Poissonian sampling rate psi oldest_age= 10 age_grid = seq(from=0,to=oldest_age,by=0.1) # choose a sufficiently fine age grid lambda = 1*exp(0.01*age_grid) # define lambda on the age grid mu = 0.2*lambda # assume similarly shaped but smaller mu # simulate deterministic HBD model # scale LTT such that it is 100 at age 1 simulation = simulate_deterministic_hbds(age_grid = age_grid, lambda = lambda, mu = mu, psi = 0.1, age0 = 1, LTT0 = 100) # plot deterministic LTT plot( x = simulation$ages, y = simulation$LTT, type='l', main='dLTT', xlab='age', ylab='lineages', xlim=c(oldest_age,0))

[Package *castor* version 1.7.0 Index]