sparr (version 0.2-2)

tolerance: Asymptotic p-value surfaces

Description

Calculates pointwise p-values based on asymptotic theory or Monte-Carlo (MC) permutations describing the extremity of risk over a given fixed or adaptive kernel-smoothed relative risk function.

Usage

tolerance(rs, pooled, test = "upper",
	method = "ASY", reduce = 1, ITER = 1000,
	exactL2 = TRUE, comment = TRUE)

Arguments

rs

An object of class "rrs" resulting from a call to risk, giving the fixed or adaptive kernel-smoothed risk function.

pooled

An object of class "bivden" resulting from a call to bivariate.density (or the component pooled from rs if it was created using raw data arguments) represent

test

A character string indicating the kind of test desired to yield the p-values. Must be one of "upper" (default - performs upper tailed tests examining heighted risk `hotspots'), "lower" (lower tailed tests examining `troughs') or

method

A character string, either "ASY" (default) or "MC" indicating which method to use for calculating the p-value surface (asymptotic and Monte-Carlo approaches respectively). The MC approach is far more computationally expensive than the asymptotic method (s

reduce

A numeric value greater than zero and less than or equal to one giving the user the option to reduce the resolution of the evaluation grid for the pointwise p-values by specifying a proportion of the size of the evaluation grid for the original density es

ITER

An integer value specifying the number of iterations to be used if method = "MC" (defaulting to 1000). Non-integer numeric values are rounded. Ignored when method = "ASY".

exactL2

Ignored if rs (and pooled) are fixed-bandwidth density estimates, or if method = "MC". A boolean value indicating whether or not to separately calculate the `L2' integral components for adaptive tolerance contours. A

comment

Boolean. Whether or not to print function progress (including starting and ending times) during execution. Defaults to TRUE.

Value

• A list with four components:
• X
  the equally spaced sequence of length ceiling(reduce*res) giving the evaluation locations on the x axis (where res is the grid resolution as specified in the calls to bivariate.density for calculation of the densities for rs and pooled)
• Y
  as above, for the y axis
• Z
  a numeric ceiling(reduce*res)*ceiling(reduce*res) matrix giving the values of the risk surface over the evaluation grid. Values corresponding to grid coordinates outside the study region are assigned NA. If method = "MC", this will be a single value of NA
• P
  a ceiling(reduce*res)*ceiling(reduce*res) matrix giving the p-values corresponding to the evaluation grid in light of the elected test

Warning

Though far less expensive computationally than calculation of Monte-Carlo p-value surfaces, the asymptotic p-value surfaces (particularly for adaptive relative risk surfaces) can still take some time to complete. The argument of reduce provides an option to reduce this computation time by decreasing the resolution of the evaluation grid. However, the accuracy and appearance of the resulting tolerance contours can be severely degraded if reduce is assigned too small a value. Care must therfore be taken and consideration given to the resolution of the original evaluation grid when altering reduce from its default value. For most practical purposes, we have found a value of reduce resulting in evaluation of a p-value surface of size 50 by 50 is adequate. The MC approach is provided as an option here for the sake of completeness only, and is coded exclusively in R. The computational cost of this approach for the adaptive risk function is enough to recommend against its use in this case, though it is faster for the fixed-bandwidth case if just comparing MC execution times between the two smoothing regimens. Comments on the issue of MC vs ASY are given in Section 3 of Hazelton and Davies (2009).

Details

This function implements developments in Hazelton and Davies (2009) (fixed) and Davies and Hazelton (2010) (adaptive) to compute pointwise p-value surfaces based on asymptotic theory of kernel-smoothed relative risk surfaces. Alternatively, the user may elect to calculate the p-value surfaces using Monte-Carlo methods (see Kelsall and Diggle, 1995). Superimposing upon a plot of the risk surface contours of these p-values at given significance levels (i.e. `tolerance contours') can be an informative way of exploring the statistical significance of the extremity of risk across the defined study region. The asymptotic approach to the p-value calculation is advantageous over a Monte-Carlo method, which can lead to excessive computation time for adaptive risk surfaces and large datasets. See the aforementioned references for further comments. Choosing different options for the argument test simply manipulates the `direction' of the p-values. That is, plotting tolerance contours at a significance level of 0.05 for a p-value surface calculated with test = "double" is equivalent to plotting tolerance contours at significance levels of 0.025 and 0.975 for test = "upper".

References

Kelsall, J.E. and Diggle, P.J. (1995), Kernel estimation of relative risk, Bernoulli, 1, 3-16. Davies, T.M. and Hazelton, M.L. (2010), Adaptive kernel estimation of spatial relative risk, Statistics in Medicine, 29(23) 2423-2437. Hazelton, M.L. and Davies, T.M. (2009), Inference based on kernel estimates of the relative risk function in geographical epidemiology, Biometrical Journal, 51(1), 98-109.

Examples

data(chorley)
ch.h <- CV.sm(chorley)

ch.pool <- bivariate.density(data = chorley, pilotH = ch.h,
 adaptive = FALSE)
ch.case <- bivariate.density(data = chorley, ID = "larynx", pilotH = ch.h,
 adaptive = FALSE)
ch.con <- bivariate.density(data = chorley, ID = "lung", pilotH = ch.h,
 adaptive = FALSE)

##Compute log-risk surface and asymptotic p-value surfaces
ch.rrs <- risk(f = ch.case, g = ch.con) 
ch.tol <- tolerance(rs = ch.rrs, pooled = ch.pool)
contour(ch.tol$X, ch.tol$Y, ch.tol$P, levels = 0.05, add = TRUE)


data(PBC)
PBC.casedata <- PBC$data[PBC$data[,3]==1,1:2]
PBC.controldata <- PBC$data[PBC$data[,3]==0,1:2]

pbc.rrs.rawdata <- risk(f = PBC.casedata, g = PBC.controldata,
 WIN = PBC$owin, adaptive = TRUE, tolerate = TRUE)
 
plot(pbc.rrs.rawdata, display = "3d", aspect = 1:2, col = heat.colors(12)[12:1],
 tolerance.matrix = pbc.rrs.rawdata$P, tol.opt = list(col = "white", raise = 0.03))