Learn R Programming

ape (version 4.0)

dbd: Probability Density Under Birth--Death Models

Description

These functions compute the probability density under some birth--death models, that is the probability of obtaining x species after a time t giving how speciation and extinction probabilities vary through time (these may be constant, or even equal to zero for extinction).

Usage

dyule(x, lambda = 0.1, t = 1, log = FALSE) dbd(x, lambda, mu, t, conditional = FALSE, log = FALSE) dbdTime(x, birth, death, t, conditional = FALSE, BIRTH = NULL, DEATH = NULL, fast = FALSE)

Arguments

x
a numeric vector of species numbers (see Details).
lambda
a numerical value giving the probability of speciation; can be a vector with several values for dyule.
mu
id. for extinction.
t
id. for the time(s).
log
a logical value specifying whether the probabilities should be returned log-transformed; the default is FALSE.
conditional
a logical specifying whether the probabilities should be computed conditional under the assumption of no extinction after time t.
birth, death
a (vectorized) function specifying how the speciation or extinction probability changes through time (see yule.time and below).
BIRTH, DEATH
a (vectorized) function giving the primitive of birth or death.
fast
a logical value specifying whether to use faster integration (see bd.time).

Value

a numeric vector.

Details

These three functions compute the probabilities to observe x species starting from a single one after time t (assumed to be continuous). The first function is a short-cut for the second one with mu = 0 and with default values for the two other arguments. dbdTime is for time-varying lambda and mu specified as R functions.

dyule is vectorized simultaneously on its three arguments x, lambda, and t, according to R's rules of recycling arguments. dbd is vectorized simultaneously x and t (to make likelihood calculations easy), and dbdTime is vectorized only on x; the other arguments are eventually shortened with a warning if necessary.

The returned value is, logically, zero for values of x out of range, i.e., negative or zero for dyule or if conditional = TRUE. However, it is not checked if the values of x are positive non-integers and the probabilities are computed and returned.

The details on the form of the arguments birth, death, BIRTH, DEATH, and fast can be found in the links below.

References

Kendall, D. G. (1948) On the generalized ``birth-and-death'' process. Annals of Mathematical Statistics, 19, 1--15.

See Also

bd.time, yule.time

Examples

Run this code
x <- 0:10
plot(x, dyule(x), type = "h", main = "Density of the Yule process")
text(7, 0.85, expression(list(lambda == 0.1, t == 1)))

y <- dbd(x, 0.1, 0.05, 10)
z <- dbd(x, 0.1, 0.05, 10, conditional = TRUE)
d <- rbind(y, z)
colnames(d) <- x
barplot(d, beside = TRUE, ylab = "Density", xlab = "Number of species",
        legend = c("unconditional", "conditional on\nno extinction"),
        args.legend = list(bty = "n"))
title("Density of the birth-death process")
text(17, 0.4, expression(list(lambda == 0.1, mu == 0.05, t == 10)))

## Not run: 
# ### generate 1000 values from a Yule process with lambda = 0.05
# x <- replicate(1e3, Ntip(rlineage(0.05, 0)))
# 
# ### the correct way to calculate the log-likelihood...:
# sum(dyule(x, 0.05, 50, log = TRUE))
# 
# ### ... and the wrong way:
# log(prod(dyule(x, 0.05, 50)))
# 
# ### a third, less preferred, way:
# sum(log(dyule(x, 0.05, 50)))
# ## End(Not run)

Run the code above in your browser using DataLab