4 \alias{skyline.coalescentIntervals}
5 \alias{skyline.collapsedIntervals}
6 \alias{find.skyline.epsilon}
8 \title{Skyline Plot Estimate of Effective Population Size}
11 \method{skyline}{phylo}(x, \dots)
12 \method{skyline}{coalescentIntervals}(x, epsilon=0, \dots)
13 \method{skyline}{collapsedIntervals}(x, old.style=FALSE, \dots)
14 find.skyline.epsilon(ci, GRID=1000, MINEPS=1e-6, \dots)
17 \item{x}{Either an ultrametric tree (i.e. an object of class
18 \code{"phylo"}), or coalescent intervals (i.e. an object of class
19 \code{"coalescentIntervals"}), or collapsed coalescent intervals
20 (i.e. an object of class \code{"collapsedIntervals"}).}
21 \item{epsilon}{collapsing parameter that controls the amount of smoothing
22 (allowed range: from \code{0} to \code{ci$total.depth}, default value: 0). This is the same parameter as in
23 \link{collapsed.intervals}.}
25 \item{old.style}{Parameter to choose between two slightly different variants of the
26 generalized skyline plot (Strimmer and Pybus, pers. comm.). The default value \code{FALSE} is
29 \item{ci}{coalescent intervals (i.e. an object of class \code{"coalescentIntervals"})}
31 \item{GRID}{Parameter for the grid search for \code{epsilon} in \code{find.skyline.epsilon}.}
33 \item{MINEPS}{Parameter for the grid search for \code{epsilon} in \code{find.skyline.epsilon}.}
35 \item{...}{Any of the above parameters.}
40 \code{skyline} computes the \emph{generalized skyline plot} estimate of effective population size
41 from an estimated phylogeny. The demographic history is approximated by
42 a step-function. The number of parameters of the skyline plot (i.e. its smoothness)
43 is controlled by a parameter \code{epsilon}.
45 \code{find.skyline.epsilon} searches for an optimal value of the \code{epsilon} parameter,
46 i.e. the value that maximizes the AICc-corrected log-likelihood (\code{logL.AICc}).
50 \code{skyline} implements the \emph{generalized skyline plot} introduced in
51 Strimmer and Pybus (2001). For \code{epsilon = 0} the
52 generalized skyline plot degenerates to the
53 \emph{classic skyline plot} described in
54 Pybus et al. (2000). The latter is in turn directly related to lineage-through-time plots
59 \code{skyline} returns an object of class \code{"skyline"} with the following entries:
61 \item{time}{ A vector with the time at the end of each coalescent
62 interval (i.e. the accumulated interval lengths from the beginning of the first interval
63 to the end of an interval)}
65 \item{interval.length}{ A vector with the length of each
68 \item{population.size}{A vector with the effective population size of each interval.}
70 \item{parameter.count}{ Number of free parameters in the skyline plot.}
71 \item{epsilon}{The value of the underlying smoothing parameter.}
73 \item{logL}{Log-likelihood of skyline plot (see Strimmer and Pybus, 2001).}
75 \item{logL.AICc}{AICc corrected log-likelihood (see Strimmer and Pybus, 2001).}
77 \code{find.skyline.epsilon} returns the value of the \code{epsilon} parameter
78 that maximizes \code{logL.AICc}.
81 \author{Korbinian Strimmer (\url{http://www.stat.uni-muenchen.de/~strimmer/})}
84 \code{\link{coalescent.intervals}}, \code{\link{collapsed.intervals}},
85 \code{\link{skylineplot}}, \code{\link{ltt.plot}}.
90 Strimmer, K. and Pybus, O. G. (2001) Exploring the demographic history
91 of DNA sequences using the generalized skyline plot. \emph{Molecular
92 Biology and Evolution}, \bold{18}, 2298--2305.
94 Pybus, O. G, Rambaut, A. and Harvey, P. H. (2000) An integrated
95 framework for the inference of viral population history from
96 reconstructed genealogies. \emph{Genetics}, \bold{155}, 1429--1437.
98 Nee, S., Holmes, E. C., Rambaut, A. and Harvey, P. H. (1995) Inferring
99 population history from molecular phylogenies. \emph{Philosophical
100 Transactions of the Royal Society of London. Series B. Biological
101 Sciences}, \bold{349}, 25--31.
106 data("hivtree.newick") # example tree in NH format
107 tree.hiv <- read.tree(text = hivtree.newick) # load tree
109 # corresponding coalescent intervals
110 ci <- coalescent.intervals(tree.hiv) # from tree
112 # collapsed intervals
113 cl1 <- collapsed.intervals(ci,0)
114 cl2 <- collapsed.intervals(ci,0.0119)
117 #### classic skyline plot ####
118 sk1 <- skyline(cl1) # from collapsed intervals
119 sk1 <- skyline(ci) # from coalescent intervals
120 sk1 <- skyline(tree.hiv) # from tree
123 plot(skyline(tree.hiv))
124 skylineplot(tree.hiv) # shortcut
126 plot(sk1, show.years=TRUE, subst.rate=0.0023, present.year = 1997)
128 #### generalized skyline plot ####
130 sk2 <- skyline(cl2) # from collapsed intervals
131 sk2 <- skyline(ci, 0.0119) # from coalescent intervals
132 sk2 <- skyline(tree.hiv, 0.0119) # from tree
138 # classic and generalized skyline plot together in one plot
139 plot(sk1, show.years=TRUE, subst.rate=0.0023, present.year = 1997, col=c(grey(.8),1))
140 lines(sk2, show.years=TRUE, subst.rate=0.0023, present.year = 1997)
141 legend(.15,500, c("classic", "generalized"), col=c(grey(.8),1),lty=1)
144 # find optimal epsilon parameter using AICc criterion
145 find.skyline.epsilon(ci)
147 sk3 <- skyline(ci, -1) # negative epsilon also triggers estimation of epsilon