twinstim_epitest: Permutation Test for Space-Time Interaction in '"twinstim"'

Description Usage Arguments Details Value Author(s) See Also Examples

Description

The function epitest takes a simple epidemic "twinstim" model (one with epidemic = ~1) and tests if the spatio-temporal interaction invoked by the epidemic model component is statistically significant.

Usage

1
2
3
4
5
6
7
8
epitest(model, data, tiles, method = "time", B = 199,
        eps.s = NULL, eps.t = NULL, fixed = NULL,
        verbose = TRUE, compress = FALSE, ...)

## S3 method for class 'epitest'
coef(object, which = c("m1", "m0"), ...)
## S3 method for class 'epitest'
plot(x, teststat = c("simpleR0", "D"), ...)

Arguments

model

a simple epidemic "twinstim" without covariates, i.e., epidemic = ~1. This is because covariate effects in the epidemic component are not well identified when there is no space-time interaction such as in the permuted data. Estimating a rich epidemic model under the null hypothesis of no space-time interaction will most likely result in singular convergence.

data

an object of class "epidataCS", the data to which the model was fitted.

tiles

(only used by method = "simulate") a "SpatialPolygons" representation of the tiles in data$stgrid.

method

one of the following character strings specifying the test method:

"LRT":

a simple likelihood ratio test of the epidemic model against the corresponding endemic-only model,

"time"/"space":

a Monte Carlo permutation test where the null distribution is obtained by relabeling time points or locations, respectively (using permute.epidataCS).

"simulate":

obtain the null distribution of the test statistic by simulations from the endemic-only model (using simulate.twinstim).

B

the number of permutations for the Monte Carlo approach. The default number is rather low; if computationally feasible, B = 999 is more appropriate. Note that this determines the “resolution” of the p-value: the smallest attainable p-value is 1/(B+1).

eps.s,eps.t

arguments for simpleR0.

fixed

optional character vector naming parameters to fix at their original value when re-fitting the model on permuted data. The special value fixed = TRUE means to fix all epidemic parameters but the intercept.

verbose

the amount of tracing in the range 0:3. Set to 0 (or FALSE) for no output, 1 (or TRUE, the default) for a progress bar, 2 for the test statistics resulting from each permutation, and to 3 for additional tracing of the log-likelihood maximization in each permutation (not useful if parallelized). Tracing does not work if permutations are parallelized using clusters. See plapply for other choices.

compress

logical indicating if the nobs-dependent elements "fitted", "fittedComponents", and "R0" should be dropped from the permutation-based model fits. Not keeping these elements saves a lot of memory especially with a large number of events. Note, however, that the returned permfits then no longer are fully valid "twinstim" objects (but most methods will still work).

...

further arguments for plapply to configure parallel operation, i.e., .parallel as well as .seed to make the results reproducible.
For the plot-method, further arguments passed to truehist.
Ignored by the coef-method.

object,x

an object of class "epitest" as returned by epitest.

which

a character string indicating either the full ("m1", default) or the endemic-only ("m0") model.

teststat

a character string determining the test statistic to plot, either "simpleR0" or "D" (twice the log-likelihood difference of the models).

Details

The test statistic is the reproduction number simpleR0. A likelihood ratio test of the supplied epidemic model against the corresponding endemic-only model is also available. The null distribution of the test statistic is obtained by a Monte Carlo permutation approach (via permute.epidataCS) similar to the knox test.

The plot-method shows a truehist of the simulated null distribution together with the observed value. The coef-method extracts the parameter estimates from the B permfits (by default for the full model which = "m1").

Value

a list (inheriting from "htest") with the following components:

method

a character string indicating the type of test performed.

data.name

a character string giving the supplied data and model arguments.

statistic

the observed test statistic.

parameter

the (effective) number of permutations used to calculate the p-value (only those with convergent fits are used).

p.value

the p-value for the test. For the methods involving resampling under the null (method != "LRT"), it is based on the subset of convergent fits only and the p-value from the simple LRT is attached as an attribute "LRT".

In addition, if method != "LRT", the result will have the following elements:

permfits

the list of model fits (endemic-only and epidemic) from the B permutations.

permstats

a data frame with B rows and the columns "l0" (log-likelihood of the endemic-only model m0), "l1" (log-likelihood of the epidemic model m1), "D" (twice their difference), "simpleR0" (the results of simpleR0(m1, eps.s, eps.t)), and "converged" (a boolean indicator if both models converged).

The plot-method invisibly returns NULL. The coef-method returns the B x length(coef(model)) matrix of parameter estimates.

Author(s)

Sebastian Meyer

See Also

permute.epidataCS, knox

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
data("imdepi")
data("imdepifit")

## test for space-time interaction of the B-cases
## assuming spatial interaction to be constant within 50 km
imdepiB50 <- update(subset(imdepi, type == "B"), eps.s = 50)
imdfitB50 <- update(imdepifit, data = imdepiB50,
                    epidemic = ~1, epilink = "identity",
                    siaf = NULL, control.siaf = NULL,
                    start = c("e.(Intercept)" = 1e-6))

## simple likelihood ratio test
epitest(imdfitB50, imdepiB50, method = "LRT")

## permutation test (only few permutations for speed, in parallel)
et <- epitest(imdfitB50, imdepiB50, B = 4 + 25*surveillance.options("allExamples"),
              verbose = 2 * (.Platform$OS.type == "unix"),
              .seed = 1, .parallel = 2)
et
plot(et)
summary(coef(et))

jimhester/surveillance documentation built on May 19, 2019, 10:33 a.m.