https://github.com/cran/ape
Tip revision: f68e51c050f9ec5cc674d158fc58a0835198d441 authored by Emmanuel Paradis on 19 December 2005, 00:00:00 UTC
version 1.8
version 1.8
Tip revision: f68e51c
skyline.Rd
\name{skyline}
\alias{skyline}
\alias{skyline.phylo}
\alias{skyline.coalescentIntervals}
\alias{skyline.collapsedIntervals}
\alias{find.skyline.epsilon}
\title{Skyline Plot Estimate of Effective Population Size}
\usage{
skyline(x, \dots)
skyline.phylo(x, \dots)
skyline.coalescentIntervals(x, epsilon=0, \dots)
skyline.collapsedIntervals(x, old.style=FALSE, \dots)
find.skyline.epsilon(ci, GRID=1000, MINEPS=1e-6, \dots)
}
\arguments{
\item{x}{Either an ultrametric tree (i.e. an object of class \code{"phylo"}),
or coalescent intervals (i.e. an object of class \code{"coalescentIntervals"}), or
collapsed coalescent intervals (i.e. an object of class \code{"collapsedIntervals"}).}
\item{epsilon}{collapsing parameter that controls the amount of smoothing
(allowed range: from \code{0} to \code{ci$total.depth}, default value: 0). This is the same parameter as in
\link{collapsed.intervals}.}
\item{old.style}{Parameter to choose between two slightly different variants of the
generalized skyline plot (Strimmer and Pybus, pers. comm.). The default value \code{FALSE} is
recommended.}
\item{ci}{coalescent intervals (i.e. an object of class \code{"coalescentIntervals"})}
\item{GRID}{Parameter for the grid search for \code{epsilon} in \code{find.skyline.epsilon}.}
\item{MINEPS}{Parameter for the grid search for \code{epsilon} in \code{find.skyline.epsilon}.}
\item{...}{Any of the above parameters.}
}
\description{
\code{skyline} computes the \emph{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 \code{epsilon}.
\code{find.skyline.epsilon} searches for an optimal value of the \code{epsilon} parameter,
i.e. the value that maximizes the AICc-corrected log-likelihood (\code{logL.AICc}).
}
\details{
\code{skyline} implements the \emph{generalized skyline plot} introduced in
Strimmer and Pybus (2001). For \code{epsilon = 0} the
generalized skyline plot degenerates to the
\emph{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).
}
\value{
\code{skyline} returns an object of class \code{"skyline"} with the following entries:
\item{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)}
\item{interval.length}{ A vector with the length of each
interval.}
\item{population.size}{A vector with the effective population size of each interval.}
\item{parameter.count}{ Number of free parameters in the skyline plot.}
\item{epsilon}{The value of the underlying smoothing parameter.}
\item{logL}{Log-likelihood of skyline plot (see Strimmer and Pybus, 2001).}
\item{logL.AICc}{AICc corrected log-likelihood (see Strimmer and Pybus, 2001).}
\code{find.skyline.epsilon} returns the value of the \code{epsilon} parameter
that maximizes \code{logL.AICc}.
}
\author{Korbinian Strimmer (\url{http://www.stat.uni-muenchen.de/~strimmer/})}
\seealso{
\code{\link{coalescent.intervals}}, \code{\link{collapsed.intervals}},
\code{\link{skylineplot}}, \code{\link{ltt.plot}}.
}
\references{
Strimmer, K. and Pybus, O. G. (2001) Exploring the demographic history
of DNA sequences using the generalized skyline plot. \emph{Molecular
Biology and Evolution}, \bold{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. \emph{Genetics}, \bold{155}, 1429--1437.
Nee, S., Holmes, E. C., Rambaut, A. and Harvey, P. H. (1995) Inferring
population history from molecular phylogenies. \emph{Philosophical
Transactions of the Royal Society of London. Series B. Biological
Sciences}, \bold{349}, 25--31.
}
\examples{
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
}
\keyword{manip}