graph.fanova: One-way graphical functional ANOVA

In myllym/GET: Global Envelopes

One-way graphical functional ANOVA

Description

One-way ANOVA tests for functional data with graphical interpretation

Usage

graph.fanova(
  nsim,
  curve_set,
  groups,
  variances = "equal",
  contrasts = FALSE,
  n.aver = 1L,
  mirror = FALSE,
  savefuns = FALSE,
  test.equality = c("mean", "var", "cov"),
  cov.lag = 1,
  ...
)

Arguments

nsim	The number of random permutations.
curve_set	The original data (an array of functions) provided as a curve_set object or a fdata object of fda.usc. The curve set should include the argument values for the functions in the component r, and the observed functions in the component obs.
groups	The original groups (a factor vector representing the assignment to groups).
variances	Either "equal" or "unequal". If "unequal", then correction for unequal variances as explained in details will be done. Only relevant for the case test.equality = "means" (default).
contrasts	Logical. FALSE and TRUE specify the two test functions as described in description part of this help file.
n.aver	If variances = "unequal", there is a possibility to use variances smoothed by appying moving average to the estimated sample variances. n.aver determines how many values on each side do contribute (incl. value itself).
mirror	The complement of the argument circular of filter. Another parameter for the moving average to estimate sample variances (see n.aver).
savefuns	Logical. If TRUE, then the functions from permutations are saved to the attribute simfuns.
test.equality	A character with possible values mean (default), var and cov. If mean, the functional ANOVA is performed to compare the means in the groups. If var, then the equality of variances of the curves in the groups is tested by performing the graphical functional ANOVA test on the functions Z_{ij}(r) = T_{ij}(r) - \bar{T}_j(r). If cov, then the equality of lag cov.lag covariance is tested by performing the fANOVA with W_{ij}(r) = \sqrt{|V_{ij}(r)|\cdot sign(V_{ij}(r))}, where V_{ij}(r) = (T_{ij}(r) - \bar{T}_j(r))((T_{ij}(r+s) - \bar{T}_j(r+s))). See Mrkvicka et al. (2020) for more details.
cov.lag	The lag of the covariance for testing the equality of covariances, see test.equality.
...	Additional parameters to be passed to global_envelope_test. For example, the type of multiple testing control, FWER or FDR must be set by typeone. And, if typeone = "fwer", the type of the global envelope can be chosen by specifying the argument type. See global_envelope_test for the defaults and available options. (The test here uses alternative="two.sided" and nstep=1 (when relevant), but all the other specifications are to be given in ....)

Details

This function can be used to perform one-way graphical functional ANOVA tests described in Mrkvička et al. (2020). Both 1d and 2d functions are allowed in curve sets.

The tests assume that there are J groups which contain n_1,\dots,n_J functions T_{ij}, i=\dots,J, j=1,\dots,n_j. The functions should be given in the argument curve_set, and the groups in the argument groups. The tests assume that T_{ij}, i=1,...,n_j is an iid sample from a stochastic process with mean function \mu_j and covariance function \gamma_j(s,t) for s,t in R and j = 1,..., J.

To test the hypothesis

H_0 : \mu_1(r) \equiv \mu_2(r)\equiv \dots \equiv \mu_J(r),

you can use the test function

\mathbf{T} = (\overline{T}_1({\bf r}), \overline{T}_2({\bf r}), \dots , \overline{T}_J({\bf r}))

where \overline{T}_i({\bf r}) is a vector of mean values of functions in the group j. This test function is used when contrasts = FALSE (default).

The hypothesis can equivalently be written as

H_0 : \mu_i(r) - \mu_j(r) = 0, i=1,\dots,J-1, j=1,\dots,J.

and, alternatively, one can use the test function (vector) taken to consist of the differences of the group averages,

\mathbf{T'} = (\overline{T}_1({\bf r})-\overline{T}_2({\bf r}), \overline{T}_1({\bf r})-\overline{T}_3({\bf r}), \dots , \overline{T}_{J-1}({\bf r})-\overline{T}_J({\bf r})).

The choice is available with the option contrasts = TRUE. This test corresponds to the post-hoc test done usually after an ANOVA test is significant, but it can be directed tested by means of the combined rank test (Mrkvička et al., 2017) with this test vector.

The test as such assumes that the variances are equal across the groups of functions. To deal with unequal variances, the differences are rescaled as the first step as follows

S_{ij}(r) = \frac{T_{ij}(r) - \overline{T}(r))}{\sqrt{Var(T_j(r))}} \sqrt{Var(T(r))} + \overline{T}(r))

where \overline{T}({\bf r}) is the overall sample mean and \sqrt{Var(T(r))} is the overall sample standard deviation. This scaling of the test functions can be obtained by giving the argument variances = "unequal".

References

Mrkvička, T., Myllymäki, M., Jilek, M. and Hahn, U. (2020) A one-way ANOVA test for functional data with graphical interpretation. Kybernetika 56 (3), 432-458. doi: 10.14736/kyb-2020-3-0432

Mrkvička, T., Myllymäki, M., and Hahn, U. (2017). Multiple Monte Carlo testing, with applications in spatial point processes. Statistics and Computing 27 (5): 1239-1255. doi:10.1007/s11222-016-9683-9

Myllymäki, M and Mrkvička, T. (2020). GET: Global envelopes in R. arXiv:1911.06583 [stat.ME]

See Also

frank.fanova

Examples

#-- NOx levels example (see for details Myllymaki and Mrkvicka, 2020)
if(require("fda.usc", quietly=TRUE)) {
  # Prepare data
  data("poblenou")
  fest <- poblenou$df$day.festive; week <- as.integer(poblenou$df$day.week)
  Type <- vector(length=length(fest))
  Type[fest == 1 | week >= 6] <- "Free"
  Type[fest == 0 & week %in% 1:4] <- "MonThu"
  Type[fest == 0 & week == 5] <- "Fri"
  Type <- factor(Type, levels = c("MonThu", "Fri", "Free"))

  # (log) Data as a curve_set
  cset <- curve_set(r = 0:23,
     obs = t(log(poblenou[['nox']][['data']])))
  # Graphical functional ANOVA
  nsim <- 2999
  
  res.c <- graph.fanova(nsim = nsim, curve_set = cset,
                        groups = Type, variances = "unequal",
                        contrasts = TRUE)
  plot(res.c) + ggplot2::labs(x = "Hour", y = "Diff.")
}

#-- Centred government expenditure centralization ratios example
# This is an example analysis of the centred GEC in Mrkvicka et al.
data("cgec")

# Number of simulations
nsim <- 2499 # increase to reduce Monte Carlo error


# Test for unequal lag 1 covariances
res.cov1 <- graph.fanova(nsim = nsim, curve_set = cgec$cgec,
                         groups = cgec$group,
                         test.equality = "cov", cov.lag = 1)
plot(res.cov1)
# Add labels
plot(res.cov1, labels = paste("Group ", 1:3, sep="")) +
  ggplot2::xlab(substitute(paste(italic(i), " (", j, ")", sep=""), list(i="r", j="Year")))
# Test for equality of variances among groups
res.var <- graph.fanova(nsim = nsim, curve_set = cgec$cgec,
                        groups = cgec$group,
                        test.equality = "var")
plot(res.var)

# Test for equality of means assuming equality of variances
# a) using 'means'
res <- graph.fanova(nsim = nsim, curve_set = cgec$cgec,
                    groups = cgec$group,
                    variances = "equal", contrasts = FALSE)
plot(res)
# b) using 'contrasts'
res2 <- graph.fanova(nsim = nsim, curve_set = cgec$cgec,
                     groups = cgec$group,
                     variances = "equal", contrasts = TRUE)
plot(res2)

# Image set examples
data("imageset3")

res <- graph.fanova(nsim = 19, # Increase nsim for serious analysis!
                    curve_set = imageset3$image_set,
                    groups = imageset3$Group)
plot(res, what = c("obs", "lo.sign", "hi.sign"), sign.type = "col")

myllym/GET documentation built on Feb. 4, 2024, 10:44 p.m.