flatness (version 0.1.4)

JP_test: Implementation of the Jolliffe-Primo flatness tests for rank histograms

Description

Given a matrix of rank histograms and an orthonormal set of deviate vector(s), this function computes the projection components, test statistics and p-values of the Jolliffe-Primo flatness test for each inputted rank histogram. See Jolliffe and Primo (2008) for details of the method.

Usage

JP_test(rkhists, deviates, ...)

Arguments

rkhists

an rkhist, a matrix, or any other object that can be coerced to a matrix. It contains the rank histogram(s) whose flatness must be tested (one in each row).

deviates

the matrix containing the deviate vectors used for testing. Each row contains a deviate vector: the vector set must be orthonormal, and each deviate vector must have its components summing to zero.

...

further arguments (currently not used).

Value

A list (with additional first class JPtest) with the following entries:

test

an array containing the result of the Jolliffe-Primo test(s). The first dimension is of length three (the projection over the deviate vectors, the test statistics and the p-values), the second and third dimensions correspond to the rank histogram and the test, respectively

deviates

the set of deviate vectors used in the test

rkhist

the tested rank histogram(s) (an rkhist object).

Details

Note that the test statistics and p-values of the projections over the residual vector (after removing all the projection on the deviates) are also computed and returned.

References

Jolliffe, Ian T., and Cristina Primo. "Evaluating rank histograms using decompositions of the chi-square test statistic." Monthly Weather Review 136.6 (2008): 2133-2139. doi:https://doi.org/10.1175/2007MWR2219.1

Examples

# NOT RUN {
require(lattice)
require(xtable)
M <- 15
N <- 100
n <- 20
fcsts <- vector("list", n)
names(fcsts) <- letters[1:n]
obs <- rnorm(N)
for (i in 1:n) {
  fcsts[[i]] <- matrix(rnorm(M*N), ncol = M)
}
rkhsts <- rkhist(fcsts, obs)
deviates <- get_deviates(M + 1)
test <- JP_test(rkhsts, deviates)
print(test)
for (what in c("projections", "statistics", "pvalues")){
  levelplot(test, what = what, main = what, rotate = what == "pvalues")
}
xtable(test$test["pvalues", ,])
xtable(t(test$test["pvalues", ,]))
# }