apTreeshape (version 1.5-0.1)

get_tree_beta: Beta parameter as a function of the proportion of remaining tips

Description

Computes the maximum likelihood estimate of the parameter beta in trees simulated by the model, as a function of the proportion of conserved species

Usage

get_tree_beta(epsilon, beta, alpha, N, sampl.frac, ntree, equal.ab = TRUE, eta = 1)

Arguments

epsilon

Minimum size of unsampled splits (see appendix 1)

beta

Imbalance index

alpha

Clade age-richness index

N

Initial tip number

sampl.frac

Vector of tips fractions for which we want to compute the beta statistic

ntree

Number of simulated trees

equal.ab

If set to TRUE, all species have the same probability to go extinct first (default to TRUE)

eta

Clade abundance-richness index (if equal.ab == FALSE)

Value

A table of size length(sample.fract)*ntree. The element in ligne i and column j is the Maximum Likelihood Estimate for the beta statistic of the j^th tree in wich a fraction sample.frac[i] has been sampled.

References

Maliet O., Gascuel F., Lambert A. (2018) Ranked tree shapes, non-random extinctions and the loss of phylogenetic diversity, bioRxiv 224295, doi: https://doi.org/10.1101/224295

Examples

# NOT RUN {
# With the field of bullets hypothesis

set.seed(813)
sampl.frac=seq(0.2,1,0.2)
Beta=get_tree_beta(epsilon=0.01, beta=-1, alpha=-1, N=20, sampl.frac=sampl.frac, ntree=3)

Beta_quantiles=sapply(1:nrow(Beta),function(x){quantile(Beta[x,],c(0.05,0.5,0.95))})

plot(1, type="n", xlab="Fraction of extinct species, p", ylab="Beta statistic", 
      ylim=c(-2,10), xlim=c(0,1))
polygon(c(1-sampl.frac, rev(1-sampl.frac)), c(Beta_quantiles[1,(1:length(sampl.frac))],
        rev(Beta_quantiles[3,(1:length(sampl.frac))])), border=NA, col=grey(0.7))
points(1-sampl.frac, Beta_quantiles[2,(1:length(sampl.frac))],t="l")


# With nonrandom extinctions

# }
# NOT RUN {
set.seed(813)
sampl.frac=seq(0.2,1,0.2)
Beta=get_tree_beta(epsilon=0.01, beta=5, alpha=2, eta=2, N=20, 
                    sampl.frac=sampl.frac, ntree=3)

Beta_quantiles=sapply(1:nrow(Beta),function(x){quantile(Beta[x,],c(0.05,0.5,0.95))})

plot(1, type="n", xlab="Fraction of extinct species, p", 
      ylab="Beta statistic", ylim=c(-2,10), xlim=c(0,1))
polygon(c(1-sampl.frac, rev(1-sampl.frac)), 
      c(Beta_quantiles[1,(1:length(sampl.frac))], 
      rev(Beta_quantiles[3,(1:length(sampl.frac))])), border=NA, col=grey(0.7))
points(1-sampl.frac, Beta_quantiles[2,(1:length(sampl.frac))],t="l")
# }
# NOT RUN {

# }