skyline                 package:ape                 R Documentation

_S_k_y_l_i_n_e _P_l_o_t _E_s_t_i_m_a_t_e _o_f _E_f_f_e_c_t_i_v_e _P_o_p_u_l_a_t_i_o_n _S_i_z_e

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

     'skyline' computes the _generalized skyline plot_ estimate of
     effective population size from an estimated phylogeny.  The
     demographic history is approximated by  a step-function.  The
     number of parameters of the skyline plot (i.e. its smoothness) is
     controlled by a parameter 'epsilon'. 

     'find.skyline.epsilon' searches for an optimal value of the
     'epsilon' parameter, i.e. the value that maximizes the
     AICc-corrected log-likelihood ('logL.AICc').

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

     skyline(x, ...)
     skyline.phylo(x, ...)
     skyline.coalescentIntervals(x, epsilon=0, ...)
     skyline.collapsedIntervals(x, old.style=FALSE, ...)
     find.skyline.epsilon(ci, GRID=1000, MINEPS=1e-6, ...)

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

       x: Either an ultrametric tree (i.e. an object of class
          '"phylo"'),  or coalescent intervals (i.e. an object of class
          '"coalescentIntervals"'), or collapsed coalescent intervals
          (i.e. an object of class '"collapsedIntervals"').

 epsilon: collapsing parameter that controls the amount of smoothing
          (allowed range: from '0' to 'ci$total.depth', default value:
          0). This is the same parameter as in collapsed.intervals.

old.style: Parameter to choose between two slightly different variants
          of the generalized skyline plot (Strimmer and Pybus, pers.
          comm.). The default value 'FALSE' is recommended.

      ci: coalescent intervals (i.e. an object of class
          '"coalescentIntervals"')

    GRID: Parameter for the grid search for 'epsilon' in
          'find.skyline.epsilon'.

  MINEPS: Parameter for the grid search for 'epsilon' in
          'find.skyline.epsilon'.

     ...: Any of the above parameters.

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

     'skyline' implements the _generalized skyline plot_  introduced in
      Strimmer and Pybus (2001).  For 'epsilon = 0' the generalized
     skyline plot degenerates to the  _classic skyline plot_ described
     in Pybus et al. (2000).  The latter is in turn directly related to
     lineage-through-time plots (Nee et al., 1995).

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

     'skyline' returns an object of class '"skyline"' with the
     following entries:

    time: A vector with the time at the end of each coalescent interval
          (i.e. the accumulated interval lengths from the beginning of
          the first interval to the end of an interval)

interval.length: A vector with the length of each  interval.

population.size: A vector with the effective population size of each
          interval.

parameter.count: Number of free parameters in the skyline plot.

 epsilon: The value of the underlying smoothing parameter.

    logL: Log-likelihood of skyline plot (see Strimmer and Pybus,
          2001).

logL.AICc: AICc corrected log-likelihood (see Strimmer and Pybus,
          2001).


     'find.skyline.epsilon' returns the value of the 'epsilon'
     parameter that maximizes 'logL.AICc'.

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

     Korbinian Strimmer (<URL:
     http://www.stat.uni-muenchen.de/~strimmer/>)

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

     Strimmer, K. and Pybus, O. G. (2001) Exploring the demographic
     history of DNA sequences using the generalized skyline plot.
     _Molecular Biology and Evolution_, *18*, 2298-2305.

     Pybus, O. G, Rambaut, A. and Harvey, P. H. (2000) An integrated
     framework for the inference of viral population history from
     reconstructed genealogies. _Genetics_, *155*, 1429-1437.

     Nee, S., Holmes, E. C., Rambaut, A. and Harvey, P. H. (1995)
     Inferring population history from molecular phylogenies.
     _Philosophical Transactions of the Royal Society of London. Series
     B. Biological Sciences_, *349*, 25-31.

_S_e_e _A_l_s_o:

     'coalescent.intervals', 'collapsed.intervals', 'skylineplot',
     'ltt.plot'.

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

     library(ape)

     # get tree
     data("hivtree.newick") # example tree in NH format
     tree.hiv <- read.tree(text = hivtree.newick) # load tree

     # corresponding coalescent intervals
     ci <- coalescent.intervals(tree.hiv) # from tree

     # collapsed intervals
     cl1 <- collapsed.intervals(ci,0)
     cl2 <- collapsed.intervals(ci,0.0119)

     #### classic skyline plot ####
     sk1 <- skyline(cl1)        # from collapsed intervals 
     sk1 <- skyline(ci)         # from coalescent intervals
     sk1 <- skyline(tree.hiv)   # from tree
     sk1

     plot(skyline(tree.hiv))
     skylineplot(tree.hiv) # shortcut

     plot(sk1, show.years=TRUE, subst.rate=0.0023, present.year = 1997)

     #### generalized skyline plot ####

     sk2 <- skyline(cl2)              # from collapsed intervals
     sk2 <- skyline(ci, 0.0119)       # from coalescent intervals
     sk2 <- skyline(tree.hiv, 0.0119) # from tree
     sk2

     plot(sk2)

     # classic and generalized skyline plot together in one plot
     plot(sk1, show.years=TRUE, subst.rate=0.0023, present.year = 1997, col=c(grey(.8),1))
     lines(sk2,  show.years=TRUE, subst.rate=0.0023, present.year = 1997)
     legend(.15,500, c("classic", "generalized"), col=c(grey(.8),1),lty=1)

     # find optimal epsilon parameter using AICc criterion
     find.skyline.epsilon(ci)

     sk3 <- skyline(ci, -1) # negative epsilon also triggers estimation of epsilon
     sk3$epsilon

