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{...}{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). In the latter case, the specification of \code{phy} and
  69   \code{model} are actually ignored: it is instead given through a
  70   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   Felsenstein, J. (1985) Phylogenies and the comparative
 112   method. \emph{American Naturalist}, \bold{125}, 1--15.
 113
 114   Martins, E. P. and Hansen, T. F. (1997) Phylogenies and the
 115   comparative method: a general approach to incorporating phylogenetic
 116   information into the analysis of interspecific data. \emph{American
 117     Naturalist}, \bold{149}, 646--667.
 118
 119   Pagel, M. (1994) Detecting correlated evolution on phylogenies: a
 120   general method for the comparative analysis of discrete
 121   characters. \emph{Proceedings of the Royal Society of London. Series
 122     B. Biological Sciences}, \bold{255}, 37--45.
 123
 124   Schluter, D., Price, T., Mooers, A. O. and Ludwig, D. (1997)
 125   Likelihood of ancestor states in adaptive radiation. \emph{Evolution},
 126   \bold{51}, 1699--1711.
 127 }
 128 \author{Emmanuel Paradis \email{Emmanuel.Paradis@mpl.ird.fr}, Ben Bolker
 129 \email{bolker@zoo.ufl.edu}}
 130 \seealso{
 131   \code{\link{corBrownian}}, \code{\link{corGrafen}},
 132   \code{\link{corMartins}}, \code{\link{compar.ou}},
 133   \code{\link[stats]{anova}}
 134 }
 135 \examples{
 136 ### Just some random data...
 137 data(bird.orders)
 138 x <- rnorm(23)
 139 ### Compare the three methods for continuous characters:
 140 ace(x, bird.orders)
 141 ace(x, bird.orders, method = "pic")
 142 ace(x, bird.orders, method = "GLS",
 143     corStruct = corBrownian(1, bird.orders))
 144 ### For discrete characters:
 145 x <- factor(c(rep(0, 5), rep(1, 18)))
 146 ans <- ace(x, bird.orders, type = "d")
 147 #### Showing the likelihoods on each node:
 148 plot(bird.orders, type = "c", FALSE, label.offset = 1)
 149 co <- c("blue", "yellow")
 150 tiplabels(pch = 22, bg = co[as.numeric(x)], cex = 2, adj = 1)
 151 nodelabels(thermo = ans$lik.anc, piecol = co, cex = 0.75)
 152 ### An example of the use of the argument `ip':
 153 tr <- character(4)
 154 tr[1] <- "((((t10:5.03,t2:5.03):2.74,(t9:4.17,"
 155 tr[2] <- "t5:4.17):3.60):2.80,(t3:4.05,t7:"
 156 tr[3] <- "4.05):6.53):2.32,((t6:4.38,t1:4.38):"
 157 tr[4] <- "2.18,(t8:2.17,t4:2.17):4.39):6.33);"
 158 tr <- read.tree(text = paste(tr, collapse = ""))
 159 y <- c(rep(1, 6), rep(2, 4))
 160 ### The default `ip = 0.1' makes ace fails:
 161 ace(y, tr, type = "d")
 162 ace(y, tr, type = "d", ip = 0.01)
 163 ### Surprisingly, using an initial value farther to the
 164 ### MLE than the default one works:
 165 ace(y, tr, type = "d", ip = 0.3)
 166 }
 167 \keyword{models}