# R/plot.evm.R In texmex: Statistical Modelling of Extreme Values

#### Documented in plot.evmOpt

#' Plots for evmOpt objects
#'
#' Various plots for \code{evmOpt} objects. These differ depending on
#' whether or not there are covariates in the model.  If there are no
#' covariates then the diagnostic plots are PP- and QQ-plots, a return
#' level plot (produced by \code{plotrl.evmSim}) and a histogram of
#' the data with superimposed density estimate.  These are all
#' calculated using the data on the original scale. If there are
#' covariates in the model then the diagnostics consist of PP- and QQ-
#' plots calculated by using the model residuals (which will be
#' standard exponential devaiates under the GPD model and standard
#' Gumbel deviates under the GEV model), and plots of residuals versus
#' fitted model parameters.
#'
#' The PP- and QQ-plots show simulated pointwise tolerance intervals.
#' The region is a \eqn{100(1 - \alpha)\%}{100(1-alpha)\%} region based
#' on \code{nsim} simulated samples.
#'
#' @param x an object of class \code{evmOpt}
#' @param main titles for diagnostic plots. Should be a vector of
#'     length 4, with values corresponding to the character strings to
#'     appear on the titles of the pp, qq, return level, and density
#'     estimate plots respectively.
#' @param xlab As for \code{main} but labels for x-axes rather than
#'     titles.
#' @param nsim The number of replicates to be simulated to produce the
#'     simulated tolerance intervals.
#' @param alpha A \eqn{100(1 - \alpha)\%}{100(1 - alpha)\%} simulation
#'     envelope is produced.
#' @param ... FIXME
#' @export
plot.evmOpt <-
function(x, main=rep(NULL,4), xlab=rep(NULL,4), nsim=1000, alpha=.05, ...){
if (!missing(main)){
if (length(main) != 1 & length(main) != 4){
stop("main should have length 1 or 4")
}
else if (length(main) == 1){ main <- rep(main, 4) }
}

plot(ppevm(x, nsim=nsim, alpha=alpha), main=main, xlab=xlab)
plot(qqevm(x, nsim=nsim, alpha=alpha), main=main, xlab=xlab)

if (all(sapply(x$data$D,ncol) == 1)){
plotrl.evmOpt(x, main=main, xlab=xlab, smooth=FALSE, ...)
plot(hist.evmOpt(x), main=main, xlab=xlab)
}
else { # Covariates in the model
np <- length(x$data$D)
lp <- predict(x,type="lp", unique.=FALSE)$obj$link
Which <- as.logical(apply(lp[,1:np],2,var)) # identifies which cols have covariates

for(i in (1:length(x$data$D))[Which]){
ParName <- names(x$data$D[i])
plot(lp[,i],resid(x),main=paste("Residuals vs fitted",ParName),xlab=paste("Fitted",ParName),ylab="Residuals")
panel.smooth(lp[,i], resid(x), col.smooth=2)
}
}

invisible()
}


## Try the texmex package in your browser

Any scripts or data that you put into this service are public.

texmex documentation built on May 2, 2019, 5:41 a.m.