https://github.com/cran/aster
Raw File
Tip revision: 303d520fe57883772999cb6e59e5ce81bb6e2741 authored by Charles J. Geyer on 23 November 2005, 00:00:00 UTC
version 0.4-1
Tip revision: 303d520
aster.Rd
\name{aster}
\alias{aster}
\alias{aster.default}
\alias{aster.formula}
\title{Aster Models}
\concept{regression}
\concept{exponential family}
\concept{graphical model}
\description{
  Fits Aster Models.
}
\usage{
aster(x, \dots)

\method{aster}{default}(x, root, pred, fam, modmat, parm,
    type = c("unconditional", "conditional"),
    method = c("trust", "nlm", "CG", "L-BFGS-B"), fscale, maxiter = 1000,
    nowarn = TRUE, newton = TRUE, optout = FALSE, coef.names, \dots)

\method{aster}{formula}(formula, pred, fam, varvar, idvar, root,
    data, parm, type = c("unconditional", "conditional"),
    method = c("trust", "nlm", "CG", "L-BFGS-B"), fscale, maxiter = 1000,
    nowarn = TRUE, newton = TRUE, optout = FALSE, \dots)
}
\arguments{

  \item{x}{an \code{nind} by \code{nnode} matrix, the data for an
    aster model.  The rows are independent and identically modeled
    random vectors.  See details below for further requirements.

    \code{aster.formula} constructs such an \code{x} from the response
    in its formula.  Hence data for \code{aster.formula} must have
    \code{nind * nnode} rows.}

  \item{root}{an object of the same shape as \code{x}, the root data.
    For \code{aster.default} an \code{nind} by \code{nnode} matrix,
    For \code{aster.formula} an \code{nind * nnode} vector.}

  \item{pred}{an integer vector of length \code{nnode} determining
    the dependence
    graph of the aster model.  \code{pred[j]} is
    the index of the predecessor of
    the node with index \code{j} unless the predecessor is a root
    node, in which case \code{pred[j] == 0}.  See details below for
    further requirements.}

  \item{fam}{an integer vector of length \code{nnode} determining
    the exponential family structure of the aster model.
    \code{fam} is an index into the vector of family names returned by
    \code{\link{families}}, so \code{families()[fam]} gives the
    one-parameter exponential families associated with each node
    of the graph.}

  \item{modmat}{an \code{nind} by \code{nnode} by \code{ncoef}
    three-dimensional array, the model matrix.

    \code{aster.formula} constructs such a \code{modmat} from 
    its formula, the data frame \verb@data@, and the variables
    in the environment of the formula.}

  \item{parm}{usually missing.  Otherwise a vector of length \code{ncoef}
    giving a starting point for the optimization.}

  \item{type}{type of model.
    The value of this argument can be abbreviated.}

  \item{method}{optimization method.  If \code{"trust"} then the
    \code{\link[trust]{trust}} function is used.  If \code{"nlm"} then the
    \code{\link{nlm}} function is used.  Otherwise the
    \code{\link{optim}} function is used with the specified \code{method}
    supplied to it.
    The value of this argument can be abbreviated.}

  \item{fscale}{an estimate of the size of the log likelihood at the maximum.
    Defaults to \code{nind}.}

  \item{maxiter}{maximum number of iterations.  Defaults to '1000'.}

  \item{nowarn}{if \code{TRUE} (the default), suppress warnings from
    the optimization routine.}

  \item{newton}{if \code{TRUE} (the default), do one Newton iteration
    on the result produced by the optimization routine, except when
    \verb@method == "trust"@ when no such Newton iteration is done,
    regardless of the value of \verb@newton@, because \code{\link[trust]{trust}}
    always terminates with a Newton iteration when it converges.}

  \item{optout}{if \code{TRUE}, save the entire result of the optimization
    routine (\code{link[trust]{trust}}, \code{link{nlm}}, or \code{link{optim}},
    as the case may be).}

  \item{coef.names}{names of the regression coefficients.  If missing,
    \code{dimnames(modmat)[[3]]} is used.  In \verb@aster.formula@ these
    are produced automatically by the R formula machinery.}

  \item{\ldots}{other arguments passed to the optimization method.}

  \item{formula}{a symbolic description of the model to be fit.  See
    \code{\link{lm}}, \code{\link{glm}}, and 
    \code{\link{formula}} for discussions of the R formula mini-language.}

  \item{varvar}{a variable of the same length as the response in
    the formula that is a factor whose levels are character strings
    treated as variable names.  The number of variable names is \code{nnode}.
    Must be of the form \code{rep(vars, each = nind)} where \code{vars} is
    a vector of variable names.  Usually found in the data frame \code{data}
    when this is produced by the \code{\link{reshape}} function.}

  \item{idvar}{a variable of the same length as the response in
    the formula that indexes individuals.  The number
    of individuals is \code{nind}.
    Must be of the form \code{rep(inds, times = nnode)} where \code{inds} is
    a vector of labels for individuals.  Usually found in the data frame
    \code{data} when this is produced by the \code{\link{reshape}} function.}

  \item{data}{an optional data frame containing the variables
    in the model.  If not found in \code{data}, the variables are taken
    from \code{environment(formula)}, typically the environment from
    which \code{aster} is called.  Usually produced by
    the \code{\link{reshape}} function.}
}
\details{
The vector \code{pred} must satisfy \code{all(pred < seq(along = pred))},
that is, each predecessor must precede in the order given in pred.
The vector \code{pred} defines a function \eqn{p}.

The joint distribution of the data matrix \code{x} is a product of conditionals
\deqn{\prod_{i \in I} \prod_{j \in J} \Pr \{ X_{i j} | X_{i p(j)} \}}{prod[i] prod[j] Pr(x[i, j], x[i, p(j)])}
When \eqn{p(j) = 0}, the notation \eqn{X_{i p(j)}}{x[i, p(j)]}
means \code{root[i, j]}.  Other elements of the matrix \code{root} are
not used.

The conditional distribution
\eqn{\Pr \{ X_{i j} | X_{i p(j)} \}}{Pr(x[i, j], x[i, p(j)])}
is the \eqn{X_{i p(j)}}{x[i, p(j)]}-fold convolution of the \eqn{j}-th family
in the vector \code{fam}, a one-parameter exponential family
(i.e., the sum of \eqn{X_{i p(j)}}{x[i, p(j)]} i.i.d. terms having
this one-parameter exponential family distribution).

For \code{type == "conditional"} the canonical parameter vector
\eqn{\theta_{i j}}{theta[i, j]} is modeled in GLM fashion as
\eqn{\theta = M \beta}{theta = M beta} where \eqn{M} is the model
matrix \code{modmat}.  Since the \dQuote{vector} \eqn{\theta}{theta} is
actually a matrix, the \dQuote{matrix} \eqn{M} must correspondingly
be a three-dimensional array.  So \eqn{\theta = M \beta}{theta = M beta}
written out in full is
\deqn{\theta_{i j} = \sum_{k \in K} m_{i j k} \beta_k}{theta[i, j] = sum[k] m[i, j, k] beta[k]}
This specifies the log likelihood.

For \code{type == "unconditional"} the canonical parameter vector
for an unconditional model is modeled in GLM fashion as
\eqn{\varphi = M \beta}{phi = M beta} (where the notation is as above).
The unconditional canonical parameters are then specified in terms of
the conditional ones by
\deqn{\varphi_{i j} = \theta_{i j} - \sum_{k \in S(j)} \psi_k(\theta_{i k})}{phi[i, j] = theta[i, j] - sum[k in S(j)] psi[k](theta[i, k])}
where \eqn{S(j)} denotes the set of successors of \eqn{j},
the \eqn{k} such that \eqn{p(k) = j}, and \eqn{\psi_k}{psi[k]} is the
cumulant function for the \eqn{k}-th exponential family.
This rather crazy looking formulation is an invertible change of parameter
and makes \eqn{\varphi}{phi}
the canonical parameter and \eqn{x} the canonical statistic of a full
flat unconditional exponential family.
Again, this specifies the log likelihood.

In order to use the R model formula mini-language we must flatten
the dimensionality, making the model matrix \code{modmat} two-dimensional
(a true matrix).  This must be done as if by
\code{matrix(modmat, ncol = ncoef)},
which imposes the requirements on \code{varvar} and \code{idvar}
given in the arguments section: they must look like \code{row(x)} and
\code{col(x)} modulo relabeling.
Then \code{x} and \code{root}
become one-dimensional, done as if by \code{as.numeric(x)}
and \code{as.numeric(root)}.

The standard way to do this in R is to use the \code{\link{reshape}}
function on a data frame in which the columns of the \code{x} matrix
are variables in the data frame.  \code{reshape} automatically puts
things in the right order and creates \code{varvar} and \code{idvar}.
}
\value{
  \code{aster} returns an object of class inheriting from \code{"aster"}.
  \code{aster.formula}, returns an object of class \code{"aster"} and
  subclass \code{"aster.formula"}.

  The function \code{\link{summary}} (i.e., \code{\link{summary.aster}}) can
  be used to obtain or print a summary of the results, the function
  \code{\link{anova}} (i.e., \code{\link{anova.aster}})
  to produce an analysis of deviance table, and the function
  \code{\link{predict}} (i.e., \code{\link{predict.aster}})
  to produce predicted values and standard errors.

  An object of class \code{"aster"} is a list containing at least the
  following components:

  \item{coefficients}{a named vector of coefficients.}
  \item{rank}{the numeric rank of the fitted generalized linear model
      part of the aster model (i.e., the rank of \code{modmat}).}
  \item{deviance}{up to a constant, minus twice the maximized
    log-likelihood.}
  \item{iter}{the number of iterations used by the optimization method.}
  \item{converged}{logical. Was the optimization algorithm judged to have
     converged?}
  \item{code}{integer.  The convergence code returned by the optimization
     method.}
  \item{gradient}{The gradient vector of minus the log likelihood at the
      fitted \code{coefficients} vector.}
  \item{hessian}{The Hessian matrix of minus the log likelihood
      (i.e., the observed Fisher information) at the
      fitted \code{coefficients} vector.
      This is also the expected Fisher information when
      \code{type == "unconditional"}.}
  \item{fisher}{Expected Fisher information at the fitted \code{coefficients}
      vector.}
  \item{optout}{The object returned by the optimization routine
    (\code{link[trust]{trust}}, \code{link{nlm}}, or \code{link{optim}}).
    Only returned when the argument \code{optout} is \code{TRUE}.}

  Calls to \code{aster.formula} return a list also containing:

  \item{call}{the matched call.}
  \item{formula}{the formula supplied.}
  \item{terms}{the \code{\link{terms}} object used.}
  \item{data}{the \code{data argument}.}
}
\seealso{
\code{\link{anova.aster}},
\code{\link{summary.aster}},
and
\code{\link{predict.aster}}
}
\examples{
### see package vignette for explanation ###
library(aster)
data(echinacea)
vars <- c("ld02", "ld03", "ld04", "fl02", "fl03", "fl04",
    "hdct02", "hdct03", "hdct04")
redata <- reshape(echinacea, varying = list(vars), direction = "long",
    timevar = "varb", times = as.factor(vars), v.names = "resp")
redata <- data.frame(redata, root = 1)
pred <- c(0, 1, 2, 1, 2, 3, 4, 5, 6)
fam <- c(1, 1, 1, 1, 1, 1, 3, 3, 3)
hdct <- grep("hdct", as.character(redata$varb))
hdct <- is.element(seq(along = redata$varb), hdct)
redata <- data.frame(redata, hdct = as.integer(hdct))
aout4 <- aster(resp ~ varb + nsloc + ewloc + pop * hdct - pop,
    pred, fam, varb, id, root, data = redata)
summary(aout4, show.graph = TRUE)
}
\keyword{models}
\keyword{regression}
back to top