Last data update: 2014.03.03

R: yule-n-rate
yule-n-rateR Documentation



Fits multi-rate variants of the pure birth (Yule) model to branching times derived from phylogenetic data. For example, the yule2rate model assumes that the clade has diversified under speciation rate r1 until some time st, at which point the speciation rate shifts to a new rate r2. The shift point(s) are found by optimizing parameters and computing likelihoods for a set of possible shift times and selecting the parameter combinations giving the maximum log-likelihood.


yule2rate(x, verbose = FALSE, ints = NULL, file = "out_yule2rate.txt")

yule3rate(x, ints = NULL, verbose = FALSE, file = "out_yule3rate.txt") 

yule4rate(x, ints = NULL)

yule5rate(x, ints = NULL)



a numeric vector of branching times


if 'verbose = TRUE', writes likelihoods and parameter estimates for all shift points considered to the specified file. Default is FALSE. Only available for yule2rate and yule3rate models


the number of intervals. See details


a filename for output if 'verbose = TRUE'


'verbose' - for yule2rate and yule3rate models, maximum log-likelihoods and parameter estimates for each shift time under consideration will be output to file. The file can then be loaded to examine the likelihood of a rate shift at different points in time.

'ints' is used in determining the number of shift points to consider. If 'ints = NULL' (the default), the model will consider only observed branching times as possible shift points. Suppose we have a small dataset with the following branching times: (100, 80, 50, 40, 30, 20, 10, 5, 2). Under the yule2rate model, we assume that the clade has diversified under some rate r1 until some time ts, at which point the rate simultaneously shifts in all lineages to a new rate r2. In this example, if 'ints = NULL', we would use the set of observed branching times only, but omitting the first and final branching times (thus, we would be considering only st = (80, 50, 40, 30, 20, 10, 5) as possible shift points. If 'ints = 100', we would consider 100 evenly spaced shift points on the interval between the 2nd branching time and the N-1 branching time (e.g., on (80, 5)). 'ints' works well for yule2rate and yule3rate models, but can result in high computational times for yule4rate and yule5rate models.


a data frame containing the following elements:


the maximum log-likelihood


the Akaike Information Criterion


the first (earliest) speciation rate giving the maximum log-likelihood


the ML estimate of the second speciation rate (in the case of yule2rate, this will be the final speciation rate)


the earliest shift point (for yule2rate, you will only have a single shift point)


speciation rates and shift points for models other than yule2rate are abbreviated by r3, st2, etc, as described above


The total number of parameters for each model is equal to the number of speciation rates and shift points subject to optimization. Thus, the yule3rate model has 3 speciation rates and 2 shift times, for a total of 5 parameters. Strictly speaking, it may be inappropriate to treat the shift time st as a free parameter, as it can only take on a limited set of values. However, in practice, it appears to work well; in many cases, using observed shift times can give higher likelihoods than when 'ints' are specified. There seems to be little improvement in the log-likelihood with 'ints' greater than 1000, at least for phylogenies with fewer than 100 tips.

Note that shift times, like branching times, are given in divergence units before present. Thus, if you have scaled a set of branching times to a basal divergence of 30 million years before present, you would interpret 'st1 = 19.5' as an inferred shift point 19.5 million years before present.


Dan Rabosky Klaus Schlep


Barraclough, T. G., and A. P. Vogler. 2002. Recent diversification rates in North American tiger beetles estimated from a dated mtDNA phylogenetic tree. Mol. Biol. Evol. 19:1706-1716.

Rabosky, D. L. 2006. Likelihood methods for inferring temporal shifts in diversification rates. Evolution 60:1152-1164.

See Also

bd , fitdAICrc , yuleWindow , pureBirth



### fitting a 3-rate Yule model to the Plethodon data:

result <- yule3rate(plethodon)

### gives data frame with maximum log-likelihood and parameter estimates 
### at the max.
### In this case, we would access individual parameters as
### result$LH (the max), result$st1 (first shift time), result$st2
### (the second shift time), and result$r1, result$r2, and result$r3
### for the speciation rates.

### Here we will use 'yule2rate' to output maximum log-likelihoods 
### for each shift point considered, then load the file and plot 
### log-likelihoods of a rate shift against 'time from basal divergence'
### to graphically explore the tempo of diversification

# result <- yule2rate(plethodon, ints = NULL, 
#          verbose = TRUE, file = 'out.txt')
# LHtable <- read.table(file = 'out.txt', header = TRUE)

### 'header = TRUE' ensures that variable names are correctly read

### rescaling shift times so that they reflect 'time from basal divergence':

# LHtable$st1 <- plethodon[1] - LHtable$st1
# plot(LHtable$LH~LHtable$st1, xlab = 'Time From Basal Divergence', 
#    ylab = 'Log-likelihood')