rtreeshape            package:apTreeshape            R Documentation

_G_e_n_e_r_a_t_e _a _l_i_s_t _o_f _r_a_n_d_o_m _b_i_n_a_r_y _t_r_e_e_s _a_c_c_o_r_d_i_n_g _t_o _a _g_i_v_e_n _m_o_d_e_l

_D_e_s_c_r_i_p_t_i_o_n:

     This function generates a tree or a list of trees of class
     '"treeshape"' according to the Yule, PDA, Biased models.
     Speciation-specified models are also allowed.

_U_s_a_g_e:

     rtreeshape(n, tip.number, p = 0.3, model="", FUN="")

_A_r_g_u_m_e_n_t_s:

       n: The number of trees to generate.

tip.number: The number of tips of the trees to generate. Can be a
          vector.

       p: Only used when 'model="biased"'. It represents the bias
          factor of the tree to generate.

   model: A character string equals to '"yule"', '"pda"',
          code{"aldous"} or '"biased"'.

     FUN: A two variables ('n' and 'i') function.

_D_e_t_a_i_l_s:

     The '"FUN"' and '"model"' arguments cannot be specified at the
     same time. An error will be returned if both arguments are
     specified. 

     If 'tip.number' is a vector, 'n' trees will be generated for each
     size contained by 'tip.number'

     'Q' enables you to build trees of class '"treeshape"' according to
     the _Markov branching_ model described by D. Aldous. 'Qn(i)' is
     the probability that the left daughter clade of an internal node
     with 'n' descendents contains 'i' tips. The 'Qn(i)' need not sum
     to one. Still, be carefull when you specify this distribution:
     computational errors may occur for complicated distributions
     and/or large trees. 

     The Yule model, also known as Markov model, can be described as
     follows. At any time, all the extant branches have the same
     probability to split into two subspecies. 
      The PDA model (Proportional to Distinguishable Arrangements) is
     not a model of growing tree. Instead, each tree with n tips has
     the same probability to be generated under this model. There is
     (2n-3)!! possible trees with n tips.
      The Biased model is a model of growing tree. When a species with
     speciation rate _r_ splits, one of its descendent species is given
     the rate _pr_ and the other is given the speciation rate (1-pr)
     where p is a probability parameter. The Biased model was
     introduced by Kirkpatrick and Slatkin (1993). The Aldous'
     Branching (AB) model is defined by the following symmetric split
     distribution q(n,i) = (1/(2*h(n-1))) * (1/(i(n-i))), where h(n) is
     the nth harmonic number. The AB model is hardly motivated by
     biological considerations.

_V_a_l_u_e:

     A list of objects of class '"treeshape"'  NULL if 'n'=0

_A_u_t_h_o_r(_s):

     Michael Blum <michael.blum@imag.fr>
       Nicolas Bortolussi <nicolas.bortolussi@imag.fr>
      Eric Durand <eric.durand@imag.fr>
      Olivier Franois <olivier.francois@imag.fr>

_R_e_f_e_r_e_n_c_e_s:

     Mooers, A. O. and Heard, S. B. (Mar., 1997), Inferring
     Evolutionnary Process from Phylogenetic Tree Shape. _The Quarterly
     Review of Biology_, *72*, 31-54, for more details about the Yule
     and PDA models.

     Aldous, D. J. (1996), _Probability Distributions on Cladograms_.
     pp.1-18 of _Random Discrete Structures_ eds D. Aldous and R.
     Pemantle, IMA Volumes Math. Appl. 76.

     Kirkpatrick, M. and Slatkin, M. (1993) Searching for evolutionary
     patterns in the shape of a phylogenetic tree. _Evolution_, *47*,
     1171 - 1181.

_E_x_a_m_p_l_e_s:

     ## Summary of a PDA tree with 100 tips:
     summary(rtreeshape(n=1, tip.number=100, model="pda")[[1]])
     ## Summary of a Yule tree with 100 tips:
     summary(rtreeshape(n=1, tip.number=100, model="yule")[[1]])
       
     ## Generate trees with different sizes
     trees=rtreeshape(n=2, tip.number=c(10,20), model="yule")
     length(trees)
     plot(trees[[1]])
     plot(trees[[2]])
       
     ## Histogram of Colless' indices for a list of 1000 PDA trees with 60 tips
     hist(sapply(rtreeshape(1000,60,model="pda"),FUN=colless,norm="pda"),freq=FALSE)

     ## Histogram of shape statistics for a list of 1000 Yule trees with 100 tips 
     ##      (takes some time to compute) 
     main="Histogram of shape statistics for a list of 1000 Yule trees"
     hist(sapply(rtreeshape(1000,100,model="yule"),FUN=shape.statistic,norm="yule"),
           freq=FALSE, main=main)
     ## It should be a gaussian with mean 0 and standard deviation 1.
     x<-seq(-4,4,by=0.01)
     lines(x,dnorm(x))       

     ## Building a tree using Markov splitting model
     Q <- function(n,i) (i==1)

     tree=rtreeshape(n=1, tip.number=10, FUN=Q)
     plot(tree[[1]])

