Graphical n sample test of correspondence of distribution functions
Graphical n sample test of correspondence of distribution functions
Compare the distributions of two (or more) samples.
GET.distrequal( x, stat ="ECDF", nsim, r = seq(min(unlist((lapply(x, min)))), max(unlist((lapply(x, max)))), length =100), tau = seq(0.1,0.9, length =100), contrasts =FALSE, GET.args =NULL, density.args =NULL, approxfun.args =NULL, rq.args =NULL, savefuns =FALSE,...)
Arguments
x: A list of numeric vectors, one for each sample.
stat: Character string indicating which test statistic to be used. See details.
nsim: The number of random permutations.
r: The sequence of argument values at which the test functions are to be compared. The default is 100 equally spaced values between the minimum and maximum over all groups.
tau: The sequence of argument values for the QR test statistic. The default values are 100 equally spaced values between 0.1 and 0.9.
contrasts: Logical. FALSE and TRUE specify the two test functions as described in description part of this help file.
GET.args: A named list of additional arguments to be passed to global_envelope_test, e.g. typeone specifies the type of multiple testing control, FWER or FDR. See global_envelope_test for the defaults and available options.
density.args: A named list of additional arguments to be passed for the estimation of the test statistic "DEN". For more details see density.
approxfun.args: A named list of additional arguments to be passed for the estimation of the the test statistic "QQ". For more details see approxfun.
rq.args: A named list of additional arguments to be passed for the estimation of the test statistic "QR". For more details see the function rq of quantreq.
savefuns: Logical. If TRUE, then the functions from permutations are saved to the attribute simfuns.
...: 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
A global envelope test can be performed to investigate whether the n distribution functions differ from each other and how do they differ. This test is a generalization of the two-sample Kolmogorov-Smirnov test with a graphical interpretation. We assume that the observations in the sample i are an i.i.d. sample from the distribution Fi(r),i=1,...,n,
and we want to test the hypothesis
F1(r)=⋯=Fn(r).F1(r)=...=Fn(r).
If contrasts = FALSE (default), then the default test statistic ("ECDF") is taken to be
T=(F^1(r),…,F^n(r))T=(F^1(r),…,F^n(r))
where F^i(r)=(F^i(r1),...,F^i(rk))
is the ecdf of the ith sample evaluated at argument values r=(r1,...,rk).
Another possibility is given by stat = "DIFF", and then the test statistic is still based on the ECDFs and constructed from all pairwise differences,
The choices contrasts = TRUE and stat = "ECDF" lead to the same test statistic. For other (or multiple) values of stat, the argument contrasts is ignored.
All the options as the test statistics are the following:
"ECDF": The ECDFs of the n-samples, as specified above
"DIFF": The pairwise differences between the ECDFs, as specified above
"DEN": The kernel estimated density functions of the n-samples as the test statistic
"QQ": The pairwise comparisons of empirical quantiles
"SHIFT" The de-trended QQ-plot (shift plot)
"QR": The quantile regression process, i.e. the beta-coefficients of the quantile regression. By default, the reference category of this test statistic is the first sample.
The test statistics are described in detail in Konstantinou et al. (2024).
The simulations under the null hypothesis that the distributions are the same are obtained by permuting the individuals of the groups. The default number of permutation, if nsim is not specified, is n∗1000−1 for the case contrasts = FALSE and (n∗(n−1)/2)∗1000−1 for the case contrasts = TRUE, where n is the length of x.
Examples
if(require("fda", quietly=TRUE)){# Heights of boys and girls at age 10 f.a <- growth$hgtf["10",]# girls at age 10 m.a <- growth$hgtm["10",]# boys at age 10# Empirical cumulative distribution functions plot(ecdf(f.a)) plot(ecdf(m.a), col='grey70', add=TRUE)# Create a list of the data fm.list <- list(Girls=f.a, Boys=m.a) res <- GET.distrequal(fm.list) plot(res)# If you want to change the labels: plot(res, scales ="free", labels = c("Girls","Boys"))# If you want to change the x-label (y-label similarly): require("ggplot2") myxlab <- substitute(paste(italic(i)," (", j,")", sep =""), list(i ="x", j ="Height in cm")) plot(res, scales ="free")+ xlab(myxlab)# Use instead the test statistics QQ and DEN res <- GET.distrequal(fm.list, stat = c("QQ","DEN")) plot(res, scales ="free")# Heights of boys and girls at age 14 f.a <- growth$hgtf["14",]# girls at age 14 m.a <- growth$hgtm["14",]# boys at age 14# Empirical cumulative distribution functions plot(ecdf(f.a)) plot(ecdf(m.a), col='grey70', add=TRUE)# Create a list of the data fm.list <- list(Girls=f.a, Boys=m.a) res <- GET.distrequal(fm.list) plot(res)+ xlab(myxlab) res <- GET.distrequal(fm.list, stat = c("QQ","DEN")) plot(res, scales ="free")+ xlab(myxlab)}if(require("datasets", quietly=TRUE)){ data("iris") virginica <- subset(iris, Species =="virginica") setosa <- subset(iris, Species =="setosa") versicolor <- subset(iris, Species =="versicolor") res <- GET.distrequal(x = list(virginica = virginica$Sepal.Length, setosa = setosa$Sepal.Length, versicolor = versicolor$Sepal.Length), stat = c("QQ","DEN")) plot(res, scales ="free", ncol =3)}
References
Konstantinou K., Mrkvička T. and Myllymäki M. (2024) Graphical n-sample tests of correspondence of distributions. arXiv:2403.01838 [stat.ME] https://doi.org/10.48550/arXiv.2403.01838