man/chronopl.Rd

   1 \name{chronopl}
   2 \alias{chronopl}
   3 \title{Molecular Dating With Penalized Likelihood}
   4 \usage{
   5 chronopl(phy, lambda, age.min = 1, age.max = NULL,
   6          node = "root", S = 1, tol = 1e-8,
   7          CV = FALSE, eval.max = 500, iter.max = 500, ...)
   8 }
   9 \arguments{
  10   \item{phy}{an object of class \code{"phylo"}.}
  11   \item{lambda}{value of the smoothing parameter.}
  12   \item{age.min}{numeric values specifying the fixed node ages (if
  13     \code{age.max = NULL}) or the youngest bound of the nodes known to
  14     be within an interval.}
  15   \item{age.max}{numeric values specifying the oldest bound of the nodes
  16     known to be within an interval.}
  17   \item{node}{the numbers of the nodes whose ages are given by
  18     \code{age.min}; \code{"root"} is a short-cut for the root.}
  19   \item{S}{the number of sites in the sequences; leave the default if
  20     branch lengths are in mean number of substitutions.}
  21   \item{tol}{the value below which branch lengths are considered
  22     effectively zero.}
  23   \item{CV}{whether to perform cross-validation.}
  24   \item{eval.max}{the maximal number of evaluations of the penalized
  25     likelihood function.}
  26   \item{iter.max}{the maximal number of iterations of the optimization
  27     algorithm.}
  28   \item{\dots}{further arguments passed to control \code{nlminb}.}
  29 }
  30 \description{
  31   This function estimates the node ages of a tree using a
  32   semi-parametric method based on penalized likelihood (Sanderson
  33   2002). The branch lengths of the input tree are interpreted as mean
  34   numbers of substitutions (i.e., per site).
  35 }
  36 \details{
  37   The idea of this method is to use a trade-off between a parametric
  38   formulation where each branch has its own rate, and a nonparametric
  39   term where changes in rates are minimized between contiguous
  40   branches. A smoothing parameter (lambda) controls this trade-off. If
  41   lambda = 0, then the parametric component dominates and rates vary as
  42   much as possible among branches, whereas for increasing values of
  43   lambda, the variation are smoother to tend to a clock-like model (same
  44   rate for all branches).
  45
  46   \code{lambda} must be given. The known ages are given in
  47   \code{age.min}, and the correponding node numbers in \code{node}.
  48   These two arguments must obviously be of the same length. By default,
  49   an age of 1 is assumed for the root, and the ages of the other nodes
  50   are estimated.
  51
  52   If \code{age.max = NULL} (the default), it is assumed that
  53   \code{age.min} gives exactly known ages. Otherwise, \code{age.max} and
  54   \code{age.min} must be of the same length and give the intervals for
  55   each node. Some node may be known exactly while the others are
  56   known within some bounds: the values will be identical in both
  57   arguments for the former (e.g., \code{age.min = c(10, 5), age.max =
  58     c(10, 6), node = c(15, 18)} means that the age of node 15 is 10
  59   units of time, and the age of node 18 is between 5 and 6).
  60
  61   The input tree may have multichotomies. If some internal branches are
  62   of zero-length, they are collapsed (with a warning), and the returned
  63   tree will have less nodes than the input one. The presence of
  64   zero-lengthed terminal branches of results in an error since it makes
  65   little sense to have zero-rate branches.
  66
  67   The cross-validation used here is different from the one proposed by
  68   Sanderson (2002). Here, each tip is dropped successively and the
  69   analysis is repeated with the reduced tree: the estimated dates for
  70   the remaining nodes are compared with the estimates from the full
  71   data. For the \eqn{i}{i}th tip the following is calculated:
  72
  73   \deqn{\sum_{j=1}^{n-2}{\frac{(t_j - t_j^{-i})^2}{t_j}}}{SUM[j = 1, ..., n-2] (tj - tj[-i])^2/tj},
  74
  75   where \eqn{t_j}{tj} is the estimated date for the \eqn{j}{j}th node
  76   with the full phylogeny, \eqn{t_j^{-i}}{tj[-i]} is the estimated date
  77   for the \eqn{j}{j}th node after removing tip \eqn{i}{i} from the tree,
  78   and \eqn{n}{n} is the number of tips.
  79
  80   The present version uses the \code{\link[stats]{nlminb}} to optimise
  81   the penalized likelihood function: see its help page for details on
  82   parameters controlling the optimisation procedure.
  83 }
  84 \value{
  85   an object of class \code{"phylo"} with branch lengths as estimated by
  86   the function. There are three or four further attributes:
  87
  88   \item{ploglik}{the maximum penalized log-likelihood.}
  89   \item{rates}{the estimated rates for each branch.}
  90   \item{message}{the message returned by \code{nlminb} indicating
  91     whether the optimisation converged.}
  92   \item{D2}{the influence of each observation on overall date
  93     estimates (if \code{CV = TRUE}).}
  94 }
  95 \references{
  96   Sanderson, M. J. (2002) Estimating absolute rates of molecular
  97   evolution and divergence times: a penalized likelihood
  98   approach. \emph{Molecular Biology and Evolution}, \bold{19},
  99   101--109.
 100 }
 101 \author{Emmanuel Paradis}
 102 \seealso{
 103   \code{\link{chronogram}}, \code{\link{ratogram}},
 104   \code{\link{NPRS.criterion}}, \code{\link{chronoMPL}}
 105 }
 106 \keyword{models}