man/ace.Rd

   1 \name{ace}
   2 \alias{ace}
   3 \alias{logLik.ace}
   4 \alias{deviance.ace}
   5 \alias{AIC.ace}
   6 \alias{anova.ace}
   7 \title{Ancestral Character Estimation}
   8 \usage{
   9 ace(x, phy, type = "continuous", method = "ML", CI = TRUE,
  10     model = if (type == "continuous") "BM" else "ER",
  11     scaled = TRUE, kappa = 1, corStruct = NULL, ip = 0.1)
  12 \method{logLik}{ace}(object, ...)
  13 \method{deviance}{ace}(object, ...)
  14 \method{AIC}{ace}(object, ..., k = 2)
  15 \method{anova}{ace}(object, ...)
  16 }
  17 \arguments{
  18   \item{x}{a vector or a factor.}
  19   \item{phy}{an object of class \code{"phylo"}.}
  20   \item{type}{the variable type; either \code{"continuous"} or
  21     \code{"discrete"} (or an abbreviation of these).}
  22   \item{method}{a character specifying the method used for
  23     estimation. Three choices are possible: \code{"ML"}, \code{"pic"},
  24     or \code{"GLS"}.}
  25   \item{CI}{a logical specifying whether to return the 95\% confidence
  26     intervals of the ancestral state estimates (for continuous
  27     characters) or the likelihood of the different states (for discrete
  28     ones).}
  29   \item{model}{a character specifying the model (ignored if \code{method
  30       = "GLS"}), or a numeric matrix if \code{type = "discrete"} (see
  31     details).}
  32   \item{scaled}{a logical specifying whether to scale the contrast
  33     estimate (used only if \code{method = "pic"}).}
  34   \item{kappa}{a positive value giving the exponent transformation of
  35     the branch lengths (see details).}
  36   \item{corStruct}{if \code{method = "GLS"}, specifies the correlation
  37     structure to be used (this also gives the assumed model).}
  38   \item{ip}{the initial value(s) used for the ML estimation procedure
  39     when \code{type == "discrete"} (possibly recycled).}
  40   \item{object}{an object of class \code{"ace"}.}
  41   \item{k}{a numeric value giving the penalty per estimated parameter;
  42     the default is \code{k = 2} which is the classical Akaike
  43     information criterion.}
  44   \item{\dots}{further arguments passed to or from other methods.}
  45 }
  46 \description{
  47   This function estimates ancestral character states, and the associated
  48   uncertainty, for continuous and discrete characters.
  49
  50   \code{logLik}, \code{deviance}, and \code{AIC} are generic functions
  51   used to extract the log-likelihood, the deviance (-2*logLik), or the
  52   Akaike information criterion of a tree. If no such values are
  53   available, \code{NULL} is returned.
  54
  55   \code{anova} is another generic function that is used to compare
  56   nested models: the significance of the additional parameter(s) is
  57   tested with likelihood ratio tests. You must ensure that the models
  58   are effectively nested (if they are not, the results will be
  59   meaningless). It is better to list the models from the smallest to the
  60   largest.
  61 }
  62 \details{
  63   If \code{type = "continuous"}, the default model is Brownian motion
  64   where characters evolve randomly following a random walk. This model
  65   can be fitted by maximum likelihood (the default, Schluter et
  66   al. 1997), least squares (\code{method = "pic"}, Felsenstein 1985), or
  67   generalized least squares (\code{method = "GLS"}, Martins and Hansen
  68   1997, Cunningham et al. 1998). In the latter case, the specification
  69   of \code{phy} and \code{model} are actually ignored: it is instead
  70   given through a correlation structure with the option \code{corStruct}.
  71
  72   For discrete characters (\code{type = "discrete"}), only maximum
  73   likelihood estimation is available (Pagel 1994). The model is
  74   specified through a numeric matrix with integer values taken as
  75   indices of the parameters. The numbers of rows and of columns of this
  76   matrix must be equal, and are taken to give the number of states of
  77   the character. For instance, \code{matrix(c(0, 1, 1, 0), 2)} will
  78   represent a model with two character states and equal rates of
  79   transition, \code{matrix(c(0, 1, 2, 0), 2)} a model with unequal
  80   rates, \code{matrix(c(0, 1, 1, 1, 0, 1, 1, 1, 0), 3)} a model with
  81   three states and equal rates of transition (the diagonal is always
  82   ignored). There are short-cuts to specify these models: \code{"ER"} is
  83   an equal-rates model (e.g., the first and third examples above),
  84   \code{"ARD"} is an all-rates-different model (the second example), and
  85   \code{"SYM"} is a symmetrical model (e.g., \code{matrix(c(0, 1, 2, 1,
  86     0, 3, 2, 3, 0), 3)}). If a short-cut is used, the number of states
  87   is determined from the data.
  88 }
  89 \value{
  90   a list with the following elements:
  91
  92   \item{ace}{if \code{type = "continuous"}, the estimates of the
  93     ancestral character values.}
  94   \item{CI95}{if \code{type = "continuous"}, the estimated 95\%
  95     confidence intervals.}
  96   \item{sigma2}{if \code{type = "continuous"}, \code{model = "BM"}, and
  97     \code{method = "ML"}, the maximum likelihood estimate of the
  98     Brownian parameter.}
  99   \item{rates}{if \code{type = "discrete"}, the maximum likelihood
 100     estimates of the transition rates.}
 101   \item{se}{if \code{type = "discrete"}, the standard-errors of
 102     estimated rates.}
 103   \item{index.matrix}{if \code{type = "discrete"}, gives the indices of
 104     the \code{rates} in the rate matrix.}
 105   \item{loglik}{if \code{method = "ML"}, the maximum log-likelihood.}
 106   \item{lik.anc}{if \code{type = "discrete"}, the scaled likelihoods of
 107     each ancestral state.}
 108   \item{call}{the function call.}
 109 }
 110 \references{
 111   Cunningham, C. W., Omland, K. E. and Oakley, T. H. (1998)
 112   Reconstructing ancestral character states: a critical
 113   reappraisal. \emph{Trends in Ecology & Evolution}, \bold{13},
 114   361--366.
 115
 116   Felsenstein, J. (1985) Phylogenies and the comparative
 117   method. \emph{American Naturalist}, \bold{125}, 1--15.
 118
 119   Martins, E. P. and Hansen, T. F. (1997) Phylogenies and the
 120   comparative method: a general approach to incorporating phylogenetic
 121   information into the analysis of interspecific data. \emph{American
 122     Naturalist}, \bold{149}, 646--667.
 123
 124   Pagel, M. (1994) Detecting correlated evolution on phylogenies: a
 125   general method for the comparative analysis of discrete
 126   characters. \emph{Proceedings of the Royal Society of London. Series
 127     B. Biological Sciences}, \bold{255}, 37--45.
 128
 129   Schluter, D., Price, T., Mooers, A. O. and Ludwig, D. (1997)
 130   Likelihood of ancestor states in adaptive radiation. \emph{Evolution},
 131   \bold{51}, 1699--1711.
 132 }
 133 \author{Emmanuel Paradis, Ben Bolker \email{bolker@zoo.ufl.edu}}
 134 \seealso{
 135   \code{\link{corBrownian}}, \code{\link{corGrafen}},
 136   \code{\link{corMartins}}, \code{\link{compar.ou}},
 137   \code{\link[stats]{anova}}
 138 }
 139 \examples{
 140 ### Just some random data...
 141 data(bird.orders)
 142 x <- rnorm(23)
 143 ### Compare the three methods for continuous characters:
 144 ace(x, bird.orders)
 145 ace(x, bird.orders, method = "pic")
 146 ace(x, bird.orders, method = "GLS",
 147     corStruct = corBrownian(1, bird.orders))
 148 ### For discrete characters:
 149 x <- factor(c(rep(0, 5), rep(1, 18)))
 150 ans <- ace(x, bird.orders, type = "d")
 151 #### Showing the likelihoods on each node:
 152 plot(bird.orders, type = "c", FALSE, label.offset = 1)
 153 co <- c("blue", "yellow")
 154 tiplabels(pch = 22, bg = co[as.numeric(x)], cex = 2, adj = 1)
 155 nodelabels(thermo = ans$lik.anc, piecol = co, cex = 0.75)
 156 ### An example of the use of the argument `ip':
 157 tr <- character(4)
 158 tr[1] <- "((((t10:5.03,t2:5.03):2.74,(t9:4.17,"
 159 tr[2] <- "t5:4.17):3.60):2.80,(t3:4.05,t7:"
 160 tr[3] <- "4.05):6.53):2.32,((t6:4.38,t1:4.38):"
 161 tr[4] <- "2.18,(t8:2.17,t4:2.17):4.39):6.33);"
 162 tr <- read.tree(text = paste(tr, collapse = ""))
 163 y <- c(rep(1, 6), rep(2, 4))
 164 ### The default `ip = 0.1' makes ace fails:
 165 ace(y, tr, type = "d")
 166 ace(y, tr, type = "d", ip = 0.01)
 167 ### Surprisingly, using an initial value farther to the
 168 ### MLE than the default one works:
 169 ace(y, tr, type = "d", ip = 0.3)
 170 }
 171 \keyword{models}