R/residuals.R

## Export: inla.ks.plot

##!\name{inla.ks.plot}
##!\alias{inla.ks.plot}
##!\alias{ks.plot}
##!\title{
##!  Kolmogorov-Smirnov Test Plots
##!}
##!\description{
##!  Illustrate a one-sample Kolmogorov-Smirnov test by plotting the
##!  empirical distribution deviation.
##!}
##!\usage{
##!inla.ks.plot(x, y, diff=TRUE, ...)
##!}
##!\arguments{
##!  \item{x}{a numeric vector of data values.}
##!  \item{y}{a cumulative distribution function such as 'pnorm'.}
##!  \item{diff}{logical, indicating if the normalised difference should be
##!    plotted.  If \code{FALSE}, the absolute distribution functions are plotted.}
##!  \item{...}{additional arguments for \code{\link{ks.test}}, ignored in
##!    the plotting.  In particular, only two-sided tests are illustrated.}
##!}
##!\details{
##!  In addition to the (normalised) empirical distribution deviation,
##!  lines for the K-S test statistic are drawn, as well as \eqn{\pm}{+/-}
##!  two standard deviations around the expectation under the null hypothesis.
##!}
##!\value{
##!  A list with class \code{"htest"}, as generated by
##!  \code{\link{ks.test}}
##!}
##!%%\references{
##!%%}
##!\author{
##!  Finn Lindgren \email{finn.lindgren@gmail.com}
##!}
##!%%\note{
##!%%}
##!
##!%% ~Make other sections like Warning with \section{Warning }{....} ~
##!
##!\seealso{
##!  \code{\link{ks.test}}
##!}
##!\examples{
##!## Check for N(0,1) data
##!data = rowSums(matrix(runif(100*12)*2-1,100,12))/2
##!inla.ks.plot(data, pnorm)
##!
##!\dontrun{
##!## Check the goodness-of-fit of cross-validated predictions
##!result = inla(..., control.predictor=list(cpo=TRUE))
##!inla.ks.plot(result$pit, punif)
##!}
##!}



## One-sample Kolmogorov-Smirnov test, plotting the normalised
## deviation between the empirical distribution function and the null
## hypothesis.
##
## Example:
## result = inla(..., control.compute=list(cpo=TRUE))
## inla.ks.plot(result$pit, punif)
inla.ks.plot = function (x, y, diff=TRUE, ...)
{
    if (any(is.na(x))) {
        x = x[!is.na(x)]
    }
    test = ks.test(x, y, ...)
    n = length(x)
    Fn = ((1:n)-0.5)/n
    F = y(sort(x))
    empirical.diff = (Fn-F)*sqrt(n)
    T = max(abs(empirical.diff))
    if (diff) {
        ylim = c(-1,1)*max(1,T)
        plot(F, empirical.diff, type='l',
             xlim=c(0,1),
             ylim=ylim,
             main=paste("K-S-test, p-value = ", test$p.value),
             ylab="(Fn-F) sqrt(n)",
             xlab="Quantile"
             )
        lines(c(0,1), c(0,0), type='l')
        lines(Fn, 2*sqrt(Fn*(1-Fn)), type='l')
        lines(Fn, -2*sqrt(Fn*(1-Fn)), type='l')
        lines(c(0,1,NA,0,1), c(T,T,NA,-T,-T), type='l')
    } else {
        plot(F, Fn, type='l',
             xlim=c(0,1),
             ylim=c(0,1),
             main=paste("K-S-test, p-value = ", test$p.value),
             ylab="Fn and F",
             xlab="Quantile"
             )
        lines(c(0,1), c(0,1), type='l')
        lines(Fn, Fn+2*sqrt(Fn*(1-Fn)/n), type='l')
        lines(Fn, Fn-2*sqrt(Fn*(1-Fn)/n), type='l')
        lines(c(0,1,NA,0,1), c(0,1,NA,0,1)+c(T,T,NA,-T,-T)/sqrt(n), type='l')
    }
    invisible(test)
}
andrewzm/INLA documentation built on May 10, 2019, 11:12 a.m.