rank_envelope: The rank envelope test
In myllym/spptest: Spatial point process testing

Description Usage Arguments Details Value References See Also Examples

The rank envelope test, p-value and global envelope

1 2	rank_envelope(curve_set, alpha = 0.05, savedevs = FALSE, alternative = c("two.sided", "less", "greater"), lexo = FALSE, ties)

`curve_set`	A curve_set (see `create_curve_set`) or an `envelope` object. If an envelope object is given, it must contain the summary functions from the simulated patterns which can be achieved by setting savefuns = TRUE when calling envelope().
`alpha`	The significance level. The 100(1-alpha)% global envelope will be calculated.
`savedevs`	Logical. Should the global rank values k_i, i=1,...,nsim+1 be returned? Default: FALSE.
`alternative`	A character string specifying the alternative hypothesis. Must be one of the following: "two.sided" (default), "less" or "greater".
`lexo`	Logical, whether or not to use rank count ordering for calculation of the p-value. See details.
`ties`	Ties method to be passed to `estimate_p_value`. Used to obtain a point estimate for the p-value. The default point estimate is the mid-rank p-value.

The rank envelope test is a completely non-parametric test, which provides a p-value interval given by the most liberal and the most conservative p-value estimate and the 100(1-alpha)% global envelope for the chosen test function T(r) on the chosen interval of distances.

Given a curve_set (or an envelope) object, the test is carried out as follows.

For each curve in the curve_set, both the data curve and the simulations, the global rank measure R is determined. If savedevs = TRUE, then the global rank values R_1, R_2, ..., R_(s+1) are returned in the component 'k', where k[1] is the value for the data.

Based on R_i, i=1, ..., s+1, the p-interval is calculated. This interval is by default plotted for the object returned by the rank_envelope function. Also a single p-value is calculated and returned in component 'p'. By default this p-value is the mid-rank p-value, but another option can be used by specifying ties argument which is passed to estimate_p_value. For options see estimate_p_value.

The 100(1-alpha)% global envelope is given by the 'k_alpha'th lower and upper envelope. For details see Myllymäki et al. (2013).

The above holds for p-value calculation if lexo == FALSE and then the test corresponds to the rank envelope test by Myllymaki et. al (2013). If lexo == TRUE, then all the pointwise ranks are used to rank the curves by rank count ordering (Myllymäki et al., 2015) and the single p-value in p is the p-value based on the rank count ordering.

The rank count ordering test allows in principle a lower number of simulations to be used, but then the test may no longer be usable as a graphical test.

An "envelope_test" object containing the following fields:

r = Distances for which the test was made.
method = The name of the envelope test.
alternative = The alternative specified in the function call.
p = A point estimate for the p-value (default is the mid-rank p-value).
ties = As the argument ties.
p_interval = The p-value interval [p_liberal, p_conservative].
k_alpha = The value of k corresponding to the 100(1-alpha)% global envelope.
k = Global rank values (k[1] is the value for the data pattern). Returned only if savedevs = TRUE.
central_curve = If the curve_set (or envelope object) contains a component 'theo', then this function is used as the central curve and returned in this component. Otherwise, the central_curve is the mean of the test functions T_i(r), i=2, ..., s+1. Used for visualization only.
data_curve = The test function for the data.
lower = The lower envelope.
upper = The upper envelope.
call = The call of the function.

Myllymäki, M., Mrkvička, T., Seijo, H., Grabarnik, P. (2013). Global envelope tests for spatial point patterns. arXiv:1307.0239 [stat.ME]

Myllymäki, M., Mrkvička, T., Grabarnik, P., Seijo, H. and Hahn, U. (2015). Global envelope tests for spatial point patterns. arXiv:1307.0239v4 [stat.ME]

random_labelling, plot.envelope_test

## Testing complete spatial randomness (CSR)
#-------------------------------------------
require(spatstat)
pp <- unmark(spruces)
# Generate nsim simulations under CSR, calculate L-function for the data and simulations
env <- envelope(pp, fun="Lest", nsim=2499, savefuns=TRUE, correction="translate")
# The rank envelope test
res <- rank_envelope(env)
# Plot the result.
# - The central curve is now obtained from env[['theo']], which is the
# value of the L-function under the null hypothesis (L(r) = r).
plot(res)
# or (requires R library ggplot2)
plot(res, use_ggplot2=TRUE)

## Advanced use:
# Choose the interval of distances [r_min, r_max] (at the same time create a curve_set from 'env')
curve_set <- crop_curves(env, r_min = 1, r_max = 7)
# For better visualisation, take the L(r)-r function
curve_set <- residual(curve_set, use_theo = TRUE)
# Do the rank envelope test
res <- rank_envelope(curve_set); plot(res, use_ggplot2=TRUE)

## Random labeling test
#----------------------
# requires library 'marksummary'
mpp <- spruces
# 1) Perform simulations under the random labelling hypothesis and calculate
# the test function T(r) for the data pattern (mpp) and each simulation.
# The command below specifies that the test function is T(r) = \hat{L}_m(r),
# which is an estimator of the mark-weighted L function, L_m(r),
# with translational edge correction (default).
# The random_labelling function returns the centred functions \hat{L}_m(r)-T_0(r),
# where T_0(r) = \hat{L}(r) is the unmarked L function.
curve_set <- random_labelling(mpp, mtf_name = 'm', nsim=2499, r_min=1.5, r_max=9.5)
# 2) Do the rank envelope test
res <- rank_envelope(curve_set)
# 3) Plot the test result
plot(res, use_ggplot2=TRUE, ylab=expression(italic(L[m](r)-L(r))))

# Make the test using instead the test function T(r) = \hat{L}_mm(r);
# which is an estimator of the mark-weighted L function, L_mm(r),
# with translational edge correction (default).
curve_set <- random_labelling(mpp, mtf_name = 'mm', nsim=2499, r_min=1.5, r_max=9.5)
res <- rank_envelope(curve_set)
plot(res, use_ggplot2=TRUE, ylab=expression(italic(L[mm](r)-L(r))))

## Goodness-of-fit test (typically conservative)
#-----------------------------------------------
pp <- unmark(spruces)
# Minimum distance between points in the pattern
min(nndist(pp))
# Fit a model
fittedmodel <- ppm(pp, interaction=Hardcore(hc=1)) # Hardcore process

## Not run: 
# Simulating Gibbs process by 'envelope' is slow, because it uses the MCMC algorithm
#env <- envelope(fittedmodel, fun="Jest", nsim=999, savefuns=TRUE,
                 correction="none", r=seq(0, 4, length=500))

# Using direct algorihm can be faster, because the perfect simulation is used here.
simulations <- NULL
for(j in 1:2499) {
   simulations[[j]] <- rHardcore(beta=exp(fittedmodel$coef[1]),
                                 R = fittedmodel$interaction$par$hc,
                                 W = pp$window);
   if(j%%10==0) cat(j, "...", sep="")
}
env <- envelope(pp, simulate=simulations, fun="Jest", nsim=length(simulations),
                savefuns=TRUE, correction="none", r=seq(0, 4, length=500))
curve_set <- crop_curves(env, r_min = 1, r_max = 3.5)
res <- rank_envelope(curve_set); plot(res, use_ggplot2=TRUE)

## End(Not run)