\name{abvnonpar}
\alias{abvnonpar}
\title{Non-parametric Estimates for Dependence Functions of the
    Bivariate Extreme Value Distribution}
\description{
  Calculate or plot non-parametric estimates for the dependence function
  \eqn{A} of the bivariate extreme value distribution.
}
\synopsis{
abvnonpar(x = 0.5, data, nsloc1 = NULL, nsloc2 = NULL,
    method = c("cfg", "pickands", "deheuvels", "hall", "tdo"),
    convex = FALSE, wf = function(t) t, kmar = NULL, plot = FALSE,
    add = FALSE, lty = 1, lwd = 1, col = 1, blty = 3, xlim = c(0, 1),
    ylim = c(0.5, 1), xlab = "", ylab = "", \dots)
}
\usage{
abvnonpar(x = 0.5, data, nsloc1 = NULL, nsloc2 = NULL,
    method = c("cfg", "pickands", "deheuvels", "hall", "tdo"),
    convex = FALSE, wf = function(t) t, plot = FALSE,
    add = FALSE, lty = 1, lwd = 1, col = 1, blty = 3, xlim = c(0, 1),
    ylim = c(0.5, 1), xlab = "", ylab = "", \dots)
}
\arguments{
  \item{x}{A vector of values at which the dependence function is
    evaluated (ignored if plot or add is \code{TRUE}). \eqn{A(1/2)}
    is returned by default since it is often a useful summary of
    dependence.}
  \item{data}{A matrix or data frame with two columns, which may contain
    missing values.}
  \item{nsloc1, nsloc2}{A data frame with the same number of rows as
    \code{data}, for linear modelling of the location
    parameter on the first/second margin.
    The data frames are treated as covariate matrices, excluding the
    intercept.
    A numeric vector can be given as an alternative to a single column
    data frame.}
  \item{method}{The estimation method (see \bold{Details}). The
    \code{"cfg"} method is used by default.}
  \item{convex}{Logical; take the convex minorant?}
  \item{wf}{The weight function used in the \code{"cfg"}
    method (see \bold{Details}). The function must be vectorized.}
  \item{plot}{Logical; if \code{TRUE} the function is plotted. The
    x and y values used to create the plot are returned invisibly.
    If \code{plot} and \code{add} are \code{FALSE} (the default),
    the arguments following \code{add} are ignored.}
  \item{add}{Logical; add to an existing plot? The existing plot
    should have been created using either \code{abvnonpar} or
    \code{\link{abvpar}}, the latter of which plots (or calculates)
    the dependence function for a number of parametric models.}
  \item{lty, blty}{Function and border line types. Set \code{blty}
    to zero to omit the border.}
  \item{lwd}{Line width.}
  \item{col}{Line colour.}
  \item{xlim, ylim}{x and y-axis limits.}
  \item{xlab, ylab}{x and y-axis labels.}
  \item{\dots}{Other high-level graphics parameters to be passed to
    \code{plot}.}
}
\details{
  The dependence function  \eqn{A(\cdot)}{A()} of the bivariate
  extreme value distribution is defined in \code{\link{abvpar}}.
  Non-parametric estimates are constructed as follows.
  Suppose \eqn{(z_{i1},z_{i2})} for \eqn{i=1,\ldots,n} are \eqn{n}
  bivariate observations that are passed using the \code{data}
  argument.
  The marginal parameters are estimated (under the assumption of
  independence) and the data is transformed using
  \deqn{y_{i1} = \{1+\hat{s}_1(z_{i1}-\hat{a}_1)/
    \hat{b}_1\}_{+}^{-1/\hat{s}_1}}{
        y_{i1} = {1 + s'_1(z_{i1}-a'_1)/b'_1}^(-1/s'_1)}
  and
  \deqn{y_{i2} = \{1+\hat{s}_2(z_{i2}-\hat{a}_2)/
    \hat{b}_2\}_{+}^{-1/\hat{s}_2}}{
        y_{i2} = {1 + s'_2(z_{i2}-a'_2)/b'_2}^(-1/s'_2)}
  for \eqn{i = 1,\ldots,n}, where
  \eqn{(\hat{a}_1,\hat{b}_1,\hat{s}_1)}{(a'_1,b'_1,s'_1)} and
  \eqn{(\hat{a}_2,\hat{b}_2,\hat{s}_2)}{(a'_2,b'_2,s'_2)}
  are the maximum likelihood estimates for the location, scale
  and shape parameters on the first and second margins.
  If \code{nsloc1} or \code{nsloc2} are given, the location
  parameters may depend on \eqn{i} (see \code{\link{fgev}}).

  Five different estimators of the dependence function can be
  implemented.
  They are defined (on \eqn{0 \leq w \leq 1}{0 <= w <= 1}) as
  follows.

  \code{method = "cfg"} (Caperaa, Fougeres and Genest, 1997)
  \deqn{A_c(w) = \exp\left\{ \{1-p(w)\} \int_{0}^{w}
    \frac{H(x) - x}{x(1-x)} \, \mbox{d}x - p(w) \int_{w}^{1}
    \frac{H(x) - x}{x(1-x)} \, \mbox{d}x \right\}}{
    A_c(w) = exp{ [1-p(w)] integral_0^w (H(x) - x)/[x(1-x)] dx -
      p(w) integral_w^1 (H(x) - x)/[x(1-x)] dx }}
  
  \code{method = "pickands"} (Pickands, 1981)
  \deqn{A_p(w) = n\left\{\sum_{i=1}^n \min\left(\frac{y_{i1}}{w},
    \frac{y_{i2}}{1-w}\right)\right\}^{-1}}{
    A_p(w) = n / {sum_{i=1}^n min[y_{i1}/w, y_{i2}/(1-w)]}}

  \code{method = "deheuvels"} (Deheuvels, 1991)
  \deqn{A_d(w) = n\left\{\sum_{i=1}^n \min\left(\frac{y_{i1}}{w},
    \frac{y_{i2}}{1-w}\right) - w\sum_{i=1}^n y_{i1} - (1-w)
    \sum_{i=1}^n y_{i2} + n\right\}^{-1}}{
        A_d(w) = n / {sum_{i=1}^n min[y_{i1}/w, y_{i2}(1-w)] -
	w sum_{i=1}^n y_{i1} - (1-w) sum_{i=1}^n y_{i2} + n}}

  \code{method = "hall"} (Hall and Tajvidi, 2000)
  \deqn{A_h(w) = n\left\{\sum_{i=1}^n \min\left(\frac{y_{i1}}{\bar{y}_1
      w},\frac{y_{i2}}{\bar{y}_2 (1-w)}\right)\right\}^{-1}}{
    A_h(w) = n (sum_{i=1}^n min[y_{i1}/(by_1 w),
    y_{i2}/(by_2 (1-w))])^{-1}}

  \code{method = "tdo"} (Tiago de Oliveira, 1997)
  \deqn{A_t(w) = 1 - \frac{1}{1 + \log n} \sum_{i=1}^n
    \min\left(\frac{w}{1+ny_{i1}},\frac{1-w}{1+ny_{i2}}\right)}{
    A_t(w) = 1 - 1/(1 + log n) sum_{i=1}^n
    min[w/(1 + n y_{i1}), (1 - w)/(1 + n y_{i2})]}
    
  In the estimator \eqn{A_h(\cdot)}{A_h()},
  \eqn{\bar{y}_j = n^{-1}\sum_{i=1}^n y_{ij}}{
    by_j = (sum_{i=1}^n y_{ij})/n} for \eqn{j = 1,2}.
  In the estimator \eqn{A_c(\cdot)}{A_c()}, \eqn{H(x)} is the
  empirical distribution function of \eqn{x_1,\ldots,x_n}, where
  \eqn{x_i = y_{i1} / (y_{i1} + y_{i2})} for \eqn{i = 1,\ldots,n},
  and \eqn{p(w)} is any bounded function on \eqn{[0,1]}, which
  can be specified using the argument \code{wf}.
  By default \code{wf} is the identity function.

  Let \eqn{A_n(\cdot)}{A_n()} be any estimator of \eqn{A(\cdot)}{A()}.
  The constraint \eqn{A_n(0) = A_n(1) = 1} is satisfied by
  \eqn{A_d(\cdot)}{A_d()}, \eqn{A_t(\cdot)}{A_t()} and
  \eqn{A_h(\cdot)}{A_h()}, and by \eqn{A_c(\cdot)}{A_c()} when
  \eqn{p(0) = 0} and \eqn{p(1) = 1}.
  None of the estimators satisfy
  \eqn{\max(w,1-w) \leq A_n(w) \leq 1}{max(w,1-w) <= A_n(w) <= 1}
  for all \eqn{0\leq w \leq1}{0 <= w <= 1}. An obvious modification is
  \deqn{A_n^{'}(w) = \min(1, \max\{A_n(w), w, 1-w\}).}{
    A'_n(w) = min(1, max{A_n(w), w, 1-w}).}
  This modification is always implemented.

  \eqn{A_t(w)}{A_t(w)} is the only estimator that is convex.
  Convex estimators can be derived from other methods by taking
  the convex minorant, which can be achieved by setting \code{convex}
  to \code{TRUE}.
}
\note{
  Appendix A of the User's Guide contains a short simulation
  study that compares the estimators defined above.
  The estimators \eqn{A_p(\cdot)}{A_p()}, \eqn{A_d(\cdot)}{A_d()}
  and \eqn{A_h(\cdot)}{A_h()} are very similar, and may not be
  distinguishable when plotted.
}

\value{
  \code{abvnonpar} calculates or plots a non-parametric estimate of
  the dependence function of the bivariate extreme value distribution.
}
\references{
  Caperaa, P. Fougeres, A.-L. and Genest, C. (1997)
  A non-parametric estimation procedure for bivariate
  extreme value copulas.
  \emph{Biometrika}, \bold{84}, 567--577.

  Deheuvels, P. (1991)
  On the limiting behaviour of the Pickands estimator 
  for bivariate extreme-value distributions.
  \emph{Statist. Probab. Letters}, \bold{12}, 429--439.

  Hall, P. and Tajvidi, N. (2000)
  Distribution and dependence-function estimation for
  bivariate extreme-value distributions.
  \emph{Bernoulli}, \bold{6}, 835--844.

  Pickands, J. (1981)
  Multivariate extreme value distributions.
  \emph{Proc. 43rd Sess. Int. Statist. Inst.}, \bold{49}, 859--878.

  Tiago de Oliveira, J. (1997)
  \emph{Statistical Analysis of Extremes.} Pendor.
}

\seealso{\code{\link{abvpar}}, \code{\link{atvnonpar}},
  \code{\link{fgev}}}

\examples{
bvdata <- rbvevd(100, dep = 0.7, model = "log")
abvnonpar(seq(0, 1, length = 10), data = bvdata, convex = TRUE)
abvnonpar(data = bvdata, method = "d", plot = TRUE)

M1 <- fitted(fbvevd(bvdata, model = "log"))
abvpar(dep = M1["dep"], model = "log", plot = TRUE)
abvnonpar(data = bvdata, add = TRUE, lty = 2)
}
\keyword{nonparametric}









\eof
\name{abvpar}
\alias{abvpar}
\title{Parametric Dependence Functions of Bivariate Extreme
    Value Models}
\description{
  Calculate or plot the dependence function \eqn{A} for
  eight parametric bivariate extreme value models.
}
\synopsis{
abvpar(x = 0.5, dep, asy = c(1,1), alpha, beta, model = c("log", "alog",
    "hr", "neglog", "aneglog", "bilog", "negbilog", "ct"),
     plot = FALSE, add = FALSE, lty = 1, lwd = 1, col = 1, blty = 3,
     xlim = c(0,1), ylim = c(0.5,1), xlab = "", ylab = "", \dots)
}
\usage{
abvpar(x = 0.5, dep, asy = c(1,1), alpha, beta, model = "log",
     plot = FALSE, add = FALSE, lty = 1, lwd = 1, col = 1, blty = 3,
     xlim = c(0,1), ylim = c(0.5,1), xlab = "", ylab = "", \dots)
}
\arguments{
  \item{x}{A vector of values at which the dependence function is
    evaluated (ignored if plot or add is \code{TRUE}). \eqn{A(1/2)}
    is returned by default since it is often a useful summary of
    dependence.}
  \item{dep}{Dependence parameter for the logistic, asymmetric
    logistic, Husler-Reiss, negative logistic and asymmetric
    negative logistic models.}
  \item{asy}{A vector of length two, containing the two asymmetry
    parameters for the asymmetric logistic and asymmetric negative
    logistic models.}
  \item{alpha, beta}{Alpha and beta parameters for the bilogistic,
    negative bilogistic and Coles-Tawn models.}
  \item{model}{The specified model; a character string. Must be
    either \code{"log"} (the default), \code{"alog"}, \code{"hr"},
    \code{"neglog"}, \code{"aneglog"}, \code{"bilog"},
    \code{"negbilog"} or \code{"ct"} (or any unique partial match),
    for the logistic, asymmetric logistic, Husler-Reiss, negative
    logistic, asymmetric negative logistic, bilogistic, negative
    bilogistic and Coles-Tawn models respectively. The definition
    of each model is given in \code{\link{rbvevd}}. If parameter
    arguments are given that do not correspond to the specified
    model those arguments are ignored, with a warning.}
  \item{plot}{Logical; if \code{TRUE} the function is plotted. The
    x and y values used to create the plot are returned invisibly.
    If \code{plot} and \code{add} are \code{FALSE} (the default),
    the arguments following \code{add} are ignored.}
  \item{add}{Logical; add to an existing plot? The existing plot
    should have been created using either \code{abvpar} or
    \code{\link{abvnonpar}}, the latter of which plots (or calculates)
    a non-parametric estimate of the dependence function.}
  \item{lty, blty}{Function and border line types. Set \code{blty}
    to zero to omit the border.}
  \item{lwd}{Line width.}
  \item{col}{Line colour.}
  \item{xlim, ylim}{x and y-axis limits.}
  \item{xlab, ylab}{x and y-axis labels.}
  \item{\dots}{Other high-level graphics parameters to be passed to
    \code{plot}.}
}
\details{
  Any bivariate extreme value distribution can be written as
  \deqn{G(z_1,z_2) = \exp\left[-(y_1+y_2)A\left(
    \frac{y_1}{y_1+y_2}\right)\right]}{
    G(z1,z2) = exp{-(y1+y2)A[y1/(y1+y2)]}}
  for some function \eqn{A(\cdot)}{A()} defined on \eqn{[0,1]}, where    
  \deqn{y_i = \{1+s_i(z_i-a_i)/b_i\}^{-1/s_i}}{
    yi = {1+si(zi-ai)/bi}^(-1/si)}
  for \eqn{1+s_i(z_i-a_i)/b_i > 0}{1+si(zi-ai)/bi > 0} and
  \eqn{i = 1,2}, with the (generalized extreme value) marginal
  parameters given by \eqn{(a_i,b_i,s_i)}{(ai,bi,si)},
  \eqn{b_i > 0}{bi > 0}.
  If \eqn{s_i = 0}{si = 0} then \eqn{y_i}{yi} is defined by
  continuity.

  \eqn{A(\cdot)}{A()} is called (by some authors) the dependence
  function.
  It follows that \eqn{A(0)=A(1)=1}, and that \eqn{A(\cdot)}{A()} is
  a convex function with \eqn{\max(x,1-x) \leq A(x)\leq 1}{max(x,1-x)
    <= A(x) <= 1} for all \eqn{0\leq x\leq1}{0 <= x <= 1}.
  The lower and upper limits of \eqn{A} are obtained under complete
  dependence and independence respectively.
  \eqn{A(\cdot)}{A()} does not depend on the marginal parameters.
}

\value{
  \code{abvpar} calculates or plots the dependence function
  for one of eight parametric bivariate extreme value models,
  at specified parameter values.
}

\seealso{\code{\link{abvnonpar}}, \code{\link{fbvevd}},
  \code{\link{rbvevd}}, \code{\link{atvpar}}}

\examples{
abvpar(dep = 2.7, model = "hr")
abvpar(seq(0,1,0.25), dep = 0.3, asy = c(.7,.9), model = "alog")
abvpar(alpha = 0.3, beta = 1.2, model = "negbi", plot = TRUE)

bvdata <- rbvevd(100, dep = 0.7, model = "log")
M1 <- fitted(fbvevd(bvdata, model = "log"))
abvpar(dep = M1["dep"], model = "log", plot = TRUE)
abvnonpar(data = bvdata, add = TRUE, lty = 2)
}
\keyword{distribution}









\eof
\name{anova.evd}
\alias{anova.evd}

\title{Compare Nested EVD Objects}
\description{
  Compute an analysis of deviance table for two or more nested
  evd objects.
}
\usage{
\method{anova}{evd}(object, object2, \dots)
}
\arguments{
  \item{object}{An object of class \code{"evd"}.}
  \item{object2}{An object of class \code{"evd"} that
    represents a model nested within \code{object}.}
  \item{\dots}{Further successively nested objects.}
}
\value{
  An object of class \code{c("anova", "data.frame")}, with one
  row for each model, and the following five columns
  \item{M.Df}{The number of parameters.}
  \item{Deviance}{The deviance.}
  \item{Df}{The number of parameters of the model in the previous row
    minus the number of parameters.}
  \item{Chisq}{The deviance minus the deviance of the model
    in the previous row.}
  \item{Pr(>chisq)}{The p-value calculated by comparing the quantile
    \code{Chisq} with a chi-squared distribution on \code{Df} degrees
    of freedom.}
}
\section{Warning}{
  Circumstances may arise such that the asymptotic distribution of the
  test statistic is not chi-squared.
  In particular, this occurs when the nested
  model is constrained at the edge of the parameter space. It is up
  to the user recognize this, and to interpret the output correctly.
}

\seealso{\code{\link{fbvevd}}, \code{\link{fextreme}},
  \code{\link{fgev}}, \code{\link{forder}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
trend <- (-49:50)/100
M1 <- fgev(uvdata, nsloc = trend)
M2 <- fgev(uvdata)
M3 <- fgev(uvdata, shape = 0)
anova(M1, M2, M3)

bvdata <- rbvevd(100, dep = 0.75, model = "log")
M1 <- fbvevd(bvdata, model = "log")
M2 <- fbvevd(bvdata, model = "log", dep = 0.75)
anova(M1, M2)
}

\keyword{models}


\eof
\name{atvnonpar}
\alias{atvnonpar}
\title{Non-parametric Estimates for Dependence Functions of the
    Trivariate Extreme Value Distribution}
\description{
  Calculate or plot non-parametric estimates for the dependence function
  \eqn{A} of the trivariate extreme value distribution.
}
\synopsis{
atvnonpar(x = rep(1/3,3), data, nsloc1 = NULL, nsloc2 = NULL, nsloc3 = NULL,
    method = c("pickands", "deheuvels", "hall"), kmar = NULL, plot = FALSE,
    col = heat.colors(12), blty = 0, grid = if(blty) 150 else 50,
    lower = 1/3, ord = 1:3, lab = as.character(1:3), lcex = 1)
}
\usage{
atvnonpar(x = rep(1/3,3), data, nsloc1 = NULL, nsloc2 = NULL, nsloc3 = NULL,
    method = c("pickands", "deheuvels", "hall"), plot = FALSE,
    col = heat.colors(12), blty = 0, grid = if(blty) 150 else 50,
    lower = 1/3, ord = 1:3, lab = as.character(1:3), lcex = 1)
}
\arguments{
  \item{x}{A vector of length three or a matrix with three columns,
    in which case the dependence function is evaluated across
    the rows (ignored if plot is \code{TRUE}). The elements/rows
    of the vector/matrix should be positive and should sum to one,
    or else they should have a positive sum, in which case the rows
    are rescaled and a warning is given. \eqn{A(1/3,1/3,1/3)} is
    returned by default since it is often a useful summary of
    dependence.}
  \item{data}{A matrix or data frame with three columns, which may
    contain missing values.}
  \item{nsloc1, nsloc2, nsloc3}{A data frame with the same number of
    rows as \code{data}, for linear modelling of the location
    parameter on the first/second/third margin.
    The data frames are treated as covariate matrices, excluding the
    intercept.
    A numeric vector can be given as an alternative to a single column
    data frame.}
  \item{method}{The estimation method; a character string. Must be
    either \code{"pickands"} (the default), \code{"deheuvels"}
    or \code{"hall"} (or any unique partial match). The three
    estimators are very similar, and may not be distinguishable
    when plotted.}
  \item{plot}{Logical; if \code{TRUE} the function is plotted. The
    minimum (evaluated) value is returned invisibly.
    If \code{FALSE} (the default), the following arguments are
    ignored.}
  \item{col}{A list of colours (see \code{\link{image}}). The first
    colours in the list represent smaller values, and hence
    stronger dependence. Each colour represents an equally spaced
    interval between \code{lower} and one.}
  \item{blty}{The border line type, for the border that surrounds
    the triangular image. By default \code{blty} is zero, so no
    border is plotted. Plotting a border leads to (by default) an
    increase in \code{grid} (and hence computation time), to ensure
    that the image fits within it.}
  \item{grid}{For plotting, the function is evaluated at \code{grid^2}
    points.}
  \item{lower}{The minimum value for which colours are plotted. By
    default \eqn{\code{lower} = 1/3} as this is the theoretical
    minimum of the dependence function of the trivariate extreme
    value distribution.}
  \item{ord}{A vector of length three, which should be a permutation
    of the set \eqn{\{1,2,3\}}{{1,2,3}}. The points \eqn{(1,0,0)},
    \eqn{(0,1,0)} and \eqn{(0,0,1)} (the vertices of the simplex)
    are depicted clockwise from the top in the order defined by
    \code{ord}.}
  \item{lab}{A character vector of length three, in which case the
    \code{i}th margin is labelled using the \code{i}th component,
    or \code{NULL}, in which case no labels are given. By default,
    \code{lab} is \code{as.character(1:3)}. The actual location of
    the margins, and hence the labels, is defined by \code{ord}.}
  \item{lcex}{A numerical value giving the amount by which the
    labels should be scaled relative to the default. Ignored
    if \code{lab} is \code{NULL}.}
}

\note{
  The rows of \code{data} that contain missing values are not used
  in the estimation of the dependence structure, but every non-missing
  value is used in estimating the generalized extreme value margins.

  The dependence function of the trivariate extreme value
  distribution is defined in \code{\link{atvpar}}.
  The function \code{\link{atvpar}} calculates and plots dependence
  functions of trivariate logistic and trivariate asymmetric
  logistic models.
}

\value{
  \code{atvnonpar} calculates or plots a non-parametric estimate of
  the dependence function of the trivariate extreme value distribution.
}

\seealso{\code{\link{atvpar}}, \code{\link{abvnonpar}},
  \code{\link{fgev}}}

\examples{
s3pts <- matrix(rexp(30), nrow = 10, ncol = 3)
s3pts <- s3pts/rowSums(s3pts)
sdat <- rmvevd(100, dep = 0.6, model = "log", d = 3)
atvnonpar(s3pts, sdat)

\dontrun{atvnonpar(data = sdat, plot = TRUE)}
\dontrun{atvnonpar(data = sdat, plot = TRUE, ord = c(2,3,1), lab = LETTERS[1:3])}
\dontrun{atvpar(dep = 0.6, model = "log", plot = TRUE)}
\dontrun{atvpar(dep = 0.6, model = "log", plot = TRUE, blty = 1)}
}
\keyword{nonparametric}









\eof
\name{atvpar}
\alias{atvpar}
\title{Parametric Dependence Functions of Trivariate Extreme
    Value Models}
\description{
  Calculate or plot the dependence function \eqn{A} for the trivariate
  logistic and trivariate asymmetric logistic models.
}
\usage{
atvpar(x = rep(1/3,3), dep, asy, model = c("log", "alog"), plot =
    FALSE, col = heat.colors(12), blty = 0, grid = if(blty) 150 else 50,
    lower = 1/3, ord = 1:3, lab = as.character(1:3), lcex = 1)
}
\arguments{
  \item{x}{A vector of length three or a matrix with three columns,
    in which case the dependence function is evaluated across
    the rows (ignored if plot is \code{TRUE}). The elements/rows
    of the vector/matrix should be positive and should sum to one,
    or else they should have a positive sum, in which case the rows
    are rescaled and a warning is given. \eqn{A(1/3,1/3,1/3)} is
    returned by default since it is often a useful summary of
    dependence.}
  \item{dep}{The dependence parameter(s). For the logistic model,
    should be a single value. For the asymmetric logistic model,
    should be a vector of length four, or a single value,
    in which case the value is used for each of the four
    parameters (see \code{\link{rmvevd}}).}
  \item{asy}{The asymmetry parameters for the asymmetric logistic
    model. Should be a list with seven vector elements; three
    of length one, three of length two and one of length three,
    containing the asymmetry parameters for each separate component
    (see \code{\link{rmvevd}} and \bold{Examples}).}
  \item{model}{The specified model; a character string. Must be
    either \code{"log"} (the default) or \code{"alog"} (or any unique
    partial match), for the logistic and asymmetric logistic models
    respectively. The definition of each model is given (for general
    dimensions) in \code{\link{rmvevd}}.}
  \item{plot}{Logical; if \code{TRUE} the function is plotted. The
    minimum (evaluated) value is returned invisibly.
    If \code{FALSE} (the default), the following arguments are
    ignored.}
  \item{col}{A list of colours (see \code{\link{image}}). The first
    colours in the list represent smaller values, and hence
    stronger dependence. Each colour represents an equally spaced
    interval between \code{lower} and one.}
  \item{blty}{The border line type, for the border that surrounds
    the triangular image. By default \code{blty} is zero, so no
    border is plotted. Plotting a border leads to (by default) an
    increase in \code{grid} (and hence computation time), to ensure
    that the image fits within it.}
  \item{grid}{For plotting, the function is evaluated at \code{grid^2}
    points.}
  \item{lower}{The minimum value for which colours are plotted. By
    defualt \eqn{\code{lower} = 1/3} as this is the theoretical
    minimum of the dependence function of the trivariate extreme
    value distribution.}
  \item{ord}{A vector of length three, which should be a permutation
    of the set \eqn{\{1,2,3\}}{{1,2,3}}. The points
    \eqn{(1,0,0)}, \eqn{(0,1,0)} and \eqn{(0,0,1)} (the vertices of
    the simplex) are depicted clockwise from the top in
    the order defined by \code{ord}.}
  \item{lab}{A character vector of length three, in which case the
    \code{i}th margin is labelled using the \code{i}th component,
    or \code{NULL}, in which case no labels are given. The actual
    location of the margins, and hence the labels, is defined by
    \code{ord}.}
  \item{lcex}{A numerical value giving the amount by which the
    labels should be scaled relative to the default. Ignored
    if \code{lab} is \code{NULL}.}
}
\details{
  Let \eqn{z = (z_1,z_2,z_3)}{z = (z1,z2,z3)} and
  \eqn{w = (w_1,w_2,w_3)}{w = (w1,w2,w3)}.
  Any trivariate extreme value distribution can be written as
  \deqn{G(z) = \exp\left\{- \left\{\sum\nolimits_{j=1}^{3} y_j
    \right\} A\left(\frac{y_1}{\sum\nolimits_{j=1}^{3} y_j},
    \frac{y_2}{\sum\nolimits_{j=1}^{3} y_j},
    \frac{y_3}{\sum\nolimits_{j=1}^{3} y_j}\right)\right\}}{
    G(z) = exp{-(y1+y2+y3) A[y1/(y1+y2+y3),y2/(y1+y2+y3),
      y3/(y1+y2+y3)]}}
  for some function \eqn{A} defined on the simplex
  \eqn{S_3 = \{w \in R^3_+ : \sum\nolimits_{j=1}^{3}
    w_j = 1\}}{S_3 = {w: w1 + w2 + w3 = 1}},
  where \deqn{y_i = \{1+s_i(z_i-a_i)/b_i\}^{-1/s_i}}{
    yi = {1+si(zi-ai)/bi}^(-1/si)}
  for \eqn{1+s_i(z_i-a_i)/b_i > 0}{1+si(zi-ai)/bi > 0} and
  \eqn{i = 1,2,3}, and where the (generalized extreme value)
  marginal parameters are given by
  \eqn{(a_i,b_i,s_i)}{(ai,bi,si)}, \eqn{b_i > 0}{bi > 0}.
  If \eqn{s_i = 0}{si = 0} then \eqn{y_i}{yi} is defined by
  continuity.

  \eqn{A} is called (by some authors) the dependence function.
  It follows that \eqn{A(1,0,0) = A(0,1,0) = A(0,0,1) = 1}, and
  that \eqn{A} is a convex function with
  \eqn{\max(w_1,w_2,w_3) \leq A(w)\leq 1}{
    max(w1,w2,w3) <= A(w) <= 1} for
  all \eqn{w} in \eqn{S_3}.
  The lower and upper limits of \eqn{A} are obtained under complete
  dependence and mutual independence respectively.
  \eqn{A} does not depend on the marginal parameters.
}

\value{
  \code{atvpar} calculates or plots the dependence function
  for the trivariate logistic and trivariate asymmetric logistic
  models, at specified parameter values.
}

\seealso{\code{\link{atvnonpar}}, \code{\link{abvpar}},
  \code{\link{rmvevd}}, \code{\link{image}}}

\examples{
atvpar(dep = 0.5, model = "log")
s3pts <- matrix(rexp(30), nrow = 10, ncol = 3)
s3pts <- s3pts/rowSums(s3pts)
atvpar(s3pts, dep = 0.5, model = "log")
\dontrun{atvpar(dep = 0.05, model = "log", plot = TRUE, blty = 1)}
atvpar(dep = 0.95, model = "log", plot = TRUE, lower = 0.94)

asy <- list(.4, .1, .6, c(.3,.2), c(.1,.1), c(.4,.1), c(.2,.3,.2))
atvpar(s3pts, dep = 0.15, asy = asy, model = "alog")
atvpar(dep = 0.15, asy = asy, model = "al", plot = TRUE, lower = 0.7)
}
\keyword{distribution}









\eof
\name{bvevd}
\alias{dbvevd}
\alias{pbvevd}
\alias{rbvevd}
\title{Parametric Bivariate Extreme Value Distributions}
\description{
  Density function, distribution function and random generation
  for eight parametric bivariate extreme value models. 
}
\synopsis{
dbvevd(x, dep, asy = c(1, 1), alpha, beta, model = c("log", "alog",
    "hr", "neglog", "aneglog", "bilog", "negbilog", "ct"),
    mar1 = c(0, 1, 0), mar2 = mar1, log = FALSE) 
pbvevd(q, dep, asy = c(1, 1), alpha, beta, model = c("log", "alog",
    "hr", "neglog", "aneglog", "bilog", "negbilog", "ct"),
    mar1 = c(0, 1, 0), mar2 = mar1, lower.tail = TRUE) 
rbvevd(n, dep, asy = c(1, 1), alpha, beta, model = c("log", "alog",
    "hr", "neglog", "aneglog", "bilog", "negbilog", "ct"),
    mar1 = c(0, 1, 0), mar2 = mar1) 
}
\usage{
dbvevd(x, dep, asy = c(1, 1), alpha, beta, model = "log",
    mar1 = c(0, 1, 0), mar2 = mar1, log = FALSE) 
pbvevd(q, dep, asy = c(1, 1), alpha, beta, model = "log",
    mar1 = c(0, 1, 0), mar2 = mar1, lower.tail = TRUE) 
rbvevd(n, dep, asy = c(1, 1), alpha, beta, model = "log",
    mar1 = c(0, 1, 0), mar2 = mar1) 
}
\arguments{
  \item{x, q}{A vector of length two or a matrix with two columns,
    in which case the density/distribution is evaluated across
    the rows.}
  \item{n}{Number of observations.}
  \item{dep}{Dependence parameter for the logistic, asymmetric
    logistic, Husler-Reiss, negative logistic and asymmetric
    negative logistic models.}
  \item{asy}{A vector of length two, containing the two asymmetry
    parameters for the asymmetric logistic and asymmetric negative
    logistic models.}
  \item{alpha, beta}{Alpha and beta parameters for the bilogistic,
    negative bilogistic and Coles-Tawn models.}
  \item{model}{The specified model; a character string. Must be
    either \code{"log"} (the default), \code{"alog"}, \code{"hr"},
    \code{"neglog"}, \code{"aneglog"}, \code{"bilog"},
    \code{"negbilog"} or \code{"ct"} (or any unique partial match),
    for the logistic, asymmetric logistic, Husler-Reiss, negative
    logistic, asymmetric negative logistic, bilogistic, negative
    bilogistic and Coles-Tawn models respectively. If parameter
    arguments are given that do not correspond to the specified
    model those arguments are ignored, with a warning.}
  \item{mar1, mar2}{Vectors of length three containing marginal
    parameters, or matrices with three columns where each
    column represents a vector of values to be passed to the
    corresponding marginal parameter.}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x].}
}
\details{
  Define
  \deqn{y_i = y_i(z_i) = \{1+s_i(z_i-a_i)/b_i\}^{-1/s_i}}{
    yi = yi(zi) = {1+si(zi-ai)/bi}^(-1/si)}
  for \eqn{1+s_i(z_i-a_i)/b_i > 0}{1+si(zi-ai)/bi > 0} and
  \eqn{i = 1,2}, where the marginal parameters are given by
  \eqn{\code{mari} = (a_i,b_i,s_i)}{\code{mari} = (ai,bi,si)},
  \eqn{b_i > 0}{bi > 0}.
  If \eqn{s_i = 0}{si = 0} then \eqn{y_i}{yi} is defined by
  continuity.
  
  In each of the bivariate distributions functions
  \eqn{G(z_1,z_2)}{G(z1,z2)} given below, the univariate margins
  are generalized extreme value, so that
  \eqn{G(z_i) = \exp(-y_i)}{G(zi) = \exp(-yi)} for \eqn{i = 1,2}.
  If \eqn{1+s_i(z_i-a_i)/b_i \leq 0}{1+si(zi-ai)/bi <= 0} for some
  \eqn{i = 1,2}, the value \eqn{z_i}{zi} is either greater than the
  upper end point (if \eqn{s_i < 0}{si < 0}), or less than the lower
  end point (if \eqn{s_i > 0}{si > 0}), of the \eqn{i}th univariate
  marginal distribution.
  
  \code{model = "log"} (Gumbel, 1960)
  
  The bivariate logistic distribution function with
  parameter \eqn{\code{dep} = r} is
  \deqn{G(z_1,z_2) = \exp\left[-(y_1^{1/r}+y_2^{1/r})^r\right]}{
    G(z1,z2) = exp{-[y1^(1/r)+y2^(1/r)]^r}}
  where \eqn{0 < r \leq 1}{0 < r <= 1}.
  This is a special case of the bivariate asymmetric logistic
  model.
  Complete dependence is obtained in the limit as
  \eqn{r} approaches zero.
  Independence is obtained when \eqn{r = 1}.

  
  \code{model = "alog"} (Tawn, 1988)
  
  The bivariate asymmetric logistic distribution function with
  parameters \eqn{\code{dep} = r} and
  \eqn{\code{asy} = (t_1,t_2)}{\code{asy} = (t1,t2)} is
  \deqn{G(z_1,z_2) = \exp\left\{-(1-t_1)y_1-(1-t_2)y_2-
    [(t_1y_1)^{1/r}+(t_2y_2)^{1/r}]^r\right\}}{
    G(z1,z2) = exp{-(1-t1)y1-(1-t2)y2-[(t1y1)^(1/r)+(t2y2)^(1/r)]^r}}
  where \eqn{0 < r \leq 1}{0 < r <= 1} and
  \eqn{0 \leq t_1,t_2 \leq 1}{0 <= t1,t2 <= 1}.
  When \eqn{t_1 = t_2 = 1}{t1 = t2 = 1} the asymmetric logistic
  model is equivalent to the logistic model.
  Independence is obtained when either \eqn{r = 1},
  \eqn{t_1 = 0}{t1 = 0} or \eqn{t_2 = 0}{t2 = 0}.
  Complete dependence is obtained in the limit when
  \eqn{t_1 = t_2 = 1}{t1 = t2 = 1} and \eqn{r}
  approaches zero.
  Different limits occur when \eqn{t_1}{t1} and \eqn{t_2}{t2}
  are fixed and \eqn{r} approaches zero.


  \code{model = "hr"} (Husler and Reiss, 1989)
  
  The Husler-Reiss distribution function with parameter
  \eqn{\code{dep} = r} is
  \deqn{G(z_1,z_2) = \exp\left(-y_1\Phi\{r^{-1}+{\textstyle\frac{1}{2}}
    r[\log(y_1/y_2)]\} - y_2\Phi\{r^{-1}+{\textstyle\frac{1}{2}}r
    [\log(y_2/y_1)]\}\right)}{
    G(z1,z2) = exp(-y1 Phi{r^{-1}+r[log(y1/y2)]/2} -
    y2 Phi{r^{-1}+r[log(y2/y1)]/2}}
  where \eqn{\Phi(\cdot)}{Phi()} is the standard normal distribution
  function and \eqn{r > 0}.
  Independence is obtained in the limit as \eqn{r} approaches zero.
  Complete dependence is obtained as \eqn{r} tends to infinity.


  \code{model = "neglog"} (Galambos, 1975)

  The bivariate negative logistic distribution function
  with parameter \eqn{\code{dep} = r} is
  \deqn{G(z_1,z_2) = \exp\left\{-y_1-y_2+
    [y_1^{-r}+y_2^{-r}]^{-1/r}\right\}}{
    G(z1,z2) = exp{-y1-y2+[y1^(-r)+y2^(-r)]^(-1/r)}}
  where \eqn{r > 0}.
  This is a special case of the bivariate asymmetric negative
  logistic model.
  Independence is obtained in the limit as \eqn{r} approaches zero.
  Complete dependence is obtained as \eqn{r} tends to infinity.
  The earliest reference to this model appears to be
  Galambos (1975, Section 4).


  \code{model = "aneglog"} (Joe, 1990)
  
  The bivariate asymmetric negative logistic distribution function
  with parameters parameters \eqn{\code{dep} = r} and
  \eqn{\code{asy} = (t_1,t_2)}{\code{asy} = (t1,t2)} is
  \deqn{G(z_1,z_2) = \exp\left\{-y_1-y_2+
    [(t_1y_1)^{-r}+(t_2y_2)^{-r}]^{-1/r}\right\}}{
    G(z1,z2) = exp{-y1-y2+[(t1y1)^(-r)+(t2y2)^(-r)]^(-1/r)}}
  where \eqn{r > 0} and \eqn{0 < t_1,t_2 \leq 1}{0 < t1,t2 <= 1}.
  When \eqn{t_1 = t_2 = 1}{t1 = t2 = 1} the asymmetric negative
  logistic model is equivalent to the negative logistic model.
  Independence is obtained in the limit as either \eqn{r},
  \eqn{t_1}{t1} or \eqn{t_2}{t2} approaches zero.
  Complete dependence is obtained in the limit when
  \eqn{t_1 = t_2 = 1}{t1 = t2 = 1} and \eqn{r}
  tends to infinity.
  Different limits occur when \eqn{t_1}{t1} and \eqn{t_2}{t2}
  are fixed and \eqn{r} tends to infinity.
  The earliest reference to this model appears to be Joe (1990),
  who introduces a multivariate extreme value distribution which
  reduces to \eqn{G(z_1,z_2)}{G(z1,z2)} in the bivariate case.


  \code{model = "bilog"} (Smith, 1990)
  
  The bilogistic distribution function with
  parameters \eqn{\code{alpha} = \alpha}{\code{alpha} = alpha}
  and \eqn{\code{beta} = \beta}{\code{beta} = beta} is
  \deqn{G(z_1,z_2) = \exp\left\{-y_1 q^{1-\alpha} -
    y_2 (1-q)^{1-\beta}\right\}}{
    G(z1,z2) = exp{- y1 q^(1-alpha) - y2 (1-q)^(1-beta)}}
  where
  \eqn{q = q(y_1,y_2;\alpha,\beta)}{q = q(y1,y2;alpha,beta)}
  is the root of the equation
  \deqn{(1-\alpha) y_1 (1-q)^\beta - (1-\beta) y_2 q^\alpha = 0,}{
    (1-alpha) y1 (1-q)^beta - (1-beta) y2 q^alpha = 0,}
  \eqn{0 < \alpha,\beta < 1}{0 < alpha,beta < 1}.
  When \eqn{\alpha = \beta}{alpha = beta} the bilogistic model
  is equivalent to the logistic model with dependence parameter
  \eqn{\code{dep} = \alpha = \beta}{\code{dep} = alpha = beta}.
  Complete dependence is obtained in the limit as
  \eqn{\alpha = \beta}{alpha = beta} approaches zero.
  Independence is obtained as
  \eqn{\alpha = \beta}{alpha = beta} approaches one, and when
  one of \eqn{\alpha,\beta}{alpha,beta} is fixed and the other
  approaches one.
  Different limits occur when one of
  \eqn{\alpha,\beta}{alpha,beta} is fixed and the other
  approaches zero.
  A bilogistic model is fitted in Smith (1990), where it appears
  to have been first introduced.


  \code{model = "negbilog"} (Coles and Tawn, 1994)

  The negative bilogistic distribution function with
  parameters \eqn{\code{alpha} = \alpha}{\code{alpha} = alpha}
  and \eqn{\code{beta} = \beta}{\code{beta} = beta} is
  \deqn{G(z_1,z_2) = \exp\left\{- y_1 - y_2 + y_1 q^{1+\alpha} +
    y_2 (1-q)^{1+\beta}\right\}}{
    G(z1,z2) = exp{- y1 - y2 + y1 q^(1+alpha) + y2 (1-q)^(1+beta)}}
  where
  \eqn{q = q(y_1,y_2;\alpha,\beta)}{q = q(y1,y2;alpha,beta)}
  is the root of the equation
  \deqn{(1+\alpha) y_1 q^\alpha - (1+\beta) y_2 (1-q)^\beta = 0,}{
    (1+alpha) y1 q^alpha - (1+beta) y2 (1-q)^beta = 0,}
  \eqn{\alpha > 0}{alpha > 0} and \eqn{\beta > 0}{beta > 0}.
  When \eqn{\alpha = \beta}{alpha = beta} the negative bilogistic
  model is equivalent to the negative logistic model with dependence
  parameter
  \eqn{\code{dep} = 1/\alpha = 1/\beta}{
    \code{dep} = 1/alpha = 1/beta}.
  Complete dependence is obtained in the limit as
  \eqn{\alpha = \beta}{alpha = beta} approaches zero.
  Independence is obtained as
  \eqn{\alpha = \beta}{alpha = beta} tends to infinity, and when
  one of \eqn{\alpha,\beta}{alpha,beta} is fixed and the other
  tends to infinity.
  Different limits occur when one of
  \eqn{\alpha,\beta}{alpha,beta} is fixed and the other
  approaches zero.


  \code{model = "ct"} (Coles and Tawn, 1991)
  
  The Coles-Tawn distribution function with
  parameters \eqn{\code{alpha} = \alpha > 0}{\code{alpha} = alpha > 0}
  and \eqn{\code{beta} = \beta > 0}{\code{beta} = beta > 0} is
  \deqn{G(z_1,z_2) =
    \exp\left\{-y_1 [1 - \mbox{Be}(q;\alpha+1,\beta)] -
    y_2 \mbox{Be}(q;\alpha,\beta+1) \right\}}{
    G(z1,z2) = exp{- y1 [1 - Be(q;alpha+1,beta)] - y2 Be(q;alpha,beta+1)}}
  where
  \eqn{q = \alpha y_2 / (\alpha y_2 + \beta y_1)}{
    q = alpha y2 / (alpha y2 + beta y1)} and
  \eqn{\mbox{Be}(q;\alpha,\beta)}{Be(q;alpha,beta)} is the beta
  distribution function evaluated at \eqn{q} with
  \eqn{\code{shape1} = \alpha}{\code{shape1} = alpha} and
  \eqn{\code{shape2} = \beta}{\code{shape2} = beta}.
  Complete dependence is obtained in the limit as
  \eqn{\alpha = \beta}{alpha = beta} tends to infinity.
  Independence is obtained as
  \eqn{\alpha = \beta}{alpha = beta} approaches zero, and when
  one of \eqn{\alpha,\beta}{alpha,beta} is fixed and the other
  approaches zero.
  Different limits occur when one of
  \eqn{\alpha,\beta}{alpha,beta} is fixed and the other
  tends to infinity.
}
\value{
  \code{dbvevd} gives the density function, \code{pbvevd} gives the
  distribution function and \code{rbvevd} generates random deviates,
  for one of eight parametric bivariate extreme value models.
}
\note{
  The logistic and asymmetric logistic models respectively are
  simulated using bivariate versions of Algorithms 1.1 and 1.2 in
  Stephenson(2003).
  All other models are simulated using a root finding algorithm
  to simulate from the conditional distributions.

  The simulation of the bilogistic and negative bilogistic models
  requires a root finding algorithm to evaluate \eqn{q}
  within the root finding algorithm used to simulate from the
  conditional distributions.
  The generation of bilogistic and negative bilogistic random
  deviates is therefore relatively slow (about 2.8 seconds per
  1000 random vectors on a 450MHz PIII, 512Mb RAM).

  The bilogistic and negative bilogistic models can be represented
  under a single model, using the integral of the maximum of two
  beta distributions (Joe, 1997).

  The Coles-Tawn model is called the Dirichelet model in Coles
  and Tawn (1991).
}
\references{
  Coles, S. G. and Tawn, J. A. (1991)
  Modelling extreme multivariate events.
  \emph{J. Roy. Statist. Soc., B}, \bold{53}, 377--392.
  
  Coles, S. G. and Tawn, J. A. (1994)
  Statistical methods for multivariate extremes: an application to
  structural design (with discussion).
  \emph{Appl. Statist.}, \bold{43}, 1--48.
  
  Galambos, J. (1975)
  Order statistics of samples from multivariate distributions.
  \emph{J. Amer. Statist. Assoc.}, \bold{70}, 674--680.
  
  Gumbel, E. J. (1960)
  Distributions des valeurs extremes en plusieurs dimensions.
  \emph{Publ. Inst. Statist. Univ. Paris}, \bold{9}, 171--173.

  Husler, J. and Reiss, R.-D. (1989)
  Maxima of normal random vectors: between independence 
  and complete dependence.
  \emph{Statist. Probab. Letters}, \bold{7}, 283--286.

  Joe, H. (1990)
  Families of min-stable multivariate exponential and multivariate
  extreme value distributions.
  \emph{Statist. Probab. Letters}, \bold{9}, 75--81.

  Joe, H. (1997)
  \emph{Multivariate Models and Dependence Concepts},
  London: Chapman & Hall.

  Smith, R. L. (1990)
  Extreme value theory. In
  \emph{Handbook of Applicable Mathematics} (ed. W. Ledermann),
  vol. 7. Chichester: John Wiley, pp. 437--471.
  
  Stephenson, A. G. (2003)
  Simulating multivariate extreme value distributions of logistic type.
  \emph{Extremes}, \bold{6}(1), 49--60.

  Tawn, J. A. (1988)
  Bivariate extreme value theory: models and estimation.
  \emph{Biometrika}, \bold{75}, 397--415.  
}

\seealso{\code{\link{abvpar}}, \code{\link{rgev}}, \code{\link{rmvevd}}}


\examples{
pbvevd(matrix(rep(0:4,2), ncol=2), dep = 0.7, model = "log")
pbvevd(c(2,2), dep = 0.7, asy = c(0.6,0.8), model = "alog")
pbvevd(c(1,1), dep = 1.7, model = "hr")

margins <- cbind(0, 1, seq(-0.5,0.5,0.1))
rbvevd(11, dep = 1.7, model = "hr", mar1 = margins)
rbvevd(10, dep = 1.2, model = "neglog", mar1 = c(10, 1, 1))
rbvevd(10, alpha = 0.7, beta = 0.52, model = "bilog")

dbvevd(c(0,0), dep = 1.2, asy = c(0.5,0.9), model = "aneglog")
dbvevd(c(0,0), alpha = 0.75, beta = 0.5, model = "ct", log = TRUE)
dbvevd(c(0,0), alpha = 0.7, beta = 1.52, model = "negbilog")
}
\keyword{distribution}





\eof
\name{chiplot}
\alias{chiplot}
\title{Dependence Measure Plots}
\description{
  Plots of estimates of the dependence measures chi and chi-bar
  for bivariate data. 
}
\usage{
chiplot(data, nq = 100, qlim = NULL, which = 1:2, conf = 0.95, lty = 1,
    cilty = 2, col = 1, cicol = 1, xlim = c(0,1), ylim1 = NULL,
    ylim2 = c(-1,1), main1 = "Chi Plot", main2 = "Chi Bar Plot", xlab =
    "Quantile", ylab1 = "Chi", ylab2 = "Chi Bar", ask = nb.fig <
    length(which) && dev.interactive(), \dots)
}
\arguments{
  \item{data}{A matrix or data frame with two columns. Rows
    (observations) with missing values are stripped from
    the data before any computations are performed.}
  \item{nq}{The number of quantiles at which the measures
    are evaluated.}
  \item{qlim}{The limits of the quantiles at which the measures
    are evaluated (see \bold{Details}).}
  \item{which}{If only one plot is required, specify \code{1}
    for chi and \code{2} for chi-bar.}
  \item{conf}{The confidence coefficient of the plotted confidence
    intervals.}
  \item{lty, cilty}{Line types for the estimates of the measures and
    for the confidence intervals respectively. Use zero to supress.}
  \item{col, cicol}{Colour types for the estimates of the measures
    and for the confidence intervals respectively.}
  \item{xlim, xlab}{Limits and labels for the x-axis; they apply
    to both plots.}
  \item{ylim1}{Limits for the y-axis of the chi plot. If this
    is \code{NULL} (the default) the upper limit is one, and the
    lower limit is the minimum of zero and the smallest plotted
    value.}
  \item{ylim2}{Limits for the y-axis of the chi-bar plot.}
  \item{main1, main2}{The plot titles for the chi and chi-bar plots
    respectively.}
  \item{ylab1, ylab2}{The y-axis labels for the chi and chi-bar plots
    respectively.} 
  \item{ask}{Logical; if \code{TRUE}, the user is asked before
    each plot.}
  \item{\dots}{Other arguments to be passed to \code{matplot}.}
}

\details{
  These measures are explained in full detail in Coles, Heffernan
  and Tawn (1999). A brief treatment is also given in Section
  8.4 of Coles(2001).
  A short summary is given as follows.
  We assume that the data are \emph{iid} random vectors with common
  bivariate distribution function \eqn{G}, and we define the random
  vector \eqn{(X,Y)} to be distributed according to \eqn{G}.
  
  The chi plot is a plot of \eqn{q} against empirical estimates of
  \deqn{\chi(q) = 2 - \log(\Pr(F_X(X) < q, F_Y(Y) < q)) / \log(q)}{
    chi(q) = 2 - log(Pr(F_X(X) < q, F_Y(Y) < q)) / log(q)}
  where \eqn{F_X} and \eqn{F_Y} are the marginal distribution
  functions, and where \eqn{q} is in the interval (0,1).
  The quantity \eqn{\chi(q)}{chi(q)} is bounded by
  \deqn{2 - \log(2u - 1)/\log(u) \leq \chi(q) \leq 1}{
    2 - log(2u - 1)/log(u) <= chi(q) <= 1}
  where the lower bound is interpreted as \code{-Inf} for
  \eqn{q \leq 1/2}{q <= 1/2} and zero for \eqn{q = 1}.
  These bounds are reflected in the corresponding estimates.

  The chi bar plot is a plot of \eqn{q} against empirical estimates of
  \deqn{\bar{\chi}(q) = 2\log(1-q)/\log(\Pr(F_X(X) > q, F_Y(Y) > q)) - 1}{
    chibar(q) = 2log(1-q)/log(Pr(F_X(X) > q, F_Y(Y) > q)) - 1}
  where \eqn{F_X} and \eqn{F_Y} are the marginal distribution
  functions, and where \eqn{q} is in the interval (0,1).
  The quantity \eqn{\bar{\chi}(q)}{chibar(q)} is bounded by
  \eqn{-1 \leq \bar{\chi}(q) \leq 1}{-1 <= chibar(q) <= 1}
  and these bounds are reflected in the corresponding estimates.

  Note that the empirical estimators for \eqn{\chi(q)}{chi(q)} and
  \eqn{\bar{\chi}(q)}{chibar(q)} are undefined near \eqn{q=0} and \eqn{q=1}. By
  default the function takes the limits of \eqn{q} so that the plots
  depicts all values at which the estimators are defined. This can be
  overridden by the argument \code{qlim}, which must represent a subset
  of the default values (and these can be determined using the
  component \code{quantile} of the invisibly returned list; see
  \bold{Value}).
  
  The confidence intervals within the plot assume that observations are
  independent, and that the marginal distributions are estimated exactly.
  The intervals are constructed using the delta method; this may
  lead to poor interval estimates near \eqn{q=0} and \eqn{q=1}.
  
  The function \eqn{\chi(q)}{chi(q)} can be interpreted as a quantile
  dependent measure of dependence. In particular, the sign of
  \eqn{\chi(q)}{chi(q)} determines whether the variables are positively
  or negatively associated at quantile level \eqn{q}.
  By definition, variables are said to be asymptotically independent
  when \eqn{\chi(1)}{chi(1)} (defined in the limit) is zero.
  For independent variables, \eqn{\chi(q) = 0}{chi(q) = 0} for all
  \eqn{q} in (0,1).
  For perfectly dependent variables, \eqn{\chi(q) = 1}{chi(q) = 1}
  for all \eqn{q} in (0,1).
  For bivariate extreme value distributions, \eqn{\chi(q) =
    2(1-A(1/2))}{chi(q) = 2(1-A(1/2))}
  for all \eqn{q} in (0,1), where \eqn{A} is the dependence function,
  as defined in \code{\link{abvpar}}. If a bivariate threshold model
  is to be fitted (using \code{\link{fbvpot}}), this plot can therefore
  act as a threshold identification plot, since e.g. the use of 95\%
  marginal quantiles as threshold values implies that \eqn{\chi(q)}{chi(q)}
  should be approximately constant above \eqn{q = 0.95}.

  The function \eqn{\bar{\chi}(q)}{chibar(q)} can again be interpreted
  as a quantile dependent measure of dependence; it is most useful
  within the class of asymptotically independent variables.
  For asymptotically dependent variables (i.e. those for which
  \eqn{\chi(1) < 1}{chi(1) < 1}), we have \eqn{\bar{\chi}(1) = 1}{
    chibar(1) = 1}, where
  \eqn{\bar{\chi}(1)}{chibar(1)} is again defined in the limit.
  For asymptotically independent variables, \eqn{\bar{\chi}(1)}{
    chibar(1)} provides a measure that increases with dependence strength.
  For independent variables \eqn{\bar{\chi}(q) = 0}{chibar(q) = 0} for
  all \eqn{q} in (0,1), and hence \eqn{\bar{\chi}(1) = 0}{chibar(1) = 0}.
}
\value{
  A list with components \code{quantile}, \code{chi} (if \code{1} is in
  \code{which}) and \code{chibar} (if \code{2} is in \code{which})
  is invisibly returned.
  The components \code{quantile} and \code{chi} contain those objects
  that were passed to the formal arguments \code{x} and \code{y} of
  \code{matplot} in order to create the chi plot.
  The components \code{quantile} and \code{chibar} contain those objects
  that were passed to the formal arguments \code{x} and \code{y} of
  \code{matplot} in order to create the chi-bar plot.
}
\references{
  Coles, S. G., Heffernan, J. and Tawn, J. A. (1999)
  Dependence measures for extreme value analyses.
  \emph{Extremes}, \bold{2}, 339--365.

  Coles, S. G. (2001)
  \emph{An Introduction to Statistical Modelling of Extreme Values},
  London: Springer--Verlag.  
}
\author{Jan Heffernan and Alec Stephenson}  

\seealso{\code{\link{fbvevd}}, \code{\link{fbvpot}},
  \code{\link{matplot}}}

\examples{
par(mfrow = c(1,2))
smdat1 <- rbvevd(1000, dep = 0.6, model = "log")
smdat2 <- rbvevd(1000, dep = 1, model = "log")
chiplot(smdat1)
chiplot(smdat2)
}
\keyword{hplot}


\eof
\name{clusters}
\alias{clusters}
\title{Identify Clusters of Exceedences}
\description{
  Identify clusters of exceedences. 
}
\usage{
clusters(data, u, r = 1, ulow = -Inf, rlow = 1, cmax = FALSE, keep.names
    = TRUE, plot = FALSE, xdata = seq(along = data), lvals = TRUE, lty =
    1, lwd = 1, pch = par("pch"), col = if(n > 250) NULL else "grey",
    xlab = "Index", ylab = "Data", ...)
}
\arguments{
  \item{data}{A numeric vector, which may contain missing values.}
  \item{u}{A single value giving the threshold, unless a time varying
    threshold is used, in which case \code{u} should be a vector of
    thresholds, typically with the same length as \code{data} (or else
    the usual recycling rules are applied).}
  \item{r}{A postive integer denoting the clustering interval
    length. By default the interval length is one.}
  \item{ulow}{A single value giving the lower threshold, unless a time
    varying lower threshold is used, in which case \code{ulow} should
    be a vector of lower thresholds, typically with the same length as
    \code{data} (or else the usual recycling rules are applied).
    By default there is no lower threshold (or equivalently, the
    lower threshold is \code{-Inf}).}
  \item{rlow}{A postive integer denoting the lower clustering interval
    length. The lower clustering interval length is only relevant if
    it is less than the clustering interval length \code{r} and if
    there exists a lower threshold (greater than \code{-Inf}).}
  \item{cmax}{Logical; if \code{FALSE} (the default), a list
    containing the clusters of exceedences is returned. If
    \code{TRUE} a numeric vector containing the cluster maxima
    is returned.}
  \item{keep.names}{Logical; if \code{FALSE}, the function makes
    no attempt to retain the names/indices of the observations
    within the returned object. If \code{data} contains a large
    number of observations, this can make the function run much
    faster. The argument is mainly designed for internal use.}
  \item{plot}{Logical; if \code{TRUE} a plot is given that depicts
    the identified clusters, and the clusters (if \code{cmax} is
    \code{FALSE}) or cluster maxima (if \code{cmax} is \code{TRUE})
    are returned invisibly. If \code{FALSE} (the default), the
    following arguments are ignored.}
  \item{xdata}{A numeric vector with the same length as \code{data},
    giving the values to be plotted on the x-axis.}
  \item{lvals}{Logical; should the values below the threshold and
    the line depicting the lower threshold be plotted?}
  \item{lty, lwd}{Line type and width for the lines depicting the
    threshold and the lower threshold.}
  \item{pch}{Plotting character.}
  \item{col}{Strips of colour \code{col} are used to identify the
    clusters. An observation is contained in the cluster if the 
    centre of the corresponding plotting character is contained
    in the coloured strip. If \code{col} is \code{NULL} the strips
    are omitted. By default the strips are coloured \code{"grey"},
    but are omitted whenever \code{data} contains more than 250
    observations.}
  \item{xlab, ylab}{Labels for the x and y axis.}
  \item{\dots}{Other graphics parameters.}
}
\details{
  The clusters of exceedences are identified as follows.
  The first exceedence of the threshold initiates the first cluster.
  The first cluster then remains active until either \code{r}
  consecutive values fall below (or are equal to) the threshold,
  or until \code{rlow} consecutive values fall below (or are equal
  to) the lower threshold.
  The next exceedence of the threshold (if it exists) then initiates
  the second cluster, and so on.
  Missing values are allowed, in which case they are treated as
  falling below (or equal to) the threshold, but falling above the
  lower threshold.
}
\value{
  If \code{cmax} is \code{FALSE} (the default), a list with one
  component for each identified cluster.
  If \code{cmax} is \code{TRUE}, a numeric vector containing the
  cluster maxima.
  In any case, the returned object has an attribute \code{acs},
  giving the average cluster size (where the cluster size is
  defined as the number of exceedences within a cluster), which
  will be \code{NaN} if there are no values above the threshold
  (and hence no clusters).
  
  If \code{plot} is \code{TRUE}, the list of clusters, or vector
  of cluster maxima, is returned invisibly.
}

\seealso{\code{\link{exi}}}

\examples{
data(portpirie)
clusters(portpirie, 4.2, 3)
clusters(portpirie, 4.2, 3, cmax = TRUE)
clusters(portpirie, 4.2, 3, 3.8, plot = TRUE)
clusters(portpirie, 4.2, 3, 3.8, plot = TRUE, lvals = FALSE)
tvu <- c(rep(4.2, 20), rep(4.1, 25), rep(4.2, 20))
clusters(portpirie, tvu, 3, plot = TRUE)
}
\keyword{manip}


\eof
\name{evd-internal}

\alias{fbvlog}
\alias{fbvalog}
\alias{fbvhr}
\alias{fbvneglog}
\alias{fbvaneglog}
\alias{fbvbilog}
\alias{fbvnegbilog}
\alias{fbvct}

\alias{pcint}
\alias{fgev.quantile}
\alias{fgev.norm}
\alias{fpot.quantile}
\alias{fpot.norm}
\alias{bvdepfn}
\alias{tvdepfn}
\alias{mvalog.check}

\alias{ccop}
\alias{mtransform}
\alias{na.vals}
\alias{bvpost.optim}
\alias{bvstart.vals}
\alias{nsloc.transform}
\alias{sep.bvdata}

\alias{dens}
\alias{pp}
\alias{qq}
\alias{rl}
\alias{dens.gev}
\alias{pp.gev}
\alias{qq.gev}
\alias{rl.gev}
\alias{dens.pot}
\alias{pp.pot}
\alias{qq.pot}
\alias{rl.pot}

\alias{bvdens}
\alias{bvcpp}
\alias{bvdp}

\alias{abvalog}
\alias{abvaneglog}
\alias{abvhr}
\alias{abvlog}
\alias{abvneglog}
\alias{abvbilog}
\alias{abvnegbilog}
\alias{abvct}

\alias{rbvlog}
\alias{rbvalog}
\alias{rbvhr}
\alias{rbvneglog}
\alias{rbvaneglog}
\alias{rbvbilog}
\alias{rbvnegbilog}
\alias{rbvct}

\alias{pbvlog}
\alias{pbvalog}
\alias{pbvhr}
\alias{pbvneglog}
\alias{pbvaneglog}
\alias{pbvbilog}
\alias{pbvnegbilog}
\alias{pbvct}

\alias{dbvlog}
\alias{dbvalog}
\alias{dbvhr}
\alias{dbvneglog}
\alias{dbvaneglog}
\alias{dbvbilog}
\alias{dbvnegbilog}
\alias{dbvct}

\alias{rmvlog}
\alias{rmvalog}
\alias{pmvlog}
\alias{pmvalog}
\alias{dmvlog}
\alias{dmvalog}
\alias{atvalog}
\alias{atvlog}

\alias{fbvcpot}
\alias{fbvppot}
\alias{bvtpost.optim}

\alias{fbvclog}
\alias{fbvcalog}
\alias{fbvcaneglog}
\alias{fbvcbilog}
\alias{fbvcct}
\alias{fbvcnegbilog}
\alias{fbvcneglog}
\alias{fbvchr}

\alias{fbvplog}
\alias{fbvpbilog}
\alias{fbvpct}
\alias{fbvpnegbilog}
\alias{fbvpneglog}

\title{Internal Functions}

\description{
  The evd package contains many internal functions that are not
  designed to be called by the user.

  The generic functions \code{dens}, \code{pp}, \code{qq} and
  \code{rl} create the diagnostic plots generated by
  \code{plot.uvevd}.
  Similarly, \code{bvdens}, \code{bvcpp} and \code{bvdp} create
  the diagnostic plots generated by \code{plot.bvevd}.

  There are internal fitting, simulation, distribution and density
  functions for each bivariate and multivariate parametric model,
  which are called from functions such as \code{rbvevd} and
  \code{rmvevd}. There also exists internal functions for the
  calculation and plotting of dependence functions of 
  bivariate and trivariate models, which are called from
  \code{abvdep} and \code{atvdep}.
  The dependence functions are ultimately plotted by the low-level
  functions \code{bvdepfn} and \code{tvdepfn}.
  
  The function \code{pcint} calculates profile confidence intervals,
  and is called from the function \code{plot.profile.evd}. The fitting
  function \code{fgev} calls the internal functions \code{fgev.quantile}
  and \code{fgev.norm} for fits under different parameterizations.
  The fitting function \code{fpot} calls the internal functions
  \code{fpot.norm} and \code{fpot.quantile}.
  Marginal transformations are executed using \code{mtransform}.
  
  The function \code{ccop} calculates condition copulas (i.e.
  conditional distributions under uniform margins) for each
  bivariate parametric model. This is needed to create the
  conditional P-P plots generated by \code{bvcpp}.

  The functions \code{nsloc.transform}, \code{na.vals},
  \code{bvpost.optim}, \code{bvstart.vals} and \code{sep.bvdata}
  are used in the fitting of bivariate models. The function
  \code{mvalog.check} checks and transforms the \code{asy} argument
  for the multivariate asymmetric model.

  For fitting bivariate threshold models, internal functions exist
  for the censored and (currently unimplemented) point process
  likelihoods, and each of these calls a further internal function
  corresponding to the specified model. The internal function
  \code{bvtpost.optim} is then used for post optimization processing.
}

\keyword{internal}

\eof
\name{evmc}
\alias{evmc}
\title{Simulate Markov Chains With Extreme Value
  Dependence Structures}
\description{
  Simulation of first order Markov chains, such that each pair
  of consecutive values has the dependence structure of one of
  eight parametric bivariate extreme value distributions. 
}
\synopsis{
evmc(n, dep, asy = c(1,1), alpha, beta, model = c("log", "alog",
    "hr", "neglog", "aneglog", "bilog", "negbilog", "ct"),
    margins = c("uniform","exponential","frechet","gumbel"))
}
\usage{
evmc(n, dep, asy = c(1,1), alpha, beta, model = "log",
    margins = "uniform")
}
\arguments{
  \item{n}{Number of observations.}
  \item{dep}{Dependence parameter for the logistic, asymmetric
    logistic, Husler-Reiss, negative logistic and asymmetric
    negative logistic models.}
  \item{asy}{A vector of length two, containing the two asymmetry
    parameters for the asymmetric logistic and asymmetric negative
    logistic models.}
  \item{alpha, beta}{Alpha and beta parameters for the bilogistic,
    negative bilogistic and Coles-Tawn models.}
  \item{model}{The specified model; a character string. Must be
    either \code{"log"} (the default), \code{"alog"}, \code{"hr"},
    \code{"neglog"}, \code{"aneglog"}, \code{"bilog"},
    \code{"negbilog"} or \code{"ct"} (or any unique partial match),
    for the logistic, asymmetric logistic, Husler-Reiss, negative
    logistic, asymmetric negative logistic, bilogistic, negative
    bilogistic and Coles-Tawn models respectively. The definition
    of each model is given in \code{\link{rbvevd}}. If parameter
    arguments are given that do not correspond to the specified
    model those arguments are ignored, with a warning.}
  \item{margins}{The marginal distribution of each value; a
    character string. Must be either \code{"uniform"} (the
    default), \code{"exponential"}, \code{"frechet"} or
    \code{"gumbel"} (or any unique partial match), for the uniform,
    standard exponential, standard Gumbel and standard Frechet
    distributions respectively.}
}
\value{
  A numeric vector of length \code{n}.
}

\seealso{\code{\link{marma}}, \code{\link{rbvevd}}}

\examples{
evmc(100, alpha = 0.1, beta = 0.1, model = "bilog")
evmc(100, dep = 10, model = "hr", margins = "exp")
}
\keyword{distribution}





\eof
\name{exi}
\alias{exi}
\title{Estimates of the Extremal Index}
\description{
  Estimates of the extremal index. 
}
\usage{
exi(data, u, r = 1, ulow = rep(-Inf, ncol(u)), rlow =
    rep(1, length(r)), dimnames = list(NULL, NULL), drop = TRUE)
}
\arguments{
  \item{data}{A numeric vector, which may contain missing values.}
  \item{u}{A numeric vector of thresholds, unless time varying
    thresholds are used, in which case \code{u} should be a matrix,
    typically with \code{length(data)} rows (or else the usual
    recycling rules are applied).}
  \item{r}{A numeric vector of postive integers denoting the
    clustering interval lengths.}
  \item{ulow}{A numeric vector of lower thresholds with length
    \code{length(u)} (if \code{u} is a vector) or \code{ncol(u)}
    (if \code{u} is a matrix), unless time varying lower
    thresholds are used, in which case \code{ulow} should be a
    matrix with \code{length(u)} or \code{ncol(u)} columns, and
    typically with \code{length(data)} rows (or else the
    usual recycling rules are applied). By default there are no
    lower thresholds (or equivalently, the lower thresholds are
    \code{-Inf}).}
  \item{rlow}{A numeric vector of postive integers, of length
    \code{length(r)}, denoting the lower clustering interval
    lengths. The lower clustering interval length is only
    relevant if it is less than the clustering interval length
    and if there exists a lower threshold (greater than
    \code{-Inf}).}
  \item{dimnames}{The \code{dimnames} of the result.}
  \item{drop}{Logical; return a vector if either \code{u} or
    \code{r} has one element/column?}
}
\details{
  The extremal index is estimated using the inverse of the
  average cluster size, using the clusters of exceedences
  derived from \code{\link{clusters}}.
  If the threshold is larger than (or equal to)
  \code{max(data)}, so that no clusters are derived, then
  the corresponding estimate is \code{NaN}.
}
\value{
  A matrix with \code{length(u)} (if \code{u} is a vector) or
  \code{ncol(u)} (if \code{u} is a matrix) rows and
  \code{length(r)} columns, such that the \code{ij}th element
  is the estimate of the extremal index using the threshold
  \code{u[i]} or \code{u[,i]}, the clustering interval length
  \code{r[j]}, the lower threshold \code{ulow[i]} or
  \code{ulow[,i]}, and the lower clustering interval length
  \code{rlow[j]}.
  If \code{drop} is \code{TRUE} (the default), a numeric vector is
  returned in preference to a matrix with one row/column.
}

\seealso{\code{\link{clusters}}}

\examples{
data(portpirie)
exi(portpirie, 4.2, 3, 3.8)
us <- seq(3.9, 4.2, len = 10)
exi(portpirie, us, 3)
exi(portpirie, us, 1:3)
tvu <- c(rep(4.2, 20), rep(4.1, 25), rep(4.2, 20))
exi(portpirie, as.matrix(tvu), 1:3)

us <- seq(3.9, 4.2, len = 100)
eis <- exi(portpirie, us, 1:3)
matplot(us, eis, type = "l", xlab = "Threshold", ylab = "Ext. Index")
}
\keyword{manip}


\eof
\name{extreme}
\alias{dextreme}
\alias{pextreme}
\alias{qextreme}
\alias{rextreme}
\title{Distributions of Maxima and Minima}
\description{
  Density function, distribution function, quantile function and
  random generation for the maximum/minimum of a given number of
  independent variables from a specified distribution. 
}
\usage{
dextreme(x, densfun, distnfun, \dots, distn, mlen = 1, largest = TRUE,
    log = FALSE)
pextreme(q, distnfun, \dots, distn, mlen = 1, largest = TRUE,
    lower.tail = TRUE) 
qextreme(p, quantfun, \dots, distn, mlen = 1, largest = TRUE,
    lower.tail = TRUE) 
rextreme(n, quantfun, \dots, distn, mlen = 1, largest = TRUE) 
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{p}{Vector of probabilities.}
  \item{n}{Number of observations.}
  \item{densfun, distnfun, quantfun}{Density, distribution and
    quantile function of the specified distribution. The density
    function must have a \code{log} argument (a simple wrapper
    can always be constructed to achieve this).}
  \item{\dots}{Parameters of the specified distribution.}
  \item{distn}{A character string, optionally given as an
    alternative to \code{densfun}, \code{distnfun} and \code{quantfun}
    such that the density, distribution and quantile functions are
    formed upon the addition of the prefixes \code{d}, \code{p} and
    \code{q} respectively.}
  \item{mlen}{The number of independent variables.}
  \item{largest}{Logical; if \code{TRUE} (default) use maxima,
    otherwise minima.}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default) probabilities
    are P[X <= x], otherwise  P[X > x].}
}
\value{
  \code{dextreme} gives the density function, \code{pextreme} gives the
  distribution function and \code{qextreme} gives the quantile function
  of the maximum/minimum of \code{mlen} independent variables from
  a specified distibution. \code{rextreme} generates random deviates.
}

\seealso{\code{\link{rgev}}, \code{\link{rorder}}}

\examples{
dextreme(2:4, dnorm, pnorm, mean = 0.5, sd = 1.2, mlen = 5)
dextreme(2:4, distn = "norm", mean = 0.5, sd = 1.2, mlen = 5)
dextreme(2:4, distn = "exp", mlen = 2)
pextreme(2:4, distn = "exp", rate = 1.2, mlen = 2)
qextreme(seq(0.9, 0.6, -0.1), distn = "exp", rate = 1.2, mlen = 2)
rextreme(5, qgamma, shape = 1, mlen = 10)
p <- (1:9)/10
pexp(qextreme(p, distn = "exp", rate = 1.2, mlen = 1), rate = 1.2)
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
}
\keyword{distribution}



\eof
\name{failure}
\alias{failure}
\title{Failure Times}
\usage{data(failure)}
\description{
  Failure times.
}
\format{A vector containing 24 observations.}
\source{
  van Montfort, M. A. J. and Otten, A. (1978)
  On testing a shape parameter in the presence of a scale
  parameter.
  \emph{Math. Operations Forsch. Statist., Ser. Statistics},
  \bold{9}, 91--104.
}
  
\keyword{datasets}

\eof
\name{fbvevd}
\alias{fbvevd}
\alias{print.bvevd}

\title{Maximum-likelihood Fitting of Bivariate Extreme Value Distributions}
\description{
  Fit models for one of eight parametric bivariate extreme value
  distributions, including linear modelling of the marginal
  location parameters, and allowing any of the parameters to be held
  fixed if desired. 
}
\synopsis{
fbvevd(x, model = c("log", "alog", "hr", "neglog", "aneglog", "bilog",
    "negbilog", "ct"), start, \dots, sym = FALSE, nsloc1 = NULL,
    nsloc2 = NULL, cshape = cscale, cscale = cloc, cloc = FALSE,
    std.err = TRUE, dsm = TRUE, corr = FALSE, method = "BFGS",
    warn.inf = TRUE)
}
\usage{
fbvevd(x, model = "log", start, \dots, sym = FALSE, nsloc1 = NULL,
    nsloc2 = NULL, cshape = cscale, cscale = cloc, cloc = FALSE,
    std.err = TRUE, dsm = TRUE, corr = FALSE, method = "BFGS",
    warn.inf = TRUE)
}
\arguments{
  \item{x}{A matrix or data frame with two columns, which may contain
    missing values.}
  \item{model}{The specified model; a character string. Must be
    either \code{"log"} (the default), \code{"alog"}, \code{"hr"},
    \code{"neglog"}, \code{"aneglog"}, \code{"bilog"},
    \code{"negbilog"} or \code{"ct"} (or any unique partial match),
    for the logistic, asymmetric logistic, Husler-Reiss, negative
    logistic, asymmetric negative logistic, bilogistic, negative
    bilogistic and Coles-Tawn models respectively. The definition
    of each model is given in \code{\link{rbvevd}}.}
  \item{start}{A named list giving the initial values for the
    parameters over which the likelihood is to be maximized.
    If \code{start} is omitted the routine attempts to find good
    starting values using marginal maximum likelihood estimators.}
  \item{\dots}{Additional parameters, either for the bivariate extreme
    value model or for the optimization function \code{optim}. If
    parameters of the model are included they will be held fixed at
    the values given (see \bold{Examples}).}
  \item{sym}{Logical; if \code{TRUE}, the dependence structure
    of the models \code{"alog"}, \code{"aneglog"} or \code{"ct"} are
    constrained to be symmetric (see \bold{Details}). For all other
    models, the argument is ignored (and a warning is given).}
  \item{nsloc1, nsloc2}{A data frame with the same number of rows as
    \code{x}, for linear modelling of the location parameter on the
    first/second margin (see \bold{Details}).
    The data frames are treated as covariate matrices, excluding the
    intercept. A numeric vector can be given as an alternative to a
    single column data frame.}
  \item{cshape}{Logical; if \code{TRUE}, a common shape parameter is
    fitted to each margin.}
  \item{cscale}{Logical; if \code{TRUE}, a common scale parameter is
    fitted to each margin, and the default value of \code{cshape}
    is then \code{TRUE}, so that under this default common scale
    and shape parameters are fitted.}
  \item{cloc}{Logical; if \code{TRUE}, a common location parameter is
    fitted to each margin, and the default values of \code{cshape}
    and \code{cscale} are then \code{TRUE}, so that under these
    defaults common marginal parameters are fitted.}
  \item{std.err}{Logical; if \code{TRUE} (the default), the standard
    errors are returned.}
  \item{dsm}{Logical; if \code{TRUE} (the default), summaries of the
    dependence structure are returned.}
  \item{corr}{Logical; if \code{TRUE}, the correlation matrix is
    returned.}
  \item{method}{The optimization method (see \code{\link{optim}} for
    details).}
  \item{warn.inf}{Logical; if \code{TRUE} (the default), a warning is
    given if the negative log-likelihood is infinite when evaluated at
    the starting values.}
}
\details{
  The dependence parameter names are one or more of \code{dep},
  \code{asy1}, \code{asy2}, \code{alpha} and \code{beta}, depending on
  the model selected (see \code{\link{rbvevd}}). The marginal parameter
  names are \code{loc1}, \code{scale1} and \code{shape1} for the first
  margin, and \code{loc2}, \code{scale2} and \code{shape2} for the
  second margin.
  If \code{nsloc1} is not \code{NULL}, so that a linear model is
  implemented for the first marginal location parameter, the parameter
  names for the first margin are \code{loc1}, \code{loc1}\emph{x1},
  \dots, \code{loc1}\emph{xn}, \code{scale} and \code{shape}, where
  \emph{x1}, \dots, \emph{xn} are the column names of \code{nsloc1},
  so that \code{loc1} is the intercept of the linear model, and
  \code{loc1}\emph{x1}, \dots, \code{loc1}\emph{xn} are the
  \code{ncol(nsloc1)} coefficients.
  When \code{nsloc2} is not \code{NULL}, the parameter names for the
  second margin are constructed similarly.
  
  It is recommended that the covariates within the linear models for
  the location parameters are (at least approximately) centered and
  scaled (i.e.\ that the columns of \code{nsloc1} and \code{nsloc2}
  are centered and scaled), particularly if automatic starting values
  are used, since the starting values for the associated parameters are
  then zero. If \code{cloc} is \code{TRUE}, both \code{nsloc1} and
  \code{nsloc2} must be identical, since a common linear model is
  then implemented on both margins.

  If \code{cshape} is true, the models are constrained so that
  \code{shape2 = shape1}. The parameter \code{shape2} is then
  taken to be specified, so that e.g. the common shape
  parameter can only be fixed at zero using \code{shape1 = 0},
  since using \code{shape2 = 0} gives an error. Similar
  comments apply for \code{cscale} and \code{cloc}.
  
  If \code{sym} is \code{TRUE}, the asymmetric logistic and
  asymmetric negative logistic models are constrained so that
  \code{asy2 = asy1}, and the Coles-Tawn model is constrained
  so that \code{beta = alpha}. The parameter \code{asy2} or
  \code{beta} is then taken to be specified, so that e.g.
  the parameters \code{asy1} and \code{asy2} can only
  be fixed at \code{0.8} using \code{asy1 = 0.8}, since
  using \code{asy2 = 0.8} gives an error.
  
  Bilogistic and negative bilogistic models constrained to
  symmetry are logistic and negative logistic models
  respectively. The mixed model (e.g. Tawn, 1998)
  is obtained by the asymmetric negative logistic model upon
  setting the dependence parameter to be one, and constraining
  the asymmetry parameters to be equal to each other. It can
  therefore be fitted using \code{model = "anegl"} with
  \code{dep = 1} and \code{sym = TRUE} (see \bold{Examples}).
  
  If \code{dsm} is \code{TRUE}, three values are returned which
  summarize the dependence structure, based on the fitted
  dependence function \eqn{A} (see \code{\link{abvpar}}).
  Two are measures of the strength of dependence.
  The first (Dependence One) is given by \eqn{2(1-A(1/2))}.
  The second (Dependence Two) is the integral of \eqn{4(1 - A(x))},
  taken over \eqn{0\leq x\leq1}{0 <= x <= 1}.
  Both measures are zero at independence and one at complete dependence.
  
  The third value (Asymmetry) is a measure of asymmetry, given by
  the integral of
  \eqn{4(A(x) - A(1-x))/(3 - 2\sqrt2)}{
  4(A(x) - A(1-x))/(3 - 2 sqrt(2))},
  taken over \eqn{0 \leq x \leq 0.5}{0 <= x <= 0.5}.
  This lies in the closed interval [-1,1], with
  larger absolute values representing stronger asymmetry.
  For the logistic, Husler-Reiss and negative logistic models
  \eqn{A(x) = A(1-x)} for all \eqn{0 \leq x \leq 0.5}{0 <= x <= 0.5},
  so the value will be zero.

  For numerical reasons the parameters of each model are subject the
  artificial constraints given in Table 1 of the User's Guide.  
}
\value{
  Returns an object of class \code{c("bvevd","evd")}.

  The generic accessor functions \code{\link{fitted}} (or
  \code{\link{fitted.values}}), \code{\link{std.errors}},
  \code{\link{deviance}}, \code{\link{logLik}} and
  \code{\link{AIC}} extract various features of the
  returned object.

  The functions \code{profile} and \code{profile2d} can be
  used to obtain deviance profiles.
  The function \code{anova} compares nested models, and the
  function \code{AIC} compares non-nested models.
  The function \code{plot} produces diagnostic plots.

  An object of class \code{c("bvevd","evd")} is a list containing
  the following components
  \item{estimate}{A vector containing the maximum likelihood estimates.}
  \item{std.err}{A vector containing the standard errors.}
  \item{fixed}{A vector containing the parameters that have been fixed
    at specific values within the optimization.}
  \item{fixed2}{A vector containing the parameters that have been set
    to be equal to other model parameters.}
  \item{param}{A vector containing all parameters (those optimized,
    those fixed to specific values, and those set to be
    equal to other model parameters).}
  \item{deviance}{The deviance at the maximum likelihood estimates.}
  \item{dep.summary}{A vector of three values, summarizing the
    dependence structure of the fitted model (see \bold{Details}).}
  \item{corr}{The correlation matrix.}
  \item{convergence, counts, message}{Components taken from the
    list returned by \code{\link{optim}}.}
  \item{data}{The data passed to the argument \code{x}.}
  \item{tdata}{The data, transformed to stationarity (for non-stationary
    models).}
  \item{nsloc1, nsloc2}{The arguments \code{nsloc1} and \code{nsloc2}.}
  \item{n}{The number of rows in \code{x}.}
  \item{sym}{The argument \code{sym}.}
  \item{cmar}{The vector \code{c(cloc, cscale, cshape)}.}
  \item{model}{The argument \code{model}.}
  \item{call}{The call of the current function.}
}
\section{Warning}{
  The standard errors and the correlation matrix in the returned
  object are taken from the observed information, calculated by a
  numerical approximation.
  They must be interpreted with caution when either of the
  marginal shape parameters are less than \eqn{-0.5}, because
  the usual asymptotic properties of maximum likelihood estimators
  do not then hold (Smith, 1985).
}
\references{
  Smith, R. L. (1985)
  Maximum likelihood estimation in a class of non-regular cases.
  \emph{Biometrika}, \bold{72}, 67--90.

  Tawn, J. A. (1988)
  Bivariate extreme value theory: models and estimation.
  \emph{Biometrika}, \bold{75}, 397--415.
}

\seealso{\code{\link{anova.evd}}, \code{\link{optim}},
  \code{\link{plot.bvevd}}, \code{\link{profile.evd}},
  \code{\link{profile2d.evd}}, \code{\link{rbvevd}}}

\examples{
bvdata <- rbvevd(100, dep = 0.6, model = "log", mar1 = c(1.2,1.4,0.4))
M1 <- fbvevd(bvdata, model = "log")
M2 <- fbvevd(bvdata, model = "log", dep = 0.75)
anova(M1, M2)
par(mfrow = c(2,2))
plot(M1)
plot(M1, mar = 1)
plot(M1, mar = 2)
\dontrun{par(mfrow = c(1,1))}
\dontrun{M1P <- profile(M1, which = "dep")}
\dontrun{plot(M1P)}

trend <- (-49:50)/100
rnd <- runif(100, min = -.5, max = .5)
fbvevd(bvdata, model = "log", nsloc1 = trend)
fbvevd(bvdata, model = "log", nsloc1 = trend, nsloc2 = data.frame(trend
= trend,  random = rnd))
fbvevd(bvdata, model = "log", nsloc1 = trend, nsloc2 = data.frame(trend
= trend, random = rnd), loc2random = 0)

bvdata <- rbvevd(100, dep = 1, asy = c(0.5,0.5), model = "anegl")
anlog <- fbvevd(bvdata, model = "anegl")
mixed <- fbvevd(bvdata, model = "anegl", dep = 1, sym = TRUE)
anova(anlog, mixed)
}

\keyword{models}


\eof
\name{fbvpot}
\alias{fbvpot}
\alias{print.bvpot}
\alias{plot.bvpot}

\title{Maximum-likelihood Fitting of Bivariate Extreme Value
  Distributions to Threshold Exceedances}
\description{
  Fit models for one of seven parametric bivariate extreme-value
  distributions using threshold exceedances, allowing any of the
  parameters to be held fixed if desired.
}
\synopsis{
fbvpot(x, threshold, model = c("log", "bilog", "alog", "neglog",
    "negbilog", "aneglog", "ct", "hr"), likelihood = c("censored",
    "poisson"), start, \dots, sym = FALSE, cshape = cscale,
    cscale = FALSE, std.err = TRUE, dsm = TRUE, corr = FALSE,
    method = "BFGS", warn.inf = TRUE)
}
\usage{
fbvpot(x, threshold, model = "log", likelihood = "censored", start,
    \dots, sym = FALSE, cshape = cscale, cscale = FALSE, std.err =
    TRUE, dsm = TRUE, corr = FALSE, method = "BFGS", warn.inf = TRUE)
}
\arguments{
  \item{x}{A matrix or data frame with two columns. If this
    contains missing values, those values are treated as if
    they fell below the corresponding marginal threshold.}
  \item{threshold}{A vector of two thresholds.}
  \item{model}{The specified model; a character string. Must be
    either \code{"log"} (the default), \code{"alog"}, \code{"hr"},
    \code{"neglog"}, \code{"aneglog"}, \code{"bilog"},
    \code{"negbilog"} or \code{"ct"} (or any unique partial match),
    for the logistic, asymmetric logistic, Husler-Reiss, negative
    logistic, asymmetric negative logistic, bilogistic, negative
    bilogistic and Coles-Tawn models respectively. The definition
    of each model is given in \code{\link{rbvevd}}.}
  \item{likelihood}{Unimplemented.}
  \item{start}{A named list giving the initial values for all of the
    parameters in the model. If \code{start} is omitted the routine
    attempts to find good starting values using marginal maximum
    likelihood estimators.}
  \item{\dots}{Additional parameters, either for the bivariate extreme
    value model or for the optimization function \code{optim}. If
    parameters of the model are included they will be held fixed at
    the values given (see \bold{Examples}).}
  \item{sym}{Logical; if \code{TRUE}, the dependence structure
    of the models \code{"alog"}, \code{"aneglog"} or \code{"ct"} are
    constrained to be symmetric (see \bold{Details}). For all other
    models, the argument is ignored (and a warning is given).}
  \item{cshape}{Logical; if \code{TRUE}, a common shape parameter is
    fitted to each margin.}
  \item{cscale}{Logical; if \code{TRUE}, a common scale parameter is
    fitted to each margin, and the default value of \code{cshape}
    is then \code{TRUE}, so that under this default common marginal
    parameters are fitted.}
  \item{std.err}{Logical; if \code{TRUE} (the default), the standard
    errors are returned.}
  \item{dsm}{Logical; if \code{TRUE} (the default), summaries of the
    dependence structure are returned.}
  \item{corr}{Logical; if \code{TRUE}, the correlation matrix is
    returned.}
  \item{method}{The optimization method (see \code{\link{optim}} for
    details).}
  \item{warn.inf}{Logical; if \code{TRUE} (the default), a warning is
    given if the negative log-likelihood is infinite when evaluated at
    the starting values.}
}
\details{
  The bivariate peaks over threshold models are fitted by
  maximizing the censored likelihood as given in e.g. Section
  8.3.1 of Coles(2001).
  
  The dependence parameter names are one or more of \code{dep},
  \code{asy1}, \code{asy2}, \code{alpha} and \code{beta}, depending on
  the model selected (see \code{\link{rbvevd}}). 
  The marginal parameter names are \code{scale1} and \code{shape1}
  for the first margin, and \code{scale2} and \code{shape2} for the
  second margin.

  If \code{cshape} is true, the models are constrained so that
  \code{shape2 = shape1}. The parameter \code{shape2} is then
  taken to be specified, so that e.g. the common shape
  parameter can only be fixed at zero using \code{shape1 = 0},
  since using \code{shape2 = 0} gives an error. Similar
  comments apply for \code{cscale}.
  
  If \code{sym} is \code{TRUE}, the asymmetric logistic and
  asymmetric negative logistic models are constrained so that
  \code{asy2 = asy1}, and the Coles-Tawn model is constrained
  so that \code{beta = alpha}. The parameter \code{asy2} or
  \code{beta} is then taken to be specified, so that e.g.
  the parameters \code{asy1} and \code{asy2} can only
  be fixed at \code{0.8} using \code{asy1 = 0.8}, since
  using \code{asy2 = 0.8} gives an error.
  
  Bilogistic and negative bilogistic models constrained to
  symmetry are logistic and negative logistic models
  respectively. The mixed model (e.g. Tawn, 1998)
  is obtained by the asymmetric negative logistic model upon
  setting the dependence parameter to be one, and constraining
  the asymmetry parameters to be equal to each other. It can
  therefore be fitted using \code{model = "anegl"} with
  \code{dep = 1} and \code{sym = TRUE}.
  
  If \code{dsm} is \code{TRUE}, three values are returned which
  summarize the dependence structure, based on the fitted
  dependence function \eqn{A} (see \code{\link{fbvevd}} for
  details).

  For numerical reasons the parameters of each model are subject the
  artificial constraints given in Table 1 of the User's Guide.  
}
\value{
  Returns an object of class \code{c("bvpot","evd")}.

  The generic accessor functions \code{\link{fitted}} (or
  \code{\link{fitted.values}}), \code{\link{std.errors}},
  \code{\link{deviance}}, \code{\link{logLik}} and
  \code{\link{AIC}} extract various features of the
  returned object.

  The functions \code{profile} and \code{profile2d} can be
  used to obtain deviance profiles.
  The function \code{anova} compares nested models, and the
  function \code{AIC} compares non-nested models.
  There is currently no plot method available.

  An object of class \code{c("bvpot","evd")} is a list containing
  the following components
  \item{estimate}{A vector containing the maximum likelihood estimates.}
  \item{std.err}{A vector containing the standard errors.}
  \item{fixed}{A vector containing the parameters that have been fixed
    at specific values within the optimization.}
  \item{fixed2}{A vector containing the parameters that have been set
    to be equal to other model parameters.}
  \item{param}{A vector containing all parameters (those optimized,
    those fixed to specific values, and those set to be
    equal to other model parameters).}
  \item{deviance}{The deviance at the maximum likelihood estimates.}
  \item{dep.summary}{A vector of three values, summarizing the
    dependence structure of the fitted model (see \bold{Details}).}
  \item{corr}{The correlation matrix.}
  \item{convergence, counts, message}{Components taken from the
    list returned by \code{\link{optim}}.}
  \item{data}{The data passed to the argument \code{x}.}
  \item{threshold}{The argument \code{threshold}.}
  \item{n}{The number of rows in \code{x}.}
  \item{nat}{The vector of length three containing the number of
    exceedances on the first, second and both margins respectively.}
  \item{sym}{The argument \code{sym}.}
  \item{cmar}{The vector \code{c(cscale, cshape)}.}
  \item{model}{The argument \code{model}.}
  \item{call}{The call of the current function.}
}
\section{Warning}{
  The standard errors and the correlation matrix in the returned
  object are taken from the observed information, calculated by a
  numerical approximation.
  They must be interpreted with caution when either of the
  marginal shape parameters are less than \eqn{-0.5}, because
  the usual asymptotic properties of maximum likelihood estimators
  do not then hold (Smith, 1985).
}
\references{
  Coles, S. G. (2001)
  \emph{An Introduction to Statistical Modelling of Extreme
  Values}, London: Springer--Verlag.

  Smith, R. L. (1985)
  Maximum likelihood estimation in a class of non-regular cases.
  \emph{Biometrika}, \bold{72}, 67--90.
}
\author{Chris Ferro and Alec Stephenson}

\seealso{\code{\link{abvpar}}, \code{\link{anova.evd}},
  \code{\link{fbvevd}}, \code{\link{optim}}, \code{\link{rbvevd}}}

\examples{
bvdata <- rbvevd(1000, dep = 0.5, model = "log")
u <- apply(bvdata, 2,  quantile, probs = 0.9)
M1 <- fbvpot(bvdata, u, model = "log")
M2 <- fbvpot(bvdata, u, "log", dep = 0.5)
anova(M1, M2)
}

\keyword{models}


\eof
\name{fextreme}
\alias{fextreme}
\title{Maximum-likelihood Fitting of Maxima and Minima}
\description{
  Maximum-likelihood fitting for the distribution of the
  maximum/minimum of a given number of independent variables from
  a specified distribution.
}
\usage{
fextreme(x, start, densfun, distnfun, \dots, distn, mlen = 1, largest =
    TRUE, std.err = TRUE, corr = FALSE, method = "Nelder-Mead")
}
\arguments{
  \item{x}{A numeric vector.}
  \item{start}{A named list giving the initial values for the
    parameters over which the likelihood is to be maximized.}
  \item{densfun, distnfun}{Density and distribution function of the
    specified distribution.}
  \item{\dots}{Additional parameters, either for the specified
    distribution or for the optimization function \code{optim}.
    If parameters of the distribution are included they will be
    held fixed at the values given (see \bold{Examples}). If
    parameters of the distribution are not included either here
    or as a named component in \code{start} they will be held
    fixed at the default values specified in the corresponding
    density and distribution functions (assuming they exist; an
    error will be generated otherwise).}
  \item{distn}{A character string, optionally specified as an alternative
    to \code{densfun} and \code{distnfun} such that the density and
    distribution functions are formed upon the addition of the
    prefixes \code{d} and \code{p} respectively.}
  \item{mlen}{The number of independent variables.}
  \item{largest}{Logical; if \code{TRUE} (default) use maxima,
    otherwise minima.}
  \item{std.err}{Logical; if \code{TRUE} (the default), the standard
    errors are returned.}
  \item{corr}{Logical; if \code{TRUE}, the correlation matrix is
    returned.}
  \item{method}{The optimization method (see \code{\link{optim}} for
    details).}
}
\details{
  Maximization of the log-likelihood is performed.
  The estimated standard errors are taken from the observed information,
  calculated by a numerical approximation. 

  If the density and distribution functions are user defined, the order
  of the arguments must mimic those in R base (i.e. data first,
  parameters second).
  Density functions must have \code{log} arguments.
}
\value{
  Returns an object of class \code{c("extreme","evd")}.

  The generic accessor functions \code{\link{fitted}} (or
  \code{\link{fitted.values}}), \code{\link{std.errors}},
  \code{\link{deviance}}, \code{\link{logLik}} and
  \code{\link{AIC}} extract various features of the
  returned object.
  The function \code{anova} compares nested models.
  
  An object of class \code{c("extreme","evd")} is a list containing
  at most the following components
  \item{estimate}{A vector containing the maximum likelihood estimates.}
  \item{std.err}{A vector containing the standard errors.}
  \item{deviance}{The deviance at the maximum likelihood estimates.}
  \item{corr}{The correlation matrix.}
  \item{convergence, counts, message}{Components taken from the
    list returned by \code{\link{optim}}.}
  \item{call}{The call of the current function.}
  \item{data}{The data passed to the argument \code{x}.}
  \item{n}{The length of \code{x}.}
}

\seealso{\code{\link{anova.evd}}, \code{\link{forder}}, \code{\link{optim}}}

\examples{
uvdata <- rextreme(100, qnorm, mean = 0.56, mlen = 365)
fextreme(uvdata, list(mean = 0, sd = 1), distn = "norm", mlen = 365)
fextreme(uvdata, list(rate = 1), distn = "exp", mlen = 365)
fextreme(uvdata, list(scale = 1), shape = 1, distn = "gamma", mlen = 365)
fextreme(uvdata, list(shape = 1, scale = 1), distn = "gamma", mlen = 365)
}
\keyword{models}



\eof
\name{fgev}
\alias{fgev}
\alias{fitted.evd}
\alias{std.errors}
\alias{std.errors.evd}
\alias{print.evd}
\alias{logLik.evd}

\title{Maximum-likelihood Fitting of the Generalized Extreme Value
  Distribution}
\description{
  Maximum-likelihood fitting for the generalized extreme value distribution,
  including linear modelling of the location parameter, and allowing any
  of the parameters to be held fixed if desired.
}
\usage{
fgev(x, start, \dots, nsloc = NULL, prob = NULL, std.err = TRUE,
    corr = FALSE, method = "BFGS", warn.inf = TRUE)
}
\arguments{
  \item{x}{A numeric vector, which may contain missing values.}
  \item{start}{A named list giving the initial values for the
    parameters over which the likelihood is to be maximized.
    If \code{start} is omitted the routine attempts to find good
    starting values using moment estimators.}
  \item{\dots}{Additional parameters, either for the GEV model
    or for the optimization function \code{optim}. If parameters
    of the model are included they will be held fixed at the
    values given (see \bold{Examples}).}
  \item{nsloc}{A data frame with the same number of rows as the
    length of \code{x}, for linear modelling of the location
    parameter.
    The data frame is treated as a covariate matrix (excluding the
    intercept).
    A numeric vector can be given as an alternative to a single column
    data frame.}
  \item{prob}{Controls the parameterization of the model (see
    \bold{Details}). Should be either \code{NULL} (the default),
    or a probability in the closed interval [0,1].}
  \item{std.err}{Logical; if \code{TRUE} (the default), the standard
    errors are returned.}
  \item{corr}{Logical; if \code{TRUE}, the correlation matrix is
    returned.}
  \item{method}{The optimization method (see \code{\link{optim}} for
    details).}
  \item{warn.inf}{Logical; if \code{TRUE} (the default), a warning is
    given if the negative log-likelihood is infinite when evaluated at
    the starting values.}
}
\details{
  If \code{prob} is \code{NULL} (the default):
  
  For stationary models the parameter names are \code{loc}, \code{scale}
  and \code{shape}, for the location, scale and shape parameters
  respectively.
  For non-stationary models, the parameter names are \code{loc},
  \code{loc}\emph{x1}, \dots, \code{loc}\emph{xn}, \code{scale} and
  \code{shape}, where \emph{x1}, \dots, \emph{xn} are the column names
  of \code{nsloc}, so that \code{loc} is the intercept of the
  linear model, and \code{loc}\emph{x1}, \dots, \code{loc}\emph{xn}
  are the \code{ncol(nsloc)} coefficients.
  If \code{nsloc} is a vector it is converted into a single column
  data frame with column name \code{trend}, and hence the associated
  trend parameter is named \code{loctrend}.

  If \eqn{\code{prob} = p} is a probability:
  
  The fit is performed using a different parameterization.
  Let \eqn{a}, \eqn{b} and \eqn{s} denote the location, scale
  and shape parameters of the GEV distribution.
  For stationary models, the distribution is parameterized
  using \eqn{(z_p,b,s)}, where
  \deqn{z_p = a - b/s (1 - (-\log(1 - p))^s)}{
        z_p = a - b/s (1 - (-log(1 - p))^s)}
  is such that \eqn{G(z_p) = 1 - p}, where \eqn{G} is the
  GEV distribution function.
  \eqn{\code{prob} = p} is therefore the probability in the upper
  tail corresponding to the quantile \eqn{z_p}.
  If \code{prob} is zero, then \eqn{z_p} is the upper end point
  \eqn{a - b/s}, and \eqn{s} is restricted to the negative
  (Weibull) axis.
  If \code{prob} is one, then \eqn{z_p} is the lower end point
  \eqn{a - b/s}, and \eqn{s} is restricted to the positive
  (Frechet) axis.
  The parameter names are \code{quantile}, \code{scale}
  and \code{shape}, for \eqn{z_p}, \eqn{b} and \eqn{s}
  respectively.
  
  For non-stationary models the parameter \eqn{z_p} is again given by
  the equation above, but \eqn{a} becomes the intercept of the linear
  model for the location parameter, so that \code{quantile} replaces
  (the intercept) \code{loc}, and hence the parameter names are
  \code{quantile}, \code{loc}\emph{x1}, \dots, \code{loc}\emph{xn},
  \code{scale} and \code{shape}, where \emph{x1}, \dots, \emph{xn} are
  the column names of \code{nsloc}.

  In either case:
  
  For non-stationary fitting it is recommended that the covariates
  within the linear model for the location parameter are (at least
  approximately) centered and scaled (i.e.\ that the columns of
  \code{nsloc} are centered and scaled), particularly if automatic
  starting values are used, since the starting values for the
  associated parameters are then zero.
}

\value{
  Returns an object of class \code{c("gev","uvevd","evd")}.

  The generic accessor functions \code{\link{fitted}} (or
  \code{\link{fitted.values}}), \code{\link{std.errors}},
  \code{\link{deviance}}, \code{\link{logLik}} and
  \code{\link{AIC}} extract various features of the
  returned object.

  The functions \code{profile} and \code{profile2d} are
  used to obtain deviance profiles for the model parameters.
  In particular, profiles of the quantile \eqn{z_p} can be
  calculated and plotted when \eqn{\code{prob} = p}.
  The function \code{anova} compares nested models.
  The function \code{plot} produces diagnostic plots.
  
  An object of class \code{c("gev","uvevd","evd")} is a list
  containing at most the following components
  \item{estimate}{A vector containing the maximum likelihood estimates.}
  \item{std.err}{A vector containing the standard errors.}
  \item{fixed}{A vector containing the parameters of the model that
    have been held fixed.}
  \item{param}{A vector containing all parameters (optimized and fixed).}
  \item{deviance}{The deviance at the maximum likelihood estimates.}
  \item{corr}{The correlation matrix.}
  \item{convergence, counts, message}{Components taken from the
    list returned by \code{\link{optim}}.}
  \item{data}{The data passed to the argument \code{x}.}
  \item{tdata}{The data, transformed to stationarity (for non-stationary
    models).}
  \item{nsloc}{The argument \code{nsloc}.}
  \item{n}{The length of \code{x}.}
  \item{prob}{The argument \code{prob}.}
  \item{loc}{The location parameter. If \code{prob} is \code{NULL}
    (the default), this will also be an element of \code{param}.}
  \item{call}{The call of the current function.}
}
\section{Warning}{
  The standard errors and the correlation matrix in the returned
  object are taken from the observed information, calculated by a
  numerical approximation.
  They must be interpreted with caution when the shape parameter
  is less than \eqn{-0.5}, because the usual asymptotic
  properties of maximum likelihood estimators do not then
  hold (Smith, 1985).
}
\references{
  Smith, R. L. (1985)
  Maximum likelihood estimation in a class of non-regular cases.
  \emph{Biometrika}, \bold{72}, 67--90.
}

\seealso{\code{\link{anova.evd}}, \code{\link{optim}},
  \code{\link{plot.uvevd}}, \code{\link{profile.evd}},
  \code{\link{profile2d.evd}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
trend <- (-49:50)/100
M1 <- fgev(uvdata, nsloc = trend, control = list(trace = 1))
M2 <- fgev(uvdata)
M3 <- fgev(uvdata, shape = 0)
M4 <- fgev(uvdata, scale = 1, shape = 0)
anova(M1, M2, M3, M4)
par(mfrow = c(2,2))
plot(M2)
\dontrun{M2P <- profile(M2)}
\dontrun{plot(M2P)}

rnd <- runif(100, min = -.5, max = .5)
fgev(uvdata, nsloc = data.frame(trend = trend, random = rnd))
fgev(uvdata, nsloc = data.frame(trend = trend, random = rnd), locrandom = 0)

uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata, prob = 0.1)
M2 <- fgev(uvdata, prob = 0.01)
\dontrun{M1P <- profile(M1, which = "quantile")}
\dontrun{M2P <- profile(M2, which = "quantile")}
\dontrun{plot(M1P)}
\dontrun{plot(M2P)}
}
\keyword{models}


\eof
\name{forder}
\alias{forder}
\title{Maximum-likelihood Fitting of Order Statistics}
\description{
  Maximum-likelihood fitting for the distribution of a selected order
  statistic of a given number of independent variables from a
  specified distribution.
}
\usage{
forder(x, start, densfun, distnfun, \dots, distn, mlen = 1, j = 1, 
largest = TRUE, std.err = TRUE, corr = FALSE, method = "Nelder-Mead")
}
\arguments{
  \item{x}{A numeric vector.}
  \item{start}{A named list giving the initial values for the
    parameters over which the likelihood is to be maximized.}
  \item{densfun, distnfun}{Density and distribution function of the
    specified distribution.}
  \item{\dots}{Additional parameters, either for the specified
    distribution or for the optimization function \code{optim}.
    If parameters of the distribution are included they will be
    held fixed at the values given (see \bold{Examples}). If
    parameters of the distribution are not included either here
    or as a named component in \code{start} they will be held
    fixed at the default values specified in the corresponding
    density and distribution functions (assuming they exist; an
    error will be generated otherwise).}
  \item{distn}{A character string, optionally specified as an alternative
    to \code{densfun} and \code{distnfun} such that the density and
    distribution and functions are formed upon the addition of the
    prefixes \code{d} and \code{p} respectively.}
  \item{mlen}{The number of independent variables.}
  \item{j}{The order statistic, taken as the \code{j}th largest
    (default) or smallest of \code{mlen}, according to the value of
    \code{largest}.}
  \item{largest}{Logical; if \code{TRUE} (default) use the \code{j}th
    largest order statistic, otherwise use the \code{j}th smallest.}
  \item{std.err}{Logical; if \code{TRUE} (the default), the standard
    errors are returned.}
  \item{corr}{Logical; if \code{TRUE}, the correlation matrix is
    returned.}
  \item{method}{The optimization method (see \code{\link{optim}} for
    details).}
}
\details{
  Maximization of the log-likelihood is performed.
  The estimated standard errors are taken from the observed information,
  calculated by a numerical approximation.

  If the density and distribution functions are user defined, the order
  of the arguments must mimic those in R base (i.e. data first,
  parameters second).
  Density functions must have \code{log} arguments.
}
\value{
  Returns an object of class \code{c("extreme","evd")}.
  This class is defined in \code{\link{fextreme}}.

  The generic accessor functions \code{\link{fitted}} (or
  \code{\link{fitted.values}}), \code{\link{std.errors}},
  \code{\link{deviance}}, \code{\link{logLik}} and
  \code{\link{AIC}} extract various features of the
  returned object.
  The function \code{anova} compares nested models.
}

\seealso{\code{\link{anova.evd}}, \code{\link{fextreme}},
  \code{\link{optim}}}

\examples{
uvd <- rorder(100, qnorm, mean = 0.56, mlen = 365, j = 2)
forder(uvd, list(mean = 0, sd = 1), distn = "norm", mlen = 365, j = 2)
forder(uvd, list(rate = 1), distn = "exp", mlen = 365, j = 2)
forder(uvd, list(scale = 1), shape = 1, distn = "gamma", mlen = 365, j = 2)
forder(uvd, list(shape = 1, scale = 1), distn = "gamma", mlen = 365, j = 2)
}
\keyword{models}



\eof
\name{fox}
\alias{fox}
\title{Maximum Annual Flood Discharges of the Fox River}
\usage{data(fox)}
\description{
  The \code{fox} data frame has 33 rows and 2 columns.
  The columns contain maximum annual flood discharges, in units
  of 1000 cubed feet per second, from the Fox River in Wisconsin,
  USA at Berlin (upstream) and Wrightstown (downstream), for the
  years 1918 to 1950.
  The row names give the years of observation.
}
\format{
  This data frame contains the following columns:
  \describe{
    \item{berlin}{A numeric vector containing maximum annual flood
      discharges at Berlin (upstream).}
    \item{wright}{A numeric vector containing maximum annual flood
      discharges at Wrightstown (downstream).}
  }
}
\source{
  Gumbel, E. J. and Mustafi, C. K. (1967)
  Some analytical properties of bivariate extremal distributions.
  \emph{J. Amer. Statist. Assoc.}, \bold{62}, 569--588.
}
\keyword{datasets}

\eof
\name{fpot}
\alias{fpot}
\alias{print.pot}
\title{Peaks Over Threshold Modelling using the Generalized Pareto
  or Point Process Representation}
\description{
  Maximum-likelihood fitting for peaks over threshold modelling,
  using the Generalized Pareto or Point Process representation,
  allowing any of the parameters to be held fixed if desired.
}
\usage{
fpot(x, threshold, model = c("gpd", "pp"), start, npp = length(x),
    cmax = FALSE, r = 1, ulow = -Inf, rlow = 1, mper = NULL, \dots,
    std.err = TRUE, corr = FALSE, method = "BFGS", warn.inf = TRUE)
}
\arguments{
  \item{x}{A numeric vector. If this contains missing values,
    those values are treated as if they fell below the
    threshold.}
  \item{threshold}{The threshold.}
  \item{model}{The model; either \code{"gpd"} (the default) or
    \code{"pp"}, for the Generalized Pareto or Point Process
    representations respectively.}
  \item{start}{A named list giving the initial values for the
    parameters over which the likelihood is to be maximized.
    If \code{start} is omitted the routine attempts to find good
    starting values using moment estimators.}
  \item{npp}{The data should contain \code{npp} observations per
    ``period'', where the return level plot produced by
    \code{plot.pot} will represent return periods in units of
    ``periods''. By default \code{npp = length(x)}, so that the
    ``period'' is the period of time over which the entire
    data set is collected. It may often be useful to change
    this default so that more sensible units are used.
    For example, if yearly periodic units are required, use
    \code{npp = 365.25} for daily data and \code{npp = 52.18}
    for weekly data. The argument only makes a difference to
    the actual fit if \code{mper} is not \code{NULL} or if
    \code{model = "pp"} (see \bold{Details}).}
  \item{cmax}{Logical; if \code{FALSE} (the default), the model
    is fitted using all exceedences over the threshold. If
    \code{TRUE}, the model is fitted using cluster maxima, using
    clusters of exceedences derived from \code{clusters}.}
  \item{r, ulow, rlow}{Arguments used for the identification of
    clusters of exceedences (see \code{\link{clusters}}).
    Ignored if \code{cmax} is \code{FALSE} (the default).}
  \item{mper}{Controls the parameterization of the generalized
    Pareto model. Should be either \code{NULL} (the default),
    or a positive number (see \bold{Details}).
    If \code{mper} is not \code{NULL} and \code{model = "pp"},
    an error is returned.}
  \item{\dots}{Additional parameters, either for the model
    or for the optimization function \code{optim}. If parameters
    of the model are included they will be held fixed at the
    values given (see \bold{Examples}).}
  \item{std.err}{Logical; if \code{TRUE} (the default), the standard
    errors are returned.}
  \item{corr}{Logical; if \code{TRUE}, the correlation matrix is
    returned.}
  \item{method}{The optimization method (see \code{\link{optim}} for
    details).}
  \item{warn.inf}{Logical; if \code{TRUE} (the default), a warning is
    given if the negative log-likelihood is infinite when evaluated at
    the starting values.}
}
\details{
  The exeedances over the threshold \code{threshold} (if \code{cmax} is
  \code{FALSE}) or the maxima of the clusters of exeedances (if
  \code{cmax} is \code{TRUE}) are (if \code{model = "gpd"}) fitted to a
  generalized Pareto distribution (GPD) with location \code{threshold}.
  If \code{model = "pp"} the exceedances are fitted to a
  non-homogeneous Poisson process (Coles, 2001).
  
  If \code{mper} is \code{NULL} (the default), the parameters of
  the model (if \code{model = "gpd"}) are \code{scale} and
  \code{shape}, for the scale and shape parameters of the GPD.
  If \code{model = "pp"} the parameters are \code{loc}, \code{scale}
  and \code{shape}. Under \code{model = "pp"} the parameters can be
  interpreted as parameters of the Generalized Extreme Value
  distribution, fitted to the maxima of \code{npp} random variables.
  In this case, the value of \code{npp} should be reasonably large.

  For both characterizations, the shape parameters are
  equivalent. The scale parameter under the generalized Pareto
  characterization is equal to \eqn{b + s(u - a)}, where \eqn{a},
  \eqn{b} and \eqn{s} are the location, scale and shape parameters
  under the Point Process characterization, and where \eqn{u} is
  the threshold.

  If \eqn{\code{mper} = m} is a positive value, then
  the generalized Pareto model is reparameterized so that the
  parameters are \code{rlevel} and \code{shape}, where
  \code{rlevel} is the \eqn{m} ``period'' return level, where
  ``period'' is defined via the argument \code{npp}.

  The \eqn{m} ``period'' return level is defined as follows.
  Let \eqn{G} be the fitted generalized Pareto distribution
  function, with location \eqn{\code{threshold} = u}, so that
  \eqn{1 - G(z)} is the fitted probability of an exceedance
  over \eqn{z > u} given an exceedance over \eqn{u}.
  The fitted probability of an exceedance over \eqn{z > u} is
  therefore \eqn{p(1 - G(z))}, where \eqn{p} is the estimated
  probabilty of exceeding \eqn{u}, which is given by the empirical
  proportion of exceedances.
  The \eqn{m} ``period'' return level \eqn{z_m} satisfies
  \eqn{p(1 - G(z_m)) = 1/(mN)}, where \eqn{N} is the number
  of points per period (multiplied by the estimate of the
  extremal index, if cluster maxima are fitted).
  In other words, \eqn{z_m} is the quantile of the fitted model
  that corresponds to the upper tail probability \eqn{1/(mN)}.
  If \code{mper} is infinite, then \eqn{z_m} is the upper end point,
  given by \code{threshold} minus \eqn{\code{scale}/\code{shape}},
  and the shape parameter is then restricted to be negative.
}

\value{
  Returns an object of class \code{c("pot","uvevd","pot")}.

  The generic accessor functions \code{\link{fitted}} (or
  \code{\link{fitted.values}}), \code{\link{std.errors}},
  \code{\link{deviance}}, \code{\link{logLik}} and
  \code{\link{AIC}} extract various features of the
  returned object.

  The function \code{profile} can be
  used to obtain deviance profiles for the model parameters.
  In particular, profiles of the \eqn{m} \code{period}
  return level \eqn{z_m} can be calculated and plotted when
  \eqn{\code{mper} = m}.
  The function \code{anova} compares nested models.
  The function \code{plot} produces diagnostic plots.
  
  An object of class \code{c("pot","uvevd","pot")} is a list containing
  the following components
  \item{estimate}{A vector containing the maximum likelihood estimates.}
  \item{std.err}{A vector containing the standard errors.}
  \item{fixed}{A vector containing the parameters of the model that
    have been held fixed.}
  \item{param}{A vector containing all parameters (optimized and fixed).}
  \item{deviance}{The deviance at the maximum likelihood estimates.}
  \item{corr}{The correlation matrix.}
  \item{convergence, counts, message}{Components taken from the
    list returned by \code{\link{optim}}.}
  \item{threshold, r, ulow, rlow, npp}{The arguments of the same
    name.}
  \item{nhigh}{The number of exceedences (if \code{cmax} is
    \code{FALSE}) or the number of clusters of exceedences (if
    \code{cmax} is \code{TRUE}).}
  \item{nat, pat}{The number and proportion of exceedences.}
  \item{extind}{The estimate of the extremal index (i.e.
    \code{nhigh} divided by \code{nat}). If \code{cmax} is
    \code{FALSE}, this is \code{NULL}.}
  \item{data}{The data passed to the argument \code{x}.}
  \item{exceedances}{The exceedences, or the maxima of the clusters
    of exceedences.}
  \item{mper}{The argument \code{mper}.}
  \item{scale}{The scale parameter for the fitted generalized Pareto
    distribution. If \code{mper} is \code{NULL} and \code{model = "gpd"}
    (the defaults), this will also be an element of \code{param}.}
  \item{call}{The call of the current function.}
}
\section{Warning}{
  The standard errors and the correlation matrix in the returned
  object are taken from the observed information, calculated by a
  numerical approximation.
  They must be interpreted with caution when the shape parameter
  is less than \eqn{-0.5}, because the usual asymptotic
  properties of maximum likelihood estimators do not then
  hold (Smith, 1985).
}
\references{
  Smith, R. L. (1985)
  Maximum likelihood estimation in a class of non-regular cases.
  \emph{Biometrika}, \bold{72}, 67--90.
}

\seealso{\code{\link{anova.evd}}, \code{\link{optim}},
  \code{\link{plot.uvevd}}, \code{\link{profile.evd}},
  \code{\link{profile2d.evd}}, \code{\link{mrlplot}},
  \code{\link{tcplot}}}

\examples{
uvdata <- rgpd(100, loc = 0, scale = 1.1, shape = 0.2)
M1 <- fpot(uvdata, 1)
M2 <- fpot(uvdata, 1, shape = 0)
anova(M1, M2)
par(mfrow = c(2,2))
plot(M1)
\dontrun{M1P <- profile(M1)}
\dontrun{plot(M1P)}

M1 <- fpot(uvdata, 1, mper = 10)
M2 <- fpot(uvdata, 1, mper = 100)
\dontrun{M1P <- profile(M1, which = "rlevel", conf=0.975, mesh=0.1)}
\dontrun{M2P <- profile(M2, which = "rlevel", conf=0.975, mesh=0.1)}
\dontrun{plot(M1P)}
\dontrun{plot(M2P)}
}
\keyword{models}


\eof
\name{frechet}
\alias{dfrechet}
\alias{pfrechet}
\alias{qfrechet}
\alias{rfrechet}
\title{The Frechet Distribution}
\description{
  Density function, distribution function, quantile function
  and random generation for the Frechet distribution with
  location, scale and shape parameters. 
}
\usage{
dfrechet(x, loc=0, scale=1, shape=1, log = FALSE) 
pfrechet(q, loc=0, scale=1, shape=1, lower.tail = TRUE) 
qfrechet(p, loc=0, scale=1, shape=1, lower.tail = TRUE)
rfrechet(n, loc=0, scale=1, shape=1)
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{p}{Vector of probabilities.}
  \item{n}{Number of observations.}
  \item{loc, scale, shape}{Location, scale and shape parameters (can be
    given as vectors).}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x]}
}
\details{
  The Frechet distribution function with parameters
  \eqn{\code{loc} = a}, \eqn{\code{scale} = b} and
  \eqn{\code{shape} = s} is
  \deqn{G(z) = \exp\left\{-\left(\frac{z-a}{b}\right)^{-s}
  \right\}}{G(x) = exp{-[(z-a)/b]^(-s)}}
  for \eqn{z > a} and zero otherwise, where \eqn{b > 0} and
  \eqn{s > 0}.
}
\value{
  \code{dfrechet} gives the density function, \code{pfrechet} gives
  the distribution function, \code{qfrechet} gives the quantile
  function, and \code{rfrechet} generates random deviates.
}

\seealso{\code{\link{rgev}}, \code{\link{rgumbel}}, \code{\link{rrweibull}}}

\examples{
dfrechet(2:4, 1, 0.5, 0.8)
pfrechet(2:4, 1, 0.5, 0.8)
qfrechet(seq(0.9, 0.6, -0.1), 2, 0.5, 0.8)
rfrechet(6, 1, 0.5, 0.8)
p <- (1:9)/10
pfrechet(qfrechet(p, 1, 2, 0.8), 1, 2, 0.8)
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
}
\keyword{distribution}


\eof
\name{gev}
\alias{dgev}
\alias{pgev}
\alias{qgev}
\alias{rgev}
\title{The Generalized Extreme Value Distribution}
\description{
  Density function, distribution function, quantile function and
  random generation for the generalized extreme value (GEV)
  distribution with location, scale and shape parameters. 
}
\usage{
dgev(x, loc=0, scale=1, shape=0, log = FALSE) 
pgev(q, loc=0, scale=1, shape=0, lower.tail = TRUE) 
qgev(p, loc=0, scale=1, shape=0, lower.tail = TRUE)
rgev(n, loc=0, scale=1, shape=0)
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{p}{Vector of probabilities.}
  \item{n}{Number of observations.}
  \item{loc, scale, shape}{Location, scale and shape parameters; the 
    \code{shape} argument cannot be a vector (must have length one).}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x]}
}
\details{
  The GEV distribution function with parameters
  \eqn{\code{loc} = a}, \eqn{\code{scale} = b} and
  \eqn{\code{shape} = s} is
  \deqn{G(z) = \exp\left[-\{1+s(z-a)/b\}^{-1/s}\right]}{
    G(x) = exp[-{1+s(z-a)/b}^(-1/s)]}
  for \eqn{1+s(z-a)/b > 0}, where \eqn{b > 0}.
  If \eqn{s = 0} the distribution is defined by continuity.
  If \eqn{1+s(z-a)/b \leq 0}{1+s(z-a)/b <= 0}, the value \eqn{z} is
  either greater than the upper end point (if \eqn{s < 0}), or less
  than the lower end point (if \eqn{s > 0}).
  
  The parametric form of the GEV encompasses that of the Gumbel,
  Frechet and reversed Weibull distributions, which are obtained
  for \eqn{s = 0}, \eqn{s > 0} and \eqn{s < 0} respectively.
  It was first introduced by Jenkinson (1955).
}
\value{
  \code{dgev} gives the density function, \code{pgev} gives the
  distribution function, \code{qgev} gives the quantile function,
  and \code{rgev} generates random deviates.
}
\references{ 
  Jenkinson, A. F. (1955)
  The frequency distribution of the annual maximum (or minimum) of
  meteorological elements.
  \emph{Quart. J. R. Met. Soc.}, \bold{81}, 158--171.
}	

\seealso{\code{\link{fgev}}, \code{\link{rfrechet}},
  \code{\link{rgumbel}}, \code{\link{rrweibull}}}

\examples{
dgev(2:4, 1, 0.5, 0.8)
pgev(2:4, 1, 0.5, 0.8)
qgev(seq(0.9, 0.6, -0.1), 2, 0.5, 0.8)
rgev(6, 1, 0.5, 0.8)
p <- (1:9)/10
pgev(qgev(p, 1, 2, 0.8), 1, 2, 0.8)
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
}
\keyword{distribution}


\eof
\name{gpd}
\alias{dgpd}
\alias{pgpd}
\alias{qgpd}
\alias{rgpd}
\title{The Generalized Pareto Distribution}
\description{
  Density function, distribution function, quantile function and
  random generation for the generalized Pareto distribution (GPD)
  with location, scale and shape parameters. 
}
\usage{
dgpd(x, loc=0, scale=1, shape=0, log = FALSE) 
pgpd(q, loc=0, scale=1, shape=0, lower.tail = TRUE) 
qgpd(p, loc=0, scale=1, shape=0, lower.tail = TRUE)
rgpd(n, loc=0, scale=1, shape=0)
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{p}{Vector of probabilities.}
  \item{n}{Number of observations.}
  \item{loc, scale, shape}{Location, scale and shape parameters; the 
    \code{shape} argument cannot be a vector (must have length one).}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x]}
}
\details{
  The generalized Pareto distribution function (Pickands, 1975) with
  parameters \eqn{\code{loc} = a}, \eqn{\code{scale} = b} and
  \eqn{\code{shape} = s} is
  \deqn{G(z) = 1 - \{1+s(z-a)/b\}^{-1/s}}{
    G(z) = 1 - {1+s(z-a)/b}^(-1/s)}
  for \eqn{1+s(z-a)/b > 0} and \eqn{z > a}, where \eqn{b > 0}.
  If \eqn{s = 0} the distribution is defined by continuity.
}
\value{
  \code{dgpd} gives the density function, \code{pgpd} gives the
  distribution function, \code{qgpd} gives the quantile function,
  and \code{rgpd} generates random deviates.
}
\references{ 
  Pickands, J. (1975)
  Statistical inference using extreme order statistics.
  \emph{Annals of Statistics}, \bold{3}, 119--131.
}	

\seealso{\code{\link{fpot}}, \code{\link{rgev}}}

\examples{
dgpd(2:4, 1, 0.5, 0.8)
pgpd(2:4, 1, 0.5, 0.8)
qgpd(seq(0.9, 0.6, -0.1), 2, 0.5, 0.8)
rgpd(6, 1, 0.5, 0.8)
p <- (1:9)/10
pgpd(qgpd(p, 1, 2, 0.8), 1, 2, 0.8)
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
}
\keyword{distribution}


\eof
\name{gumbel}
\alias{dgumbel}
\alias{pgumbel}
\alias{qgumbel}
\alias{rgumbel}
\title{The Gumbel Distribution}
\description{
  Density function, distribution function, quantile function and
  random generation for the Gumbel distribution with location and
  scale parameters. 
}
\usage{
dgumbel(x, loc=0, scale=1, log = FALSE) 
pgumbel(q, loc=0, scale=1, lower.tail = TRUE) 
qgumbel(p, loc=0, scale=1, lower.tail = TRUE)
rgumbel(n, loc=0, scale=1)
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{p}{Vector of probabilities.}
  \item{n}{Number of observations.}
  \item{loc, scale}{Location and scale parameters (can be given
    as vectors).}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x]}
}
\details{
  The Gumbel distribution function with parameters
  \eqn{\code{loc} = a} and \eqn{\code{scale} = b} is
  \deqn{G(z) = \exp\left\{-\exp\left[-\left(\frac{z-a}{b}\right)
    \right]\right\}}{G(x) = exp{-exp[-(z-a)/b]}}
  for all real \eqn{z}, where \eqn{b > 0}.
}
\value{
  \code{dgumbel} gives the density function, \code{pgumbel} gives
  the distribution function, \code{qgumbel} gives the quantile
  function, and \code{rgumbel} generates random deviates.
}
\section{Interesting Fact}{
  Gumbel is the name of my teddy bear.
  And a famous statistician.
}
\seealso{\code{\link{rfrechet}}, \code{\link{rgev}}, \code{\link{rrweibull}}}

\examples{
dgumbel(-1:2, -1, 0.5)
pgumbel(-1:2, -1, 0.5)
qgumbel(seq(0.9, 0.6, -0.1), 2, 0.5)
rgumbel(6, -1, 0.5)
p <- (1:9)/10
pgumbel(qgumbel(p, -1, 2), -1, 2)
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
}
\keyword{distribution}



\eof
\name{lisbon}
\alias{lisbon}
\title{Annual Maximum Wind Speeds at Lisbon}
\usage{data(lisbon)}
\description{
  A numeric vector containing annual maximum wind speeds, in
  kilometers per hour, from 1941 to 1970 at Lisbon, Portugal.
}
\format{A vector containing 30 observations.}
\source{
  Tiago de Oliveira, J. (1997)
  \emph{Statistical Analysis of Extremes.}
  Pendor.
}
\keyword{datasets}

\eof
\name{marma}
\alias{marma}
\alias{mar}
\alias{mma}
\title{Simulate MARMA(p,q) Processes}
\description{
  Simulation of MARMA(p,q) processes. 
}
\usage{
marma(n, p = 0, q = 0, psi, theta, init = rep(0, p), n.start = p,
    rand.gen = rfrechet, \dots)
mar(n, p = 1, psi, init = rep(0, p), n.start = p, rand.gen =
    rfrechet, \dots)
mma(n, q = 1, theta, rand.gen = rfrechet, \dots)
}
\arguments{
  \item{n}{The number of observations.}
  \item{p}{The AR order of the MARMA process.}
  \item{q}{The MA order of the MARMA process.}
  \item{psi}{A vector of non-negative parameters, of length
    \code{p}. Can be omitted if \code{p} is zero.}
  \item{theta}{A vector of non-negative parameters, of length
    \code{q}. Can be omitted if \code{q} is zero.}
  \item{init}{A vector of non-negative starting values, of
    length \code{p}.}
  \item{n.start}{A non-negative value denoting the length of the
    burn-in period. If \code{n.start} is less than \code{p}, then
    \code{p} minus \code{n.start} starting values will be included
    in the output series.}
  \item{rand.gen}{A simulation function to generate the
    innovations.}
  \item{\dots}{Additional arguments for \code{rand.gen}.  Most
    usefully, the scale and shape parameters of the innovations
    generated by \code{rfrechet} can be specified by \code{scale}
    and \code{shape} respectively.}
}
\details{
  A max autoregressive moving average process \eqn{\{X_k\}}{{X_k}},
  denoted by MARMA(p,q), satisfies
  \deqn{X_k = \max\{\phi_1 X_{k-1}, \ldots, \phi_p X_{k-p}, \epsilon_k,
    \theta_1 \epsilon_{k-1}, \ldots, \theta_q \epsilon_{k-q}\}}{
    X_k = max[phi_1 X_{k-1}, \ldots, phi_p X_{k-p}, epsilon_k,
    theta_1 epsilon_{k-1}, \ldots, theta_q epsilon_{k-q}]}
  where \eqn{\code{phi} = (\phi_1, \ldots, \phi_p)}{
    \code{phi} = (phi_1, \ldots, phi_p)}
  and \eqn{\code{theta} = (\theta_1, \ldots, \theta_q)}{
    \code{theta} = (theta_1, \ldots, theta_q)}
  are non-negative vectors of parameters, and where
  \eqn{\{\epsilon_k\}}{{epsilon_k}} is a series of \emph{iid}
  random variables with a common distribution defined by
  \code{rand.gen}.

  The functions \code{mar} and \code{mma} generate MAR(p) and
  MMA(q) processes respectively.
  A MAR(p) process \eqn{\{X_k\}}{{X_k}} is equivalent to a
  MARMA(p, 0) process, so that
  \deqn{X_k = \max\{\phi_1 X_{k-1}, \ldots, \phi_p X_{k-p},
    \epsilon_k\}.}{X_k = max[phi_1 X_{k-1}, \ldots, phi_p X_{k-p},
    epsilon_k].}
  A MMA(q) process \eqn{\{X_k\}}{{X_k}} is equivalent to a
  MARMA(0, q) process, so that
  \deqn{X_k = \max\{\epsilon_k, \theta_1 \epsilon_{k-1}, \ldots,
    \theta_q \epsilon_{k-q}\}.}{X_k = max[epsilon_k, theta_1
    epsilon_{k-1}, \ldots, theta_q epsilon_{k-q}].}
}
\value{
  A numeric vector of length \code{n}.
}

\seealso{\code{\link{evmc}}}

\examples{
marma(100, p = 1, q = 1, psi = 0.75, theta = 0.65)
mar(100, psi = 0.85, n.start = 20)
mma(100, q = 2, theta = c(0.75, 0.8))
}
\keyword{distribution}


\eof
\name{mrlplot}
\alias{mrlplot}
\title{Empirical Mean Residual Life Plot}
\description{
  The empirical mean residual life plot. 
}
\usage{
mrlplot(data, tlim, nt = max(100, length(data)), lty = c(2,1,2),
    col = 1, conf = 0.95, main = "Mean Residual Life Plot", xlab =
    "Threshold", ylab = "Mean Excess", \dots)
}
\arguments{
  \item{data}{A numeric vector.}
  \item{tlim}{A numeric vector of length two, giving the limits for
    the thresholds at which the mean residual life plot is
    evaluated. If \code{tlim} is not given, sensible defaults
    are used.}
  \item{nt}{The number of thresholds at which the mean residual life
    plot is evaluated.}
  \item{lty, col}{Arguments passed to \code{matplot}. The first
    and last elements of \code{lty} correspond to the lower and
    upper confidence limits respectively. Use zero to supress.}
  \item{conf}{The (pointwise) confidence coefficient for the
    plotted confidence intervals.} 
  \item{main}{Plot title.}     
  \item{xlab, ylab}{x and y axis labels.}
  \item{\dots}{Other arguments to be passed to \code{matplot}.}
}
\details{
  The empirical mean residual life plot is the locus of points
  \deqn{\left(u,\frac{1}{n_u} \sum\nolimits_{i=1}^{n_u}
    (x_{(i)} - u) \right)}{{u,1/n_u \sum_{i=1}^{n_u} (x(i) - u)}}
  where \eqn{x_{(1)}, \dots, x_{(n_u)}}{x(1), \dots, x(n_u)} are
  the \eqn{n_u} observations that exceed the threshold \eqn{u}.
  If the exceedances of a threshold \eqn{u_0}{u0}
  are generalized Pareto, the empirical mean residual life plot
  should be approximately linear for \eqn{u > u_0}{u > u0}.
  
  The confidence intervals within the plot are symmetric intervals
  based on the approximate normality of sample means.
}
\value{
  A list with components \code{x} and \code{y} is invisibly returned.
  The components contain those objects that were passed to the formal
  arguments \code{x} and \code{y} of \code{matplot} in order to create
  the mean residual life plot.
}
\author{Stuart Coles and Alec Stephenson}  

\seealso{\code{\link{fpot}}, \code{\link{matplot}},
    \code{\link{tcplot}}}

\examples{
data(portpirie)
mrlplot(portpirie)
}
\keyword{hplot}


\eof
\name{mvevd}
\alias{pmvevd}
\alias{rmvevd}
\alias{dmvevd}
\title{Parametric Multivariate Extreme Value Distributions}
\description{
  Density function, distribution function and random generation
  for the multivariate logistic and multivariate asymmetric
  logistic models.
}
\usage{
pmvevd(q, dep, asy, model = c("log", "alog"), d = 2, mar = c(0,1,0),
    lower.tail = TRUE)
rmvevd(n, dep, asy, model = c("log", "alog"), d = 2, mar = c(0,1,0))
dmvevd(x, dep, asy, model = c("log", "alog"), d = 2, mar = c(0,1,0),
    log = FALSE) 
}
\arguments{
  \item{x, q}{A vector of length \code{d} or a matrix with \code{d}
    columns, in which case the density/distribution is evaluated
    across the rows.}
  \item{n}{Number of observations.}
  \item{dep}{The dependence parameter(s). For the logistic model,
    should be a single value. For the asymmetric logistic model,
    should be a vector of length \eqn{2^d-d-1}, or a single value,
    in which case the value is used for each of the \eqn{2^d-d-1}
    parameters (see \bold{Details}).}
  \item{asy}{The asymmetry parameters for the asymmetric logistic
    model. Should be a list with \eqn{2^d-1} vector elements
    containing the asymmetry parameters for each separate component
    (see \bold{Details}).}
  \item{model}{The specified model; a character string. Must be either
    \code{"log"} (the default) or \code{"alog"} (or any unique partial
    match), for the logistic and asymmetric logistic models respectively.}
  \item{d}{The dimension.}
  \item{mar}{A vector of length three containing marginal parameters
    for every univariate margin, or a matrix with three columns where
    each column represents a vector of values to be passed to the
    corresponding marginal parameter. It can also be a list with
    \code{d} elements, such that each element is either a vector of
    length three or a matrix with three columns, in which case the
    \eqn{i}th element represents the marginal parameters on the
    \eqn{i}th margin.}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x]}
}
\details{
  Define
  \deqn{y_i = y_i(z_i) = \{1+s_i(z_i-a_i)/b_i\}^{-1/s_i}}{
    yi = yi(zi) = {1+si(zi-ai)/bi}^(-1/si)}
  for \eqn{1+s_i(z_i-a_i)/b_i > 0}{1+si(zi-ai)/bi > 0} and
  \eqn{i = 1,\ldots,d}{i = 1,\ldots,d}, where the marginal
  parameters are given by
  \eqn{(a_i,b_i,s_i)}{(ai,bi,si)}, \eqn{b_i > 0}{bi > 0}.
  If \eqn{s_i = 0}{si = 0} then \eqn{y_i}{yi} is defined by
  continuity.
  Let \eqn{z = (z_1,z_2,\ldots,z_d)}{z = (z1,z2,\ldots,zd)}.
  In each of the multivariate distributions functions
  \eqn{G(z)} given below, the
  univariate margins are generalized extreme value, so that
  \eqn{G(z_i) = \exp(-y_i)}{G(zi) = \exp(-yi)} for
  \eqn{i = 1,\ldots,d}{i = 1,\ldots,d}.
  If \eqn{1+s_i(z_i-a_i)/b_i \leq 0}{1+si(zi-ai)/bi <= 0} for some
  \eqn{i = 1,\ldots,d}{i = 1,\ldots,d}, the value \eqn{z_i}{zi} is
  either greater than the upper end point (if \eqn{s_i < 0}{si < 0}),
  or less than the lower end point (if \eqn{s_i > 0}{si > 0}), of the
  \eqn{i}th univariate marginal distribution.
  
  \code{model = "log"} (Gumbel, 1960)
  
  The \code{d} dimensional multivariate logistic distribution
  function with parameter \eqn{\code{dep} = r} is
  \deqn{G(z) = \exp\left\{-\left(\sum\nolimits_{i = 1}^{d}
    y_i^{1/r}\right)^r\right\}}{
    G(z) = exp{-[sum_{i=1}^d yi^(1/r)]^r}}
  where \eqn{0 < r \leq 1}{0 < r <= 1}.
  This is a special case of the multivariate asymmetric logistic
  model.

  \code{model = "alog"} (Tawn, 1990)
  
  Let \eqn{B} be the set of all non-empty subsets of
  \eqn{\{1,\ldots,d\}}{{1,\ldots,d}}, let
  \eqn{B_1=\{b \in B:|b|=1\}}{B1={b in B:|b|=1}}, where \eqn{|b|}
  denotes the number of elements in the set \eqn{b}, and let
  \eqn{B_{(i)}=\{b \in B:i \in b\}}{B(i)={b in B:i in b}}.
  The \code{d} dimensional multivariate asymmetric logistic distribution
  function is
  \deqn{G(z)=\exp\left\{-\sum\nolimits_{b \in B} \left[\sum\nolimits_
    {i\in b}(t_{i,b}y_i)^{1/r_b}\right]^{r_b}\right\},}{G(z) =
    exp{-sum{b in B} [sum{i in b}(t{i,b}yi)^(1/r{b})]^r{b}},}
  where the dependence parameters \eqn{r_b\in(0,1]}{r{b} in (0,1]} for
  all \eqn{b\in B \setminus B_1}{b in B\B1}, and the asymmetry parameters
  \eqn{t_{i,b}\in[0,1]}{t{i,b} in [0,1]} for all \eqn{b\in B}{b in B}
  and \eqn{i\in b}{i in b}.
  The constraints
  \eqn{\sum_{b \in B_{(i)}}t_{i,b}=1}{sum{b in B(i)} t{i,b}=1} for
  \eqn{i = 1,\ldots,d}
  ensure that the marginal distributions are generalized extreme value.
  Further constraints arise from the possible redundancy of asymmetry
  parameters in the expansion of the distribution form.
  Let \eqn{b_{-i_0} = \{i \in b:i \neq i_0\}}{
    b_{-i0} = {i in b:i is not i_0}}.
  If \eqn{r_b = 1}{r{b} = 1} for some
  \eqn{b\in B \setminus B_1}{b in B\B1} then
  \eqn{t_{i,b} = 0}{t{i,b} = 0} for all \eqn{i\in b}{i in b}.
  Furthermore, if for some
  \eqn{b\in B \setminus B_1}{b in B\B1},
  \eqn{t_{i,b} = 0}{t{i,b} = 0} for all
  \eqn{i\in b_{-i_0}}{i in b_{-i0}}, then
  \eqn{t_{i_0,b} = 0}{t{i0,b} = 0}.

  \code{dep} should be a vector of length \eqn{2^d-d-1} which contains
  \eqn{\{r_b:b\in B \setminus B_1\}}{{r{b}:b in B\B1}}, with
  the order defined by the natural set ordering on the index.
  For example, for the trivariate model,
  \eqn{\code{dep} = (r_{12},r_{13},r_{23},r_{123})}{
    \code{dep} = (r{12},r{13},r{23},r{123})}.
  \code{asy} should be a list with \eqn{2^d-1} elements.
  Each element is a vector which corresponds to a set
  \eqn{b\in B}{b in B}, containing \eqn{t_{i,b}}{t{i,b}} for
  every integer \eqn{i\in b}{i in b}.
  The elements should be given using the natural set ordering on the
  \eqn{b\in B}{b in B}, so that the first \eqn{d} elements are vectors
  of length one corresponding to the sets
  \eqn{\{1\},\ldots,\{d\}}{{1},\ldots,{d}}, and the last element is a
  a vector of length \eqn{d}, corresponding to the set
  \eqn{\{1,\ldots,d\}}{{1,\ldots,d}}.
  \code{asy} must be constructed to ensure that all constraints are
  satisfied or an error will occur.
}
\value{
  \code{pmvevd} gives the distribution function, \code{dmvevd} gives
  the density function and \code{rmvevd} generates random deviates, for
  the multivariate logistic or multivariate asymmetric logistic model.
}
\note{
  Multivariate extensions of other bivariate models are more complex.
  A multivariate extension of the Husler-Reiss model exists, involving
  a multidimensional integral and one parameter for each bivariate
  margin.
  Multivariate extensions for the negative logistic model can be derived
  but are considerably more complex and appear to be less flexible.
  The ``multivariate negative logistic model'' often presented in the
  literature (e.g. Kotz \emph{et al}, 2000) is not a valid distribution
  function and should not be used.

  The logistic and asymmetric logistic models respectively are
  simulated using Algorithms 2.1 and 2.2 in Stephenson(2003b).

  The density function of the logistic model is evaluated using
  the representation of Shi(1995).
  The density function of the asymmetric logistic model is evaluated
  using the representation given in Stephenson(2003a).
}
\references{
  Gumbel, E. J. (1960)
  Distributions des valeurs extremes en plusieurs dimensions.
  \emph{Publ. Inst. Statist. Univ. Paris}, \bold{9}, 171--173.
  
  Kotz, S. and Balakrishnan, N. and Johnson, N. L. (2000)
  \emph{Continuous Multivariate Distributions}, vol. 1.
  New York: John Wiley & Sons, 2nd edn.

  Shi, D. (1995)
  Fisher information for a multivariate extreme value distribution.
  \emph{Biometrika}, \bold{82}(3), 644--649.

  Stephenson, A. G. (2003a)
  \emph{Extreme Value Distributions and their Application}.
  Ph.D. Thesis, Lancaster University, Lancaster, UK.
  
  Stephenson, A. G. (2003b)
  Simulating multivariate extreme value distributions of logistic type.
  \emph{Extremes}, \bold{6}(1), 49--60.

  Tawn, J. A. (1990)
  Modelling multivariate extreme value distributions.
  \emph{Biometrika}, \bold{77}, 245--253.
}


\seealso{\code{\link{rbvevd}}, \code{\link{rgev}}}

\examples{
pmvevd(matrix(rep(0:4,5), ncol=5), dep = .7, model = "log", d = 5)
pmvevd(rep(4,5), dep = .7, model = "log", d = 5)
rmvevd(10, dep = .7, model = "log", d = 5)
dmvevd(rep(-1,20), dep = .7, model = "log", d = 20, log = TRUE)

asy <- list(.4, .1, .6, c(.3,.2), c(.1,.1), c(.4,.1), c(.2,.3,.2))
pmvevd(rep(2,3), dep = c(.6,.5,.8,.3), asy = asy, model = "alog", d = 3)
asy <- list(.4, .0, .6, c(.3,.2), c(.1,.1), c(.4,.1), c(.2,.4,.2))
rmvevd(10, dep = c(.6,.5,.8,.3), asy = asy, model = "alog", d = 3)
dmvevd(rep(0,3), dep = c(.6,.5,.8,.3), asy = asy, model = "alog", d = 3)

asy <- list(0, 0, 0, 0, c(0,0), c(0,0), c(0,0), c(0,0), c(0,0), c(0,0),
  c(.2,.1,.2), c(.1,.1,.2), c(.3,.4,.1), c(.2,.2,.2), c(.4,.6,.2,.5))
rmvevd(10, dep = .7, asy = asy, model = "alog", d = 4)
rmvevd(10, dep = c(rep(1,6), rep(.7,5)), asy = asy, model = "alog", d = 4)
}
\keyword{distribution}


\eof
\name{ocmulgee}
\alias{ocmulgee}
\title{Maximum Annual Flood Discharges of the Ocmulgee River}
\usage{data(ocmulgee)}
\description{
  The \code{ocmulgee} data frame has 40 rows and 2 columns.
  The columns contain maximum annual flood discharges, in units
  of 1000 cubed feet per second, from the Ocmulgee River in
  Georgia, USA at Hawkinsville (upstream) and Macon (downstream),
  for the years 1910 to 1949.
  The row names give the years of observation.
}
\format{
  This data frame contains the following columns:
  \describe{
    \item{hawk}{A numeric vector containing maximum annual flood
      discharges at Hawkinsville (upstream).}
    \item{macon}{A numeric vector containing maximum annual flood
      discharges at Macon (downstream).}
  }
}
\source{
  Gumbel, E. J. and Goldstein, N. (1964)
  Analysis of empirical bivariate extremal distributions.
  \emph{J. Amer. Statist. Assoc.}, \bold{59}, 794--816.
}
\keyword{datasets}

\eof
\name{oldage}
\alias{oldage}
\title{Oldest Ages for Swedish Males and Females}
\usage{data(oldage)}
\description{
  The \code{oldage} data frame has 66 rows and 2 columns.
  The columns contain the oldest ages at death for men and women in
  Sweden, for the period 1905--1970.
  The row names give the years of observation.
}
\format{
  This data frame contains the following columns:
  \describe{
    \item{men}{A numeric vector containing the oldest ages at death
      for men.}
    \item{women}{A numeric vector containing the oldest ages at death
      for women.}
  }
}
\source{
  Fransen, A. and Tiago de Oliveira, J. (1984)
  Statistical choice of univariate extreme models, part II, in
  \emph{Statistical Extremes and Applications},
  J. Tiago de Oliveira ed., 373--394, D. Reidel, Dordrect.
}
\keyword{datasets}

\eof
\name{order}
\alias{dorder}
\alias{porder}
\alias{rorder}
\title{Distributions of Order Statistics}
\description{
  Density function, distribution function and random generation for
  a selected order statistic of a given number of independent variables
  from a specified distribution. 
}
\usage{
dorder(x, densfun, distnfun, \dots, distn, mlen = 1, j = 1,
      largest = TRUE, log = FALSE)
porder(q, distnfun, \dots, distn, mlen = 1, j = 1, largest = TRUE,
      lower.tail = TRUE)  
rorder(n, quantfun, \dots, distn, mlen = 1, j = 1, largest = TRUE) 
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{n}{Number of observations.}
  \item{densfun, distnfun, quantfun}{Density, distribution and
    quantile function of the specified distribution. The density
    function must have a \code{log} argument (a simple wrapper
    can always be constructed to achieve this).}
  \item{\dots}{Parameters of the specified distribution.}
  \item{distn}{A character string, optionally specified as an
    alternative to \code{densfun}, \code{distnfun} and \code{quantfun}
    such that the density, distribution and quantile functions are
    formed upon the addition of the prefixes \code{d}, \code{p} and
    \code{q} respectively.}
  \item{mlen}{The number of independent variables.}
  \item{j}{The order statistic, taken as the \code{j}th largest
    (default) or smallest of \code{mlen}, according to the value of
    \code{largest}.}
  \item{largest}{Logical; if \code{TRUE} (default) use the \code{j}th
    largest order statistic, otherwise use the \code{j}th smallest.}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default) probabilities
    are P[X <= x], otherwise  P[X > x].}
}
\value{
  \code{dorder} gives the density function, \code{porder} gives the
  distribution function and \code{qorder} gives the quantile function
  of a selected order statistic from a sample of size \code{mlen},
  from a specified distibution. \code{rorder} generates random deviates.
}

\seealso{\code{\link{rextreme}}, \code{\link{rgev}}}

\examples{
dorder(2:4, dnorm, pnorm, mean = 0.5, sd = 1.2, mlen = 5, j = 2)
dorder(2:4, distn = "norm", mean = 0.5, sd = 1.2, mlen = 5, j = 2)
dorder(2:4, distn = "exp", mlen = 2, j = 2)
porder(2:4, distn = "exp", rate = 1.2, mlen = 2, j = 2)
rorder(5, qgamma, shape = 1, mlen = 10, j = 2)
}
\keyword{distribution}


\eof
\name{oxford}
\alias{oxford}
\title{Annual Maximum Temperatures at Oxford}
\usage{data(oxford)}
\description{
  A numeric vector containing annual maximum temperatures, in degrees
  Fahrenheit, from 1901 to 1980 at Oxford, England.
}
\format{A vector containing 80 observations.}
\source{
  Tabony, R. C. (1983)
  Extreme value analysis in meteorology.
  \emph{The Meteorological Magazine} \bold{112}, 77--98.
}  
\keyword{datasets}

\eof
\name{plot.bvevd}
\alias{plot.bvevd}

\title{Plot Diagnostics for a Bivariate EVD Object}
\description{
  Four plots (selectable by \code{which}) are currently provided:
  two conditional P-P plots (conditioning on each margin),
  a density plot and a dependence function plot.
  Plot diagnostics for the generalized extreme value margins
  (selectable by \code{mar} and \code{which}) are also available.
}
\usage{
\method{plot}{bvevd}(x, mar = 0, which = 1:4, main = c("Conditional Plot One", 
    "Conditional Plot Two", "Density Plot", "Dependence Function"), 
    ask = nb.fig < length(which) && dev.interactive(), ci = TRUE, 
    jitter = FALSE, grid = 50, nplty = 2, blty = 3, method = "cfg",
    convex = FALSE, wf = function(t) t, \dots)
}
\arguments{
  \item{x}{An object of class \code{"bvevd"}.}
  \item{mar}{If \code{mar = 1} or \code{mar = 2} diagnostics
    are given for the first or second genereralized extreme
    value margin respectively. The values of the remaining
    parameters are then passed to the plot method
    \code{\link{plot.uvevd}}.}
  \item{which}{If a subset of the plots is required, specify a
    subset of the numbers \code{1:4}.}
  \item{main}{Title of each plot.}
  \item{ask}{Logical; if \code{TRUE}, the user is asked before
    each plot.}
  \item{ci}{Logical; if \code{TRUE} (the default), plot simulated
    95\% confidence intervals for the conditional P-P plots.}
  \item{jitter, grid}{Arguments for the density plot. The (possibly
    transformed) data is plotted with a contour plot of the
    bivariate density of the fitted model. The density is evaluated
    at \code{grid^2} points. If \code{jitter} is \code{TRUE}, the
    data are jittered. This need only be used if the data contains
    repeated values.}
  \item{nplty, blty, method, convex, wf}{Arguments to the dependence
    function plot. The dependence function for the fitted model is
    plotted and (optionally) compared to a non-parameteric estimate.
    See \code{\link{abvnonpar}} for a definition of the dependence
    function, and for a description of the arguments \code{method},
    \code{modify} and \code{wf}, which alter the behaviour of the
    non-parametric estimator. \code{nplty} is the line type of the
    non-parametric estimate. To omit the non-parametric estimate set
    \code{nplty} to zero. \code{blty} is the line type of the
    triangular border. To omit the border estimate set \code{blty}
    to zero.}
  \item{\dots}{Other arguments to be passed through to plotting
    functions.}
}
\details{
  The following discussion assumes that the fitted model is
  stationary. For non-stationary models the data are transformed
  to stationarity. The plot then corresponds to the distribution
  obtained when all covariates are zero.
  
  A conditional P-P plot is a P-P plot for the condition
  distribution function of a bivariate evd object.
  Let \eqn{G(.|.)} be the conditional distribution of
  the first margin given the second, under the fitted model.
  Let \eqn{z_1,\ldots,z_m} be the data used in the fitted model,
  where \eqn{z_j = (z_{1j}, z_{2j})} for \eqn{j = 1,\ldots,m}.
  The plot that (by default) is labelled Conditional Plot Two,
  conditioning on the second margin, consists of the points
  \deqn{\{(p_i, c_i), i = 1,\ldots,m\}}{
    {(p_i, c_i), i = 1,\ldots,m}}
  where \eqn{p_1,\ldots,p_m} are plotting points defined by
  \code{\link{ppoints}} and \eqn{c_i} is the \eqn{i}th largest
  value from the sample
  \eqn{\{G(z_{j1}|z_{j2}), j = 1,\ldots,m\}.}{
    {G(z_{j1}|z_{j2}), j = 1,\ldots,m}.}
  The margins are reversed for Conditional Plot One, so that
  \eqn{G(.|.)} is the conditional distribution of the second
  margin given the first. 
}

\seealso{\code{\link{plot.uvevd}}, \code{\link{contour}},
  \code{\link{jitter}}, \code{\link{abvnonpar}}}

\examples{
bvdata <- rbvevd(100, dep = 0.6, model = "log")
M1 <- fbvevd(bvdata, model = "log")
\dontrun{par(mfrow = c(2,2))}
\dontrun{plot(M1)}
\dontrun{plot(M1, mar = 1)}
\dontrun{plot(M1, mar = 2)}
}
\keyword{hplot}

\eof
\name{plot.profile.evd}
\alias{plot.profile.evd}

\title{Plot Profile Deviances and Calculate Profile Confidence
  Intervals}
\description{
  Displays profile deviances from a model profiled with
  \code{\link{profile.evd}}, and invisibly returns profile
  confidence intervals.
}
\usage{
\method{plot}{profile.evd}(x, which = names(x), main = NULL,
    ask = nb.fig < length(which) && dev.interactive(), ci = 0.95,
    clty = 2, \dots)
}
\arguments{
  \item{x}{An object of class \code{"profile.evd"}.}
  \item{which}{A character vector giving the parameters for which
    the profile deviance is plotted, and for which profile confidence
    intervals are calculated. By default all profiled parameters in
    \code{x} are used.}
  \item{main}{Title of each plot; a character vector, the
    same length as \code{which}.}
  \item{ask}{Logical; if \code{TRUE}, the user is asked before
    each plot.}
  \item{ci}{A numeric vector. For each parameter in \code{which}
    profile confidence intervals are calculated, for each confidence
    coefficient in \code{ci} (but see \bold{Warning}).
    The intervals are returned invisibly as a list of
    vectors/matrices. Each plot then (by default) includes
    horizonal lines that represent each interval.}
  \item{clty}{The line type of the horizontal lines that represent
    the profile confidence intervals. To omit the lines set
    \code{clty} to zero.}
  \item{\dots}{Other graphics parameters.}
}
\value{
  Profile devainces are plotted for each parameter in
  \code{which}.
  A list with one element for each parameter in \code{which}
  is also returned invisibly.
  Each element is a vector of length two or a matrix with two
  columns, giving the lower and upper limits of each confidence
  interval.
}
\section{Warning}{
  The profile confidence intervals may not have confidence coefficient
  \code{ci}, because the usual asymptotic properties of maximum
  likelihood estimators may not hold.
  For the GEV model, the usual asymptotic properties hold when the
  shape parameter is greater than \eqn{-0.5} (Smith, 1985).
  Fortunately, this is usually the case.
}
\references{
  Smith, R. L. (1985)
  Maximum likelihood estimation in a class of non-regular cases.
  \emph{Biometrika}, \bold{72}, 67--90.
}

\seealso{\code{\link{plot.profile2d.evd}}, \code{\link{profile.evd}},
  \code{\link{profile2d.evd}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata)
\dontrun{M1P <- profile(M1)}
\dontrun{par(mfrow = c(2,2))}
\dontrun{cint <- plot(M1P, ci = c(0.95, 0.99))}
\dontrun{cint}
}
\keyword{hplot}

\eof
\name{plot.profile2d.evd}
\alias{plot.profile2d.evd}

\title{Plot Joint Profile Devainces}
\description{
  Displays an image plot of the joint profile deviance from a
  model profiled with \code{\link{profile.evd}} and
  \code{\link{profile2d.evd}}.
}
\usage{
\method{plot}{profile2d.evd}(x, main = NULL,
    ci = c(0.5, 0.8, 0.9, 0.95, 0.975, 0.99, 0.995),
    col = heat.colors(8), intpts = 75, xaxs = "r", yaxs = "r", \dots)
}
\arguments{
  \item{x}{An object of class \code{"profile2d.evd"}.}
  \item{main}{Title of plot; a character string.}
  \item{ci}{A numeric vector whose length is one less than the
    length of \code{col}. The colours of the image plot,
    excluding the background colour, represent confidence sets
    with confidence coefficients \code{ci} (but see \bold{Warning}).}
  \item{col}{A list of colors such as that generated by
    \code{rainbow}, \code{heat.colors}, \code{topo.colors},
    \code{terrain.colors} or similar functions.}
  \item{intpts}{If the package \bold{akima} is available,
    interpolation is performed using \code{intpts} points
    for each parameter. The function is interpolated at
    \code{intpts^2} points in total.}
  \item{xaxs,yaxs}{Graphics parameters (see \code{\link{par}}).
    The default, \code{"r"}, overrides the default set by
    \code{image}.}
  \item{\dots}{Other parameters to be passed to \code{image}.}
}
\section{Warning}{
  The sets represented by different colours may not be
  confidence sets with confidence coefficients \code{ci}, because
  the usual asymptotic properties of maximum likelihood estimators
  may not hold.
  For the GEV model, the usual asymptotic properties hold when the
  shape parameter is greater than \eqn{-0.5} (Smith, 1985).
  Fortunately, this is usually the case.
}
\references{
  Smith, R. L. (1985)
  Maximum likelihood estimation in a class of non-regular cases.
  \emph{Biometrika}, \bold{72}, 67--90.
}

\seealso{\code{\link{plot.profile.evd}}, \code{\link{profile.evd}},
  \code{\link{profile2d.evd}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata)
\dontrun{M1P <- profile(M1)}
\dontrun{M1JP <- profile2d(M1, M1P, which = c("scale", "shape"))}
\dontrun{plot(M1JP)}
}
\keyword{hplot}

\eof
\name{plot.uvevd}
\alias{plot.uvevd}

\title{Plot Diagnostics for a Univariate EVD Object}
\description{
  Four plots (selectable by \code{which}) are currently provided:
  a P-P plot, a Q-Q plot, a density plot and a return level plot.
}
\usage{
\method{plot}{uvevd}(x, which = 1:4, main = c("Probability Plot",
    "Quantile Plot", "Density Plot", "Return Level Plot"),
    ask = nb.fig < length(which) && dev.interactive(),
    ci = TRUE, adjust = 1, jitter = FALSE, nplty = 2, \dots)
}
\arguments{
  \item{x}{An object that inherits from class \code{"uvevd"}.}
  \item{which}{If a subset of the plots is required, specify a
    subset of the numbers \code{1:4}.}
  \item{main}{Title of each plot.}
  \item{ask}{Logical; if \code{TRUE}, the user is asked before
    each plot.}
  \item{ci}{Logical; if \code{TRUE} (the default), plot simulated
    95\% confidence intervals for the P-P, Q-Q and return level
    plots.}
  \item{adjust, jitter, nplty}{Arguments to the density plot.
    The density of the fitted model is plotted with a rug plot and
    (optionally) a non-parameteric estimate. The argument
    \code{adjust} controls the smoothing bandwidth for the
    non-parametric estimate (see \code{\link{density}}).
    \code{jitter} is logical; if \code{TRUE}, the (possibly
    transformed) data are jittered to produce the rug plot.
    This need only be used if the data contains repeated
    values. \code{nplty} is the line type of the non-parametric
    estimate. To omit the non-parametric estimate set \code{nplty}
    to zero.}
  \item{\dots}{Other parameters to be passed through to plotting
    functions.}
}
\details{
  The following discussion assumes that the fitted model is
  stationary. For non-stationary generalized extreme value
  models the data are transformed to stationarity. The plot
  then corresponds to the distribution obtained when all
  covariates are zero.

  The P-P plot consists of the points
  \deqn{\{(G_n(z_i), G(z_i)), i = 1,\ldots,m\}}{
    {(G_n(z_i), G(z_i)), i = 1,\ldots,m}}
  where \eqn{G_n} is the empirical distribution function
  (defined using \code{\link{ppoints}}), G is the model based
  estimate of the distribution (generalized extreme value
  or generalized Pareto), and \eqn{z_1,\ldots,z_m} are the data
  used in the fitted model, sorted into ascending order.

  The Q-Q plot consists of the points
  \deqn{\{(G^{-1}(p_i), z_i), i = 1,\ldots,m\}}{
    {(G^{-1}(p_i), z_i), i = 1,\ldots,m}}
  where \eqn{G^{-1}} is the model based estimate of the quantile
  function (generalized extreme value or generalized Pareto),
  \eqn{p_1,\ldots,p_m} are plotting points defined by
  \code{\link{ppoints}}, and \eqn{z_1,\ldots,z_m} are the data
  used in the fitted model, sorted into ascending order.

  The return level plot for generalized extreme value models
  is defined as follows.
  
  Let \eqn{G} be the generalized extreme value distribution
  function, with location, scale and shape parameters \eqn{a},
  \eqn{b} and \eqn{s} respectively.
  Let \eqn{z_t} be defined by \eqn{G(z_t) = 1 - 1/t}.
  In common terminology, \eqn{z_t} is the return level
  associated with the return period \eqn{t}.

  Let \eqn{y_t = -1/\log(1 - 1/t)}{y_t = -1/log(1 - 1/t)}.
  It follows that
  \deqn{z_t = a + b(y_t^s - 1)/s.}{
    z_t = a + b((y_t)^s - 1)/s.}
  When \eqn{s = 0}, \eqn{z_t} is defined by continuity, so that
  \deqn{z_t = a + b\log(y_t).}{
    z_t = a + b log(y_t).}
  The curve within the return level plot is \eqn{z_t} plotted
  against \eqn{y_t} on a logarithmic scale, using maximum likelihood
  estimates of \eqn{(a,b,s)}. If the estimate of \eqn{s} is zero, the
  curve will be linear.
  For large values of \eqn{t}, \eqn{y_t} is approximately equal
  to the return period \eqn{t}. It is usual practice to label the
  x-axis as the return period.  

  The points on the plot are  
  \deqn{\{(-1/\log(p_i), z_i), i = 1,\ldots,m\}}{
    {(-1/log(p_i), z_i), i = 1,\ldots,m}}
  where \eqn{p_1,\ldots,p_m} are plotting points defined by
  \code{\link{ppoints}}, and \eqn{z_1,\ldots,z_m} are the data
  used in the fitted model, sorted into ascending order.
  For a good fit the points should lie ``close'' to the curve.

  The return level plot for peaks over threshold models is
  defined as follows.
  
  Let \eqn{G} be the generalized Pareto distribution function,
  with location, scale and shape parameters \eqn{u}, \eqn{b}
  and \eqn{s} respectively, where \eqn{u} is the model threshold.
  Let \eqn{z_m} denote the \eqn{m} period return level
  (see \code{\link{fpot}} and the notation therein).
  It follows that 
  \deqn{z_m = u + b((pmN)^s - 1)/s.}{
    z_m = u + b((pmN)^s - 1)/s.}
  When \eqn{s = 0}, \eqn{z_m} is defined by continuity, so that
  \deqn{z_m = u + b\log(pmN).}{
    z_m = u + b log(pmN).}
  The curve within the return level plot is \eqn{z_m} plotted
  against \eqn{m} on a logarithmic scale, using maximum likelihood
  estimates of \eqn{(b,s,p)}. If the estimate of \eqn{s} is zero,
  the curve will be linear.  

  The points on the plot are  
  \deqn{\{(1/(pN(1-p_i)), z_i), i = 1,\ldots,m\}}{
    {(1/(pN(1-p_i)), z_i), i = 1,\ldots,m}}
  where \eqn{p_1,\ldots,p_m} are plotting points defined by
  \code{\link{ppoints}}, and \eqn{z_1,\ldots,z_m} are the data
  used in the fitted model, sorted into ascending order.
  For a good fit the points should lie ``close'' to the curve.
}

\seealso{\code{\link{plot.bvevd}}, \code{\link{density}},
  \code{\link{jitter}}, \code{\link{rug}}, \code{\link{ppoints}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata)
\dontrun{par(mfrow = c(2,2))}
\dontrun{plot(M1)}

uvdata <- rgpd(100, loc = 0, scale = 1.1, shape = 0.2)
M1 <- fpot(uvdata, 1)
\dontrun{par(mfrow = c(2,2))}
\dontrun{plot(M1)}
}
\keyword{hplot}

\eof
\name{portpirie}
\alias{portpirie}
\title{Annual Maximum Sea Levels at Port Pirie}
\usage{data(portpirie)}
\description{
  A numeric vector containing annual maximum sea levels, in metres,
  from 1923 to 1987 at Port Pirie, South Australia.
}
\format{A vector containing 65 observations.}
\source{
  Tawn, J. A. (1993)
  Extreme sea-levels, in \emph{Statistics in the Environment},
  243--263, eds. V. Barnett and F. Turkman, Wiley.
}
\references{
  Coles, S. G. (2001)
  \emph{An Introduction to Statistical Modeling of Extreme Values}.
  London: Springer-Verlag.
}
  
\keyword{datasets}

\eof
\name{profile.evd}
\alias{profile.evd}

\title{Method for Profiling EVD Objects}
\description{
  Calculate profile traces for fitted models.
}
\usage{
\method{profile}{evd}(fitted, which = names(fitted$estimate), conf = 0.999,
    mesh = fitted$std.err[which]/4, xmin = rep(-Inf, length(which)),
    xmax = rep(Inf, length(which)), convergence = FALSE, method = "BFGS",
    control = list(maxit = 500), \dots)
}
\arguments{
  \item{fitted}{An object of class \code{"evd"}.}
  \item{which}{A character vector giving the model parameters that
    are to be profiled. By default, all parameters are profiled.}
  \item{conf}{Controls the range over which the parameters are profiled.
    The profile trace is constructed so that (assuming the usual
    asymptotic properties hold) profile confidence intervals with
    confidence coefficients \code{conf} or less can be derived from it.}
  \item{mesh}{A numeric vector containing one value for each
    parameter in \code{which}. The values represent the
    distance between the points profiled. By default \code{mesh} is
    one quarter of the standard errors. If the fitted object does not
    contain standard errors the argument must be specified.
    The argument should also be specified when an estimator is
    on or close to a parameter boundary, since the approximated
    ``standard error'' will then be close to zero.}
  \item{xmin, xmax}{Numeric vectors containing one value for each
    parameter in \code{which}. Each value represents the theoretical
    lower/upper bound of the corresponding parameter.
    The arguments are needed only when a parameter has a
    lower/upper bound at which the likelihood is non-zero. Do not
    use these arguments to specify plotting ranges in a subsequent
    plot (as they are used in the calculation of profile confidence
    intervals); to do this use \code{xlim} in the call to \code{plot}.}
  \item{convergence}{Logical; print convergence code after each
    optimization? (A warning is given for each non-zero convergence
    code, irrespective of the value of \code{convergence}.)}
  \item{method}{The optimization method.}
  \item{control}{Passed to \code{optim}. See \code{\link{optim}} for
    details.}
  \item{\dots}{Ignored.}
}
\value{
  An object of class \code{"profile.evd"}, which is a list with an
  element for each parameter being profiled. The elements are
  matrices. The first column contains the values of the profiled
  parameter. The second column contains profile deviances. The
  remaining columns contain the constrained maximum likelihood
  estimates for the remaining model parameters.
}

\seealso{\code{\link{profile2d.evd}}, \code{\link{plot.profile.evd}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata)
\dontrun{M1P <- profile(M1)}
\dontrun{par(mfrow = c(2,2))}
\dontrun{cint <- plot(M1P)}
\dontrun{cint}
}
\keyword{models}




\eof
\name{profile2d.evd}
\alias{profile2d}
\alias{profile2d.evd}

\title{Method for Profiling EVD Objects}
\description{
  Calculate joint profile traces for fitted models.
}
\usage{
\method{profile2d}{evd}(fitted, prof, which, pts = 20, convergence =
    FALSE, method = "Nelder-Mead", control = list(maxit = 5000), \dots)
}
\arguments{
  \item{fitted}{An object of class \code{"evd"}.}
  \item{prof}{An object of class \code{"profile.evd"}, created using
    \code{\link{profile.evd}} with argument \code{fitted}.
    The object must contain the (marginal) profile traces for the two
    parameters specified in \code{which}.}
  \item{which}{A character vector of length two containing the original 
     model parameters that are to be jointly profiled.}
   \item{pts}{The number of distinct values used for each profiled
     parameter in \code{which}.
    There are \code{pts^2} optimizations performed in total.}
  \item{convergence}{Logical; print convergence code after each
    optimization? (A warning is given for each non-zero convergence
    code, irrespective of the value of \code{convergence}.)}
  \item{method}{The optimization method.}
  \item{control}{Passed to \code{optim}. See \code{\link{optim}} for
    details.}
  \item{\dots}{Ignored.}
}
\value{
  An object of class \code{"profile2d.evd"}, which is a list with three
  elements.
  The first element, a matrix named \code{trace}, has the same structure
  as the elements of an object of class \code{"profile.evd"}.
  The last two elements give the distinct values used for each profiled
  parameter in \code{which}.
}

\seealso{\code{\link{profile.evd}}, \code{\link{plot.profile2d.evd}}}

\examples{
uvdata <- rgev(100, loc = 0.13, scale = 1.1, shape = 0.2)
M1 <- fgev(uvdata)
\dontrun{M1P <- profile(M1)}
\dontrun{M1JP <- profile2d(M1, M1P, which = c("scale", "shape"))}
\dontrun{plot(M1JP)}
}
\keyword{models}

\eof
\name{rweibull}
\alias{drweibull}
\alias{prweibull}
\alias{qrweibull}
\alias{rrweibull}
\title{The Reversed Weibull Distribution}
\description{
  Density function, distribution function, quantile function and
  random generation for the reversed Weibull distribution with
  location, scale and shape parameters. 
}
\usage{
drweibull(x, loc=0, scale=1, shape=1, log = FALSE) 
prweibull(q, loc=0, scale=1, shape=1, lower.tail = TRUE) 
qrweibull(p, loc=0, scale=1, shape=1, lower.tail = TRUE)
rrweibull(n, loc=0, scale=1, shape=1)
}
\arguments{
  \item{x, q}{Vector of quantiles.}
  \item{p}{Vector of probabilities.}
  \item{n}{Number of observations.}
  \item{loc, scale, shape}{Location, scale and shape parameters (can be
    given as vectors).}
  \item{log}{Logical; if \code{TRUE}, the log density is returned.}
  \item{lower.tail}{Logical; if \code{TRUE} (default), probabilities
    are P[X <= x], otherwise, P[X > x]}
}
\details{
  The reversed Weibull distribution function with parameters
  \eqn{\code{loc} = a}, \eqn{\code{scale} = b} and
  \eqn{\code{shape} = s} is
  \deqn{G(z) = \exp\left\{-\left[-\left(\frac{z-a}{b}\right)
    \right]^s\right\}}{G(x) = exp{-[-(z-a)/b]^s}}
  for \eqn{z < a} and one otherwise, where \eqn{b > 0} and
  \eqn{s > 0}.
}
\note{
  Within extreme value theory the reversed Weibull distibution is
  usually referred to as the Weibull distribution.
  I make a distinction to avoid confusion with the three-parameter
  distribution used in survival analysis, which is related by a
  change of sign to the distribution given above.
}
\value{
  \code{drweibull} gives the density function, \code{prweibull}
  gives the distribution function, \code{qrweibull} gives the
  quantile function, and \code{rrweibull} generates random deviates.
}

\seealso{\code{\link{rfrechet}}, \code{\link{rgev}}, \code{\link{rgumbel}}}

\examples{
drweibull(-5:-3, -1, 0.5, 0.8)
prweibull(-5:-3, -1, 0.5, 0.8)
qrweibull(seq(0.9, 0.6, -0.1), 2, 0.5, 0.8)
rrweibull(6, -1, 0.5, 0.8)
p <- (1:9)/10
prweibull(qrweibull(p, -1, 2, 0.8), -1, 2, 0.8)
## [1] 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
}
\keyword{distribution}





\eof
\name{sask}
\alias{sask}
\title{Maximum Annual Flood Discharges of the North Saskachevan River}
\usage{data(sask)}
\description{
  A numeric vector containing maximum annual flood discharges, in
  units of 1000 cubic feet per second, of the North Saskachevan
  River at Edmonton, over a period of 47 years. Unfortunately,
  the data are ordered from largest to smallest.
}
\format{A vector containing 47 observations.}
\source{
  van Montfort, M. A. J. (1970)
  On testing that the distribution is of type I when type II
  is the alternative.
  \emph{J. Hydrology}, \bold{11}, 421--427.
}

  
\keyword{datasets}

\eof
\name{sealevel}
\alias{sealevel}
\title{Annual Sea Level Maxima at Dover and Harwich}
\usage{data(sealevel)}
\description{
  The \code{sealevel} data frame has 81 rows and 2 columns.
  The columns contain annual sea level maxima from 1912 to 1992 at
  Dover and Harwich respectively, two sites on the coast of Britain.
  The row names give the years of observation.
  There are 39 missing values.
}
\format{
  This data frame contains the following columns:
  \describe{
    \item{dover}{A numeric vector containing annual sea level maxima
    at Dover, including 9 missing values.}
    \item{harwich}{A numeric vector containing sea annual level maxima
    at Harwich, including 30 missing values.}
  }
}
\source{
  Coles, S. G. and Tawn, J. A. (1990)
  Statistics of coastal flood prevention.
  \emph{Phil. Trans. R. Soc. Lond., A} \bold{332}, 457--476.
}
\keyword{datasets}

\eof
\name{tcplot}
\alias{tcplot}
\title{Threshold Choice Plot}
\description{
  Plots of parameter estimates at various thresholds for peaks
  over threshold modelling, using the Generalized Pareto or Point
  Process representation. 
}
\usage{
tcplot(data, tlim, model = c("gpd","pp"), cmax = FALSE, r = 1,
    ulow = -Inf, rlow = 1, nt = 25, which = 1:npar, conf = 0.95,
    lty = 1, lwd = 1, type = "b", cilty = 1, ask = nb.fig <
    length(which) && dev.interactive(), \dots)
}
\arguments{
  \item{data}{A numeric vector.}
  \item{tlim}{A numeric vector of length two, giving the limits for
    the thresholds at which the model is fitted.}
  \item{model}{The model; either \code{"gpd"} (the default) or
    \code{"pp"}, for the Generalized Pareto or Point Process
    representations respectively.}
  \item{cmax}{Logical; if \code{FALSE} (the default), the models
    are fitted using all exceedences over the thresholds. If
    \code{TRUE}, the models are fitted using cluster maxima, using
    clusters of exceedences derived from \code{clusters}.}
  \item{r, ulow, rlow}{Arguments used for the identification of
    clusters of exceedences (see \code{\link{clusters}}).
    Ignored if \code{cmax} is \code{FALSE} (the default).}
  \item{nt}{The number of thresholds at which the model is fitted.}
  \item{which}{If a subset of the plots is required, specify a
    subset of the numbers \code{1:npar}, where \code{npar} is
    the number of parameters, so that \code{npar = 2} when
    \code{model = "gpd"} (the default) and \code{npar = 3}
    when \code{model = "pp"}.}
  \item{conf}{The (pointwise) confidence coefficient for the
    plotted confidence intervals. Use zero to suppress.}
  \item{lty, lwd}{The line type and width of the line connecting
    the parameter estimates.}
  \item{type}{The form taken by the line connecting the parameter
    estimates and the points denoting these estimates. Possible
    values include \code{"b"} (the default) for points joined by
    lines, \code{"o"} for overplotted points and lines, and
    \code{"l"} for an unbroken line with no points.}
  \item{cilty}{The line type of the lines depicting the
    confidence intervals.}
  \item{ask}{Logical; if \code{TRUE}, the user is asked before
    each plot.}
  \item{\dots}{Other arguments to be passed to the model fit
    function \code{fpot}.}
}
\details{
  For each of the \code{nt} thresholds a peaks over threshold model
  is fitted using the function \code{fpot}. When \code{model = "gpd"}
  (the default), the maximum likelihood estimates for the shape and the
  modified scale parameter (modified by subtracting the shape multiplied
  by the threshold) are plotted against the thresholds.
  When \code{model = "pp"} the maximum likelihood estimates for the
  location, scale and shape parameters are plotted against the
  thresholds. (The modified scale parameter in the \code{"gpd"} case
  is equivalent to the scale parameter in the \code{"pp"} case.)
  If the threshold \code{u} is a valid threshold to be used for peaks
  over threshold modelling, the parameter estimates depicted should
  be approximately constant above \code{u}.
}
\value{
  A list is invisibly returned. Each component is a matrix with
  three columns giving parameter estimates and confidence limits.
}
\author{Stuart Coles and Alec Stephenson}

\seealso{\code{\link{fpot}}, \code{\link{mrlplot}},
    \code{\link{clusters}}}

\examples{
data(portpirie)
tlim <- c(3.6, 4.2)
\dontrun{tcplot(portpirie, tlim)}
\dontrun{tcplot(portpirie, tlim, nt = 100, lwd = 3, type = "l")}
\dontrun{tcplot(portpirie, tlim, model = "pp")}
}
\keyword{hplot}


\eof
\name{uccle}
\alias{uccle}
\title{Rainfall Maxima at Uccle, Belguim}
\usage{data(uccle)}
\description{
  The \code{uccle} data frame has 35 rows and 4 columns.
  The columns contain annual rainfall maxima (in milimetres) from
  1938 to 1972 at Uccle, Belguim, over the durations of one day,
  one hour, ten minutes and one minute.
  The row names give the years of observation.
}
\format{
  This data frame contains the following columns:
  \describe{
    \item{day}{Annual daily rainfall maxima.}
    \item{hour}{Annual hourly rainfall maxima.}
    \item{tmin}{Annual rainfall maxima over ten minute durations.}
    \item{min}{Annual rainfall maxima over one minute durations.}
  }
}
\source{
  Sneyers, R. (1977)
  L'intensit\'{e} maximale des pr\'{e}cipitations en Belgique.
  \emph{Inst. Royal M\'{e}t\'{e}or. Belgique, B} \bold{86}.
}
\keyword{datasets}

\eof
\name{venice}
\alias{venice}
\title{Largest Sea Levels in Venice}
\usage{data(venice)}
\description{
  The \code{venice} data frame has 51 rows and 10 columns.
  The jth column contains the jth largest sea levels in Venice,
  for the years 1931--1981.
  Only the largest six measurements are available for the year 1935;
  the corresponding row contains four missing values.
  The years for each set of measurements are given as row names. 
}
\format{A data frame with 51 rows and 10 columns.}
\source{
  Smith, R. L. (1986)
  Extreme value theory based on the \eqn{r} largest annual events.
  \emph{Journal of Hydrology}, \bold{86}, 27--43.
}
\references{
  Coles, S. G. (2001)
  \emph{An Introduction to Statistical Modeling of Extreme Values}.
  London: Springer-Verlag.
}
\keyword{datasets}

\eof
