Type:	Package
Title:	Calculate Pairwise Multiple Comparisons of Mean Rank Sums Extended
Version:	1.9.12
Date:	2024-09-07
Description:	For one-way layout experiments the one-way ANOVA can be performed as an omnibus test. All-pairs multiple comparisons tests (Tukey-Kramer test, Scheffe test, LSD-test) and many-to-one tests (Dunnett test) for normally distributed residuals and equal within variance are available. Furthermore, all-pairs tests (Games-Howell test, Tamhane's T2 test, Dunnett T3 test, Ury-Wiggins-Hochberg test) and many-to-one (Tamhane-Dunnett Test) for normally distributed residuals and heterogeneous variances are provided. Van der Waerden's normal scores test for omnibus, all-pairs and many-to-one tests is provided for non-normally distributed residuals and homogeneous variances. The Kruskal-Wallis, BWS and Anderson-Darling omnibus test and all-pairs tests (Nemenyi test, Dunn test, Conover test, Dwass-Steele-Critchlow- Fligner test) as well as many-to-one (Nemenyi test, Dunn test, U-test) are given for the analysis of variance by ranks. Non-parametric trend tests (Jonckheere test, Cuzick test, Johnson-Mehrotra test, Spearman test) are included. In addition, a Friedman-test for one-way ANOVA with repeated measures on ranks (CRBD) and Skillings-Mack test for unbalanced CRBD is provided with consequent all-pairs tests (Nemenyi test, Siegel test, Miller test, Conover test, Exact test) and many-to-one tests (Nemenyi test, Demsar test, Exact test). A trend can be tested with Pages's test. Durbin's test for a two-way balanced incomplete block design (BIBD) is given in this package as well as Gore's test for CRBD with multiple observations per cell is given. Outlier tests, Mandel's k- and h statistic as well as functions for Type I error and Power analysis as well as generic summary, print and plot methods are provided.
Depends:	R (≥ 3.5.0)
Imports:	mvtnorm (≥ 1.0), multcompView, gmp, Rmpfr, SuppDists, kSamples (≥ 1.2.7), BWStest (≥ 0.2.1), MASS, stats
Suggests:	xtable, graphics, knitr, rmarkdown, car, e1071, multcomp, pwr, NSM3
SystemRequirements:	gmp (>= 4.2.3), mpfr (>= 3.0.0) | file README.md
SystemRequirementsNote:	see >> README.md
SysDataCompression:	gzip
VignetteBuilder:	knitr, rmarkdown
Classification/MSC-2010:	62J15, 62J10, 62G10, 62F03, 62G30
NeedsCompilation:	yes
Encoding:	UTF-8
LazyData:	true
RoxygenNote:	7.3.1
License:	GPL (≥ 3)
Packaged:	2024-09-08 09:18:42 UTC; thorsten
Author:	Thorsten Pohlert [aut, cre]
Maintainer:	Thorsten Pohlert <thorsten.pohlert@gmx.de>
Repository:	CRAN
Date/Publication:	2024-09-08 10:10:03 UTC

Cochran's distribution

Description

Distribution function and quantile function for Cochran's distribution.

Usage

qcochran(p, k, n, lower.tail = TRUE, log.p = FALSE)

pcochran(q, k, n, lower.tail = TRUE, log.p = FALSE)

Arguments

Value

pcochran gives the distribution function and qcochran gives the quantile function.

References

Cochran, W.G. (1941) The distribution of the largest of a set of estimated variances as a fraction of their total. Ann. Eugen. 11, 47–52.

Wilrich, P.-T. (2011) Critical values of Mandel's h and k, Grubbs and the Cochran test statistic. Adv. Stat. Anal.. doi:10.1007/s10182-011-0185-y.

See Also

FDist

Examples

qcochran(0.05, 7, 3)

Grubbs D* distribution

Description

Distribution function for Grubbs D* distribution.

Usage

pdgrubbs(q, n, m = 10000, lower.tail = TRUE, log.p = FALSE)

Arguments

Value

pgrubbs gives the distribution function

References

Grubbs, F.E. (1950) Sample criteria for testing outlying observations, Ann. Math. Stat. 21, 27–58.

Wilrich, P.-T. (2011) Critical values of Mandel's h and k, Grubbs and the Cochran test statistic, Adv. Stat. Anal.. doi:10.1007/s10182-011-0185-y.

See Also

Grubbs

Examples

pdgrubbs(0.62, 7, 1E4)

Generalized Siegel-Tukey Test of Homogeneity of Scales

Description

Performs a Siegel-Tukey k-sample rank dispersion test.

Usage

GSTTest(x, ...)

## Default S3 method:
GSTTest(x, g, dist = c("Chisquare", "KruskalWallis"), ...)

## S3 method for class 'formula'
GSTTest(
  formula,
  data,
  subset,
  na.action,
  dist = c("Chisquare", "KruskalWallis"),
  ...
)

Arguments

Details

Meyer-Bahlburg (1970) has proposed a generalized Siegel-Tukey rank dispersion test for the k-sample case. Likewise to the fligner.test, this test is a nonparametric test for testing the homogegeneity of scales in several groups. Let \theta_i, and \lambda_i denote location and scale parameter of the ith group, then for the two-tailed case, the null hypothesis H: \lambda_i / \lambda_j = 1 | \theta_i = \theta_j, ~ i \ne j is tested against the alternative, A: \lambda_i / \lambda_j \ne 1 with at least one inequality beeing strict.

The data are combinedly ranked according to Siegel-Tukey. The ranking is done by alternate extremes (rank 1 is lowest, 2 and 3 are the two highest, 4 and 5 are the two next lowest, etc.).

Meyer-Bahlburg (1970) showed, that the Kruskal-Wallis H-test can be employed on the Siegel-Tukey ranks. The H-statistic is assymptotically chi-squared distributed with v = k - 1 degree of freedom, the default test distribution is consequently dist = "Chisquare". If dist = "KruskalWallis" is selected, an incomplete beta approximation is used for the calculation of p-values as implemented in the function pKruskalWallis of the package SuppDists.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

If ties are present, a tie correction is performed and a warning message is given. The GSTTest is sensitive to median differences, likewise to the Siegel-Tukey test. It is thus appropriate to apply this test on the residuals of a one-way ANOVA, rather than on the original data (see example).

References

H.F.L. Meyer-Bahlburg (1970), A nonparametric test for relative spread in k unpaired samples, Metrika 15, 23–29.

See Also

fligner.test, pKruskalWallis, Chisquare, fligner.test

Examples

GSTTest(count ~ spray, data = InsectSprays)

## as means/medians differ, apply the test to residuals
## of one-way ANOVA
ans <- aov(count ~ spray, data = InsectSprays)
GSTTest( residuals( ans) ~ spray, data =InsectSprays)

Grubbs distribution

Description

Distribution function and quantile function for Grubbs distribution.

Usage

qgrubbs(p, n)

pgrubbs(q, n, lower.tail = TRUE)

Arguments

Value

pgrubbs gives the distribution function and qgrubbs gives the quantile function.

References

Grubbs, F. E. (1950) Sample criteria for testing outlying observations. Ann. Math. Stat. 21, 27–58.

Wilrich, P.-T. (2011) Critical values of Mandel's h and k, Grubbs and the Cochran test statistic. Adv. Stat. Anal.. doi:10.1007/s10182-011-0185-y.

See Also

TDist

Examples

qgrubbs(0.05, 7)

Extended One-Sided Studentised Range Test

Description

Performs Nashimoto-Wright's extended one-sided studentised range test against an ordered alternative for normal data with equal variances.

Usage

MTest(x, ...)

## Default S3 method:
MTest(x, g, alternative = c("greater", "less"), ...)

## S3 method for class 'formula'
MTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  ...
)

## S3 method for class 'aov'
MTest(x, alternative = c("greater", "less"), ...)

Arguments

Details

The procedure uses the property of a simple order, \theta_m' - \mu_m \le \mu_j - \mu_i \le \mu_l' - \mu_l \qquad (l \le i \le m~\mathrm{and}~ m' \le j \le l'). The null hypothesis H_{ij}: \mu_i = \mu_j is tested against the alternative A_{ij}: \mu_i < \mu_j for any 1 \le i < j \le k.

The all-pairs comparisons test statistics for a balanced design are

\hat{h}_{ij} = \max_{i \le m < m' \le j} \frac{\left(\bar{x}_{m'} - \bar{x}_m \right)} {s_{\mathrm{in}} / \sqrt{n}},

with n = n_i; ~ N = \sum_i^k n_i ~~ (1 \le i \le k), \bar{x}_i the arithmetic mean of the ith group, and s_{\mathrm{in}}^2 the within ANOVA variance. The null hypothesis is rejected, if \hat{h} > h_{k,\alpha,v}, with v = N - k degree of freedom.

For the unbalanced case with moderate imbalance the test statistic is

\hat{h}_{ij} = \max_{i \le m < m' \le j} \frac{\left(\bar{x}_{m'} - \bar{x}_m \right)} {s_{\mathrm{in}} \left(1/n_m + 1/n_{m'}\right)^{1/2}},

The null hypothesis is rejected, if \hat{h}_{ij} > h_{k,\alpha,v} / \sqrt{2}.

The function does not return p-values. Instead the critical h-values as given in the tables of Hayter (1990) for \alpha = 0.05 (one-sided) are looked up according to the number of groups (k) and the degree of freedoms (v).

Value

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated statistic(s)

crit.value

critical values for \alpha = 0.05.

alternative

a character string describing the alternative hypothesis.

parameter

the parameter(s) of the test distribution.

dist

a string that denotes the test distribution.

There are print and summary methods available.

Note

The function will give a warning for the unbalanced case and returns the critical value h_{k,\alpha,\infty} / \sqrt{2}.

References

Hayter, A. J.(1990) A One-Sided Studentised Range Test for Testing Against a Simple Ordered Alternative, Journal of the American Statistical Association 85, 778–785.

Nashimoto, K., Wright, F.T., (2005) Multiple comparison procedures for detecting differences in simply ordered means. Comput. Statist. Data Anal. 48, 291–306.

See Also

osrtTest, NPMTest

Examples

##
md <- aov(weight ~ group, PlantGrowth)
anova(md)
osrtTest(md)
MTest(md)

Mandel's h Distribution

Description

Distribution function and quantile function for Mandel's h distribution.

Usage

qmandelh(p, k, lower.tail = TRUE, log.p = FALSE)

pmandelh(q, k, lower.tail = TRUE, log.p = FALSE)

Arguments

Value

pmandelh gives the distribution function and qmandelh gives the quantile function.

Source

The code for pmandelh was taken from:
Stephen L R Ellison. (2017). metRology: Support for Metrological Applications. R package version 0.9-26-2. https://CRAN.R-project.org/package=metRology

References

Practice E 691 (2005) Standard Practice for Conducting an Interlaboratory Study to Determine the Precision of a Test Method, ASTM International.

See Also

mandelhTest

Examples

## We need a two-sided upper-tail quantile
qmandelh(p = 0.005/2, k = 7, lower.tail=FALSE)

Mandel's k Distribution

Description

Distribution function and quantile function for Mandel's k distribution.

Usage

qmandelk(p, k, n, lower.tail = TRUE, log.p = FALSE)

pmandelk(q, k, n, lower.tail = TRUE, log.p = FALSE)

Arguments

Value

pmandelk gives the distribution function and qmandelk gives the quantile function.

Source

The code for pmandelk was taken from:
Stephen L R Ellison. (2017). metRology: Support for Metrological Applications. R package version 0.9-26-2. https://CRAN.R-project.org/package=metRology

Note

The functions are only appropriate for balanced designs.

References

Practice E 691 (2005) Standard Practice for Conducting an Interlaboratory Study to Determine the Precision of a Test Method, ASTM International.

See Also

mandelkTest

pmandelh, qmandelh

Examples

qmandelk(0.005, 7, 3, lower.tail=FALSE)

All-Pairs Comparisons for Simply Ordered Mean Ranksums

Description

Performs Nashimoto and Wright's all-pairs comparison procedure for simply ordered mean ranksums.

Usage

NPMTest(x, ...)

## Default S3 method:
NPMTest(
  x,
  g,
  alternative = c("greater", "less"),
  method = c("look-up", "boot", "asympt"),
  nperm = 10000,
  ...
)

## S3 method for class 'formula'
NPMTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  method = c("look-up", "boot", "asympt"),
  nperm = 10000,
  ...
)

Arguments

Details

The procedure uses the property of a simple order, \theta_m' - \theta_m \le \theta_j - \theta_i \le \theta_l' - \theta_l \qquad (l \le i \le m~\mathrm{and}~ m' \le j \le l'). The null hypothesis H_{ij}: \theta_i = \theta_j is tested against the alternative A_{ij}: \theta_i < \theta_j for any 1 \le i < j \le k.

The all-pairs comparisons test statistics for a balanced design are

\hat{h}_{ij} = \max_{i \le m < m' \le j} \frac{\left(\bar{R}_{m'} - \bar{R}_m \right)}{\sigma_a / \sqrt{n}},

with n = n_i; ~ N = \sum_i^k n_i ~~ (1 \le i \le k), \bar{R}_i the mean rank for the ith group, and \sigma_a = \sqrt{N \left(N + 1 \right) / 12}. The null hypothesis is rejected, if h_{ij} > h_{k,\alpha,\infty}.

For the unbalanced case with moderate imbalance the test statistic is

\hat{h}_{ij} = \max_{i \le m < m' \le j} \frac{\left(\bar{R}_{m'} - \bar{R}_m \right)} {\sigma_a \left(1/n_m + 1/n_{m'}\right)^{1/2}},

The null hypothesis is rejected, if \hat{h}_{ij} > h_{k,\alpha,\infty} / \sqrt{2}.

If method = "look-up" the function will not return p-values. Instead the critical h-values as given in the tables of Hayter (1990) for \alpha = 0.05 (one-sided) are looked up according to the number of groups (k) and the degree of freedoms (v = \infty).

If method = "boot" an asymetric permutation test is conducted and p-values is returned.

If method = "asympt" is selected the asymptotic p-value is estimated as implemented in the function pHayStonLSA of the package NSM3.

Value

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated statistic(s)

crit.value

critical values for \alpha = 0.05.

alternative

a character string describing the alternative hypothesis.

parameter

the parameter(s) of the test distribution.

dist

a string that denotes the test distribution.

There are print and summary methods available.

Either a list of class "PMCMR" or a list with class "osrt" that contains the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated statistic(s)

crit.value

critical values for \alpha = 0.05.

alternative

a character string describing the alternative hypothesis.

parameter

the parameter(s) of the test distribution.

dist

a string that denotes the test distribution.

There are print and summary methods available.

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

The function will give a warning for the unbalanced case and returns the critical value h_{k,\alpha,\infty} / \sqrt{2}.

Source

If method = "asympt" is selected, this function calls an internal probability function pHS. The GPL-2 code for this function was taken from pHayStonLSA of the the package NSM3:

Grant Schneider, Eric Chicken and Rachel Becvarik (2020) NSM3: Functions and Datasets to Accompany Hollander, Wolfe, and Chicken - Nonparametric Statistical Methods, Third Edition. R package version 1.15. https://CRAN.R-project.org/package=NSM3

References

Hayter, A. J.(1990) A One-Sided Studentised Range Test for Testing Against a Simple Ordered Alternative, Journal of the American Statistical Association 85, 778–785.

Nashimoto, K., Wright, F.T. (2007) Nonparametric Multiple-Comparison Methods for Simply Ordered Medians. Comput Stat Data Anal 51, 5068–5076.

See Also

MTest

Examples

## Example from Shirley (1977)
## Reaction times of mice to stimuli to their tails.
x <- c(2.4, 3, 3, 2.2, 2.2, 2.2, 2.2, 2.8, 2, 3,
 2.8, 2.2, 3.8, 9.4, 8.4, 3, 3.2, 4.4, 3.2, 7.4, 9.8, 3.2, 5.8,
 7.8, 2.6, 2.2, 6.2, 9.4, 7.8, 3.4, 7, 9.8, 9.4, 8.8, 8.8, 3.4,
 9, 8.4, 2.4, 7.8)
g <- gl(4, 10)

## Shirley's test
## one-sided test using look-up table
shirleyWilliamsTest(x ~ g, alternative = "greater")

## Chacko's global hypothesis test for 'greater'
chackoTest(x , g)

## post-hoc test, default is standard normal distribution (NPT'-test)
summary(chaAllPairsNashimotoTest(x, g, p.adjust.method = "none"))

## same but h-distribution (NPY'-test)
chaAllPairsNashimotoTest(x, g, dist = "h")

## NPM-test
NPMTest(x, g)

## Hayter-Stone test
hayterStoneTest(x, g)

## all-pairs comparisons
hsAllPairsTest(x, g)

Pentosan Dataset

Description

A benchmark dataset of an interlaboratory study for determining the precision of a test method on several levels of the material Pentosan.

Format

A data frame with 189 obs. of 3 variables:

value

numeric, test result (no unit specified)

lab

factor, identifier of the lab (1–7)

material

factor, identifier of the level of the material (A–I)

Source

Tab. 8, Practice E 691, 2005, Standard Practice for Conducting an Interlaboratory Study to Determine the Precision of a Test Method, ASTM International.

Anderson-Darling All-Pairs Comparison Test

Description

Performs Anderson-Darling all-pairs comparison test.

Usage

adAllPairsTest(x, ...)

## Default S3 method:
adAllPairsTest(x, g, p.adjust.method = p.adjust.methods, ...)

## S3 method for class 'formula'
adAllPairsTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Anderson-Darling's all-pairs comparison test can be used. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: F_i(x) = F_j(x) is tested in the two-tailed test against the alternative A_{ij}: F_i(x) \ne F_j(x), ~~ i \ne j.

This function is a wrapper function that sequentially calls adKSampleTest for each pair. The calculated p-values for Pr(>|T2N|) can be adjusted to account for Type I error multiplicity using any method as implemented in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Scholz, F.W., Stephens, M.A. (1987) K-Sample Anderson-Darling Tests. Journal of the American Statistical Association 82, 918–924.

See Also

adKSampleTest, adManyOneTest, ad.pval.

Examples

adKSampleTest(count ~ spray, InsectSprays)

out <- adAllPairsTest(count ~ spray, InsectSprays, p.adjust="holm")
summary(out)
summaryGroup(out)

Anderson-Darling k-Sample Test

Description

Performs Anderson-Darling k-sample test.

Usage

adKSampleTest(x, ...)

## Default S3 method:
adKSampleTest(x, g, ...)

## S3 method for class 'formula'
adKSampleTest(formula, data, subset, na.action, ...)

Arguments

Details

The null hypothesis, H_0: F_1 = F_2 = \ldots = F_k is tested against the alternative, H_\mathrm{A}: F_i \ne F_j ~~(i \ne j), with at least one unequality beeing strict.

This function only evaluates version 1 of the k-sample Anderson-Darling test (i.e. Eq. 6) of Scholz and Stephens (1987). The p-values are estimated with the extended empirical function as implemented in ad.pval of the package kSamples.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Scholz, F.W., Stephens, M.A. (1987) K-Sample Anderson-Darling Tests. Journal of the American Statistical Association 82, 918–924.

See Also

adAllPairsTest, adManyOneTest, ad.pval.

Examples

## Hollander & Wolfe (1973), 116.
## Mucociliary efficiency from the rate of removal of dust in normal
## subjects, subjects with obstructive airway disease, and subjects
## with asbestosis.
x <- c(2.9, 3.0, 2.5, 2.6, 3.2) # normal subjects
y <- c(3.8, 2.7, 4.0, 2.4)      # with obstructive airway disease
z <- c(2.8, 3.4, 3.7, 2.2, 2.0) # with asbestosis
g <- factor(x = c(rep(1, length(x)),
                   rep(2, length(y)),
                   rep(3, length(z))),
             labels = c("ns", "oad", "a"))
dat <- data.frame(
   g = g,
   x = c(x, y, z))

## AD-Test
adKSampleTest(x ~ g, data = dat)

## BWS-Test
bwsKSampleTest(x ~ g, data = dat)

## Kruskal-Test
## Using incomplete beta approximation
kruskalTest(x ~ g, dat, dist="KruskalWallis")
## Using chisquare distribution
kruskalTest(x ~ g, dat, dist="Chisquare")

## Not run: 
## Check with kruskal.test from R stats
kruskal.test(x ~ g, dat)

## End(Not run)
## Using Conover's F
kruskalTest(x ~ g, dat, dist="FDist")

## Not run: 
## Check with aov on ranks
anova(aov(rank(x) ~ g, dat))
## Check with oneway.test
oneway.test(rank(x) ~ g, dat, var.equal = TRUE)

## End(Not run)

## Median Test asymptotic
medianTest(x ~ g, dat)

## Median Test with simulated p-values
set.seed(112)
medianTest(x ~ g, dat, simulate.p.value = TRUE)

Anderson-Darling Many-To-One Comparison Test

Description

Performs Anderson-Darling many-to-one comparison test.

Usage

adManyOneTest(x, ...)

## Default S3 method:
adManyOneTest(x, g, p.adjust.method = p.adjust.methods, ...)

## S3 method for class 'formula'
adManyOneTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in an one-factorial layout with non-normally distributed residuals Anderson-Darling's non-parametric test can be performed. Let there be k groups including the control, then the number of treatment levels is m = k - 1. Then m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: F_0 = F_i is tested in the two-tailed case against A_i: F_0 \ne F_i, ~~ (1 \le i \le m).

This function is a wrapper function that sequentially calls adKSampleTest for each pair. The calculated p-values for Pr(>|T2N|) can be adjusted to account for Type I error inflation using any method as implemented in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Scholz, F.W., Stephens, M.A. (1987) K-Sample Anderson-Darling Tests. Journal of the American Statistical Association 82, 918–924.

See Also

adKSampleTest, adAllPairsTest, ad.pval.

Examples

## Data set PlantGrowth
## Global test
adKSampleTest(weight ~ group, data = PlantGrowth)

##
ans <- adManyOneTest(weight ~ group,
                             data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

Algae Growth Inhibition Data Set

Description

A dose-response experiment was conducted using Atrazine at 9 different dose-levels including the zero-dose control and the biomass of algae (Selenastrum capricornutum) as the response variable. Three replicates were measured at day 0, 1 and 2. The fluorescence method (Mayer et al. 1997) was applied to measure biomass.

Format

A data frame with 22 observations on the following 10 variables.

concentration

a numeric vector of dose value in mg / L

Day.0

a numeric vector, total biomass

Day.0.1

a numeric vector, total biomass

Day.0.2

a numeric vector, total biomass

Day.1

a numeric vector, total biomass

Day.1.1

a numeric vector, total biomass

Day.1.2

a numeric vector, total biomass

Day.2

a numeric vector, total biomass

Day.2.1

a numeric vector, total biomass

Day.2.2

a numeric vector, total biomass

Source

ENV/JM/MONO(2006)18/ANN, page 24.

References

OECD (ed. 2006) Current approaches in the statistical analysis of ecotoxicity data: A guidance to application - Annexes, OECD Series on testing and assessment, No. 54, (ENV/JM/MONO(2006)18/ANN).

See Also

demo(algae)

Plotting PMCMR Objects

Description

Plots a bar-plot for objects of class "PMCMR".

Usage

barPlot(x, alpha = 0.05, ...)

Arguments

Value

A barplot where the height of the bars corresponds to the arithmetic mean. The extend of the whiskers are \pm z_{(1-\alpha/2)} \times s_{\mathrm{E},i}, where the latter denotes the standard error of the ith group. Symbolic letters are depicted on top of the bars, whereas different letters indicate significant differences between groups for the selected level of alpha.

Note

The barplot is strictly spoken only valid for normal data, as the depicted significance intervall implies symetry.

Examples

## data set chickwts
ans <- tukeyTest(weight ~ feed, data = chickwts)
barPlot(ans)

BWS All-Pairs Comparison Test

Description

Performs Baumgartner-Weiß-Schindler all-pairs comparison test.

Usage

bwsAllPairsTest(x, ...)

## Default S3 method:
bwsAllPairsTest(
  x,
  g,
  method = c("BWS", "Murakami"),
  p.adjust.method = p.adjust.methods,
  ...
)

## S3 method for class 'formula'
bwsAllPairsTest(
  formula,
  data,
  subset,
  na.action,
  method = c("BWS", "Murakami"),
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Baumgartner-Weiß-Schindler all-pairs comparison test can be used. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: F_i(x) = F_j(x) is tested in the two-tailed test against the alternative A_{ij}: F_i(x) \ne F_j(x), ~~ i \ne j.

This function is a wrapper function that sequentially calls bws_test for each pair. The default test method ("BWS") is the original Baumgartner-Weiß-Schindler test statistic B. For method == "Murakami" it is the modified BWS statistic denoted B*. The calculated p-values for Pr(>|B|) or Pr(>|B*|) can be adjusted to account for Type I error inflation using any method as implemented in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Baumgartner, W., Weiss, P., Schindler, H. (1998) A nonparametric test for the general two-sample problem, Biometrics 54, 1129–1135.

Murakami, H. (2006) K-sample rank test based on modified Baumgartner statistic and its power comparison, J. Jpn. Comp. Statist. 19, 1–13.

See Also

bws_test.

Examples

out <- bwsAllPairsTest(count ~ spray, InsectSprays, p.adjust="holm")
summary(out)
summaryGroup(out)

Murakami's k-Sample BWS Test

Description

Performs Murakami's k-Sample BWS Test.

Usage

bwsKSampleTest(x, ...)

## Default S3 method:
bwsKSampleTest(x, g, nperm = 1000, ...)

## S3 method for class 'formula'
bwsKSampleTest(formula, data, subset, na.action, nperm = 1000, ...)

Arguments

Details

Let X_{ij} ~ (1 \le i \le k,~ 1 \le 1 \le n_i) denote an identically and independently distributed variable that is obtained from an unknown continuous distribution F_i(x). Let R_{ij} be the rank of X_{ij}, where X_{ij} is jointly ranked from 1 to N, ~ N = \sum_{i=1}^k n_i. In the k-sample test the null hypothesis, H: F_i = F_j is tested against the alternative, A: F_i \ne F_j ~~(i \ne j) with at least one inequality beeing strict. Murakami (2006) has generalized the two-sample Baumgartner-Weiß-Schindler test (Baumgartner et al. 1998) and proposed a modified statistic B_k^* defined by

B_{k}^* = \frac{1}{k}\sum_{i=1}^k \left\{\frac{1}{n_i} \sum_{j=1}^{n_i} \frac{(R_{ij} - \mathsf{E}[R_{ij}])^2} {\mathsf{Var}[R_{ij}]}\right\},

where

\mathsf{E}[R_{ij}] = \frac{N + 1}{n_i + 1} j

and

\mathsf{Var}[R_{ij}] = \frac{j}{n_i + 1} \left(1 - \frac{j}{n_i + 1}\right) \frac{\left(N-n_i\right)\left(N+1\right)}{n_i + 2}.

The p-values are estimated via an assymptotic boot-strap method. It should be noted that the B_k^* detects both differences in the unknown location parameters and / or differences in the unknown scale parameters of the k-samples.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

One may increase the number of permutations to e.g. nperm = 10000 in order to get more precise p-values. However, this will be on the expense of computational time.

References

Baumgartner, W., Weiss, P., Schindler, H. (1998) A nonparametric test for the general two-sample problem, Biometrics 54, 1129–1135.

Murakami, H. (2006) K-sample rank test based on modified Baumgartner statistic and its power comparison, J. Jpn. Comp. Statist. 19, 1–13.

See Also

sample, bwsAllPairsTest, bwsManyOneTest.

Examples

## Hollander & Wolfe (1973), 116.
## Mucociliary efficiency from the rate of removal of dust in normal
## subjects, subjects with obstructive airway disease, and subjects
## with asbestosis.
x <- c(2.9, 3.0, 2.5, 2.6, 3.2) # normal subjects
y <- c(3.8, 2.7, 4.0, 2.4)      # with obstructive airway disease
z <- c(2.8, 3.4, 3.7, 2.2, 2.0) # with asbestosis
g <- factor(x = c(rep(1, length(x)),
                   rep(2, length(y)),
                   rep(3, length(z))),
             labels = c("ns", "oad", "a"))
dat <- data.frame(
   g = g,
   x = c(x, y, z))

## AD-Test
adKSampleTest(x ~ g, data = dat)

## BWS-Test
bwsKSampleTest(x ~ g, data = dat)

## Kruskal-Test
## Using incomplete beta approximation
kruskalTest(x ~ g, dat, dist="KruskalWallis")
## Using chisquare distribution
kruskalTest(x ~ g, dat, dist="Chisquare")

## Not run: 
## Check with kruskal.test from R stats
kruskal.test(x ~ g, dat)

## End(Not run)
## Using Conover's F
kruskalTest(x ~ g, dat, dist="FDist")

## Not run: 
## Check with aov on ranks
anova(aov(rank(x) ~ g, dat))
## Check with oneway.test
oneway.test(rank(x) ~ g, dat, var.equal = TRUE)

## End(Not run)

## Median Test asymptotic
medianTest(x ~ g, dat)

## Median Test with simulated p-values
set.seed(112)
medianTest(x ~ g, dat, simulate.p.value = TRUE)

BWS Many-To-One Comparison Test

Description

Performs Baumgartner-Weiß-Schindler many-to-one comparison test.

Usage

bwsManyOneTest(x, ...)

## Default S3 method:
bwsManyOneTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  method = c("BWS", "Murakami", "Neuhauser"),
  p.adjust.method = p.adjust.methods,
  ...
)

## S3 method for class 'formula'
bwsManyOneTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  method = c("BWS", "Murakami", "Neuhauser"),
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in an one-factorial layout with non-normally distributed residuals Baumgartner-Weiß-Schindler's non-parametric test can be performed. Let there be k groups including the control, then the number of treatment levels is m = k - 1. Then m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: F_0 = F_i is tested in the two-tailed case against A_i: F_0 \ne F_i, ~~ (1 \le i \le m).

This function is a wrapper function that sequentially calls bws_stat and bws_cdf for each pair. For the default test method ("BWS") the original Baumgartner-Weiß-Schindler test statistic B and its corresponding Pr(>|B|) is calculated. For method == "BWS" only a two-sided test is possible.

For method == "Murakami" the modified BWS statistic denoted B* and its corresponding Pr(>|B*|) is computed by sequentially calling murakami_stat and murakami_cdf. For method == "Murakami" only a two-sided test is possible.

If alternative == "greater" then the alternative, if one population is stochastically larger than the other is tested: H_i: F_0 = F_i against A_i: F_0 \ge F_i, ~~ (1 \le i \le m). The modified test-statistic B* according to Neuhäuser (2001) and its corresponding Pr(>B*) or Pr(<B*) is computed by sequentally calling murakami_stat and murakami_cdf with flavor = 2.

The p-values can be adjusted to account for Type I error inflation using any method as implemented in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Baumgartner, W., Weiss, P., Schindler, H. (1998) A nonparametric test for the general two-sample problem, Biometrics 54, 1129–1135.

Murakami, H. (2006) K-sample rank test based on modified Baumgartner statistic and its power comparison, J Jpn Comp Statist 19, 1–13.

Neuhäuser, M. (2001) One-Side Two-Sample and Trend Tests Based on a Modified Baumgartner-Weiss-Schindler Statistic. J Nonparametric Stat 13, 729–739.

See Also

murakami_stat, murakami_cdf, bws_stat, bws_cdf.

Examples

out <- bwsManyOneTest(weight ~ group, PlantGrowth, p.adjust="holm")
summary(out)

## A two-sample test
set.seed(1245)
x <- c(rnorm(20), rnorm(20,0.3))
g <- gl(2, 20)
summary(bwsManyOneTest(x ~ g, alternative = "less", p.adjust="none"))
summary(bwsManyOneTest(x ~ g, alternative = "greater", p.adjust="none"))

## Not run: 
## Check with the implementation in package BWStest
BWStest::bws_test(x=x[g==1], y=x[g==2], alternative = "less")
BWStest::bws_test(x=x[g==1], y=x[g==2], alternative = "greater")

## End(Not run)

Testing against Ordered Alternatives (Murakami's BWS Trend Test)

Description

Performs Murakami's modified Baumgartner-Weiß-Schindler test for testing against ordered alternatives.

Usage

bwsTrendTest(x, ...)

## Default S3 method:
bwsTrendTest(x, g, nperm = 1000, ...)

## S3 method for class 'formula'
bwsTrendTest(formula, data, subset, na.action, nperm = 1000, ...)

Arguments

Details

The null hypothesis, H_0: F_1(u) = F_2(u) = \ldots = F_k(u) ~~ u \in R is tested against a simple order hypothesis, H_\mathrm{A}: F_1(u) \le F_2(u) \le \ldots \le F_k(u),~F_1(u) < F_k(u), ~~ u \in R.

The p-values are estimated through an assymptotic boot-strap method using the function sample.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

One may increase the number of permutations to e.g. nperm = 10000 in order to get more precise p-values. However, this will be on the expense of computational time.

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Baumgartner, W., Weiss, P., Schindler, H. (1998) A nonparametric test for the general two-sample problem, Biometrics 54, 1129–1135.

Murakami, H. (2006) K-sample rank test based on modified Baumgartner statistic and its power comparison, J Jpn Comp Statist 19, 1–13.

Neuhäuser, M. (2001) One-Side Two-Sample and Trend Tests Based on a Modified Baumgartner-Weiss-Schindler Statistic. J Nonparametric Stat 13, 729–739.

See Also

sample, bwsAllPairsTest, bwsManyOneTest.

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

All-Pairs Comparisons for Simply Ordered Mean Ranksums

Description

Performs Nashimoto and Wright's all-pairs comparison procedure for simply ordered mean ranksums (NPT'-test and NPY'-test).

According to the authors, the procedure shall only be applied after Chacko's test (see chackoTest) indicates global significance.

Usage

chaAllPairsNashimotoTest(x, ...)

## Default S3 method:
chaAllPairsNashimotoTest(
  x,
  g,
  p.adjust.method = c(p.adjust.methods),
  alternative = c("greater", "less"),
  dist = c("Normal", "h"),
  ...
)

## S3 method for class 'formula'
chaAllPairsNashimotoTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = c(p.adjust.methods),
  alternative = c("greater", "less"),
  dist = c("Normal", "h"),
  ...
)

Arguments

Details

The modified procedure uses the property of a simple order, \theta_m' - \theta_m \le \theta_j - \theta_i \le \theta_l' - \theta_l \qquad (l \le i \le m~\mathrm{and}~ m' \le j \le l'). The null hypothesis H_{ij}: \theta_i = \theta_j is tested against the alternative A_{ij}: \theta_i < \theta_j for any 1 \le i < j \le k.

Let R_{ij} be the rank of X_{ij}, where X_{ij} is jointly ranked from \left\{1, 2, \ldots, N \right\}, ~~ N = \sum_{i=1}^k n_i, then the test statistics for all-pairs comparisons and a balanced design is calculated as

\hat{T}_{ij} = \max_{i \le m < m' \le j} \frac{\left(\bar{R}_{m'} - \bar{R}_m \right)} {\sigma_a / \sqrt{n}},

with n = n_i; ~ N = \sum_i^k n_i ~~ (1 \le i \le k), \bar{R}_i the mean rank for the ith group, and the expected variance (without ties) \sigma_a^2 = N \left(N + 1 \right) / 12.

For the NPY'-test (dist = "h"), if T_{ij} > h_{k-1,\alpha,\infty}.

For the unbalanced case with moderate imbalance the test statistic is

\hat{T}_{ij} = \max_{i \le m < m' \le j} \frac{\left(\bar{R}_{m'} - \bar{R}_m \right)} {\sigma_a \left(1/n_m + 1/n_{m'}\right)^{1/2}},

For the NPY'-test (dist="h") the null hypothesis is rejected in an unbalanced design, if \hat{T}_{ij} > h_{k,\alpha,\infty} / \sqrt{2}. In case of a NPY'-test, the function does not return p-values. Instead the critical h-values as given in the tables of Hayter (1990) for \alpha = 0.05 (one-sided) are looked up according to the number of groups (k-1) and the degree of freedoms (v = \infty).

For the NPT'-test (dist = "Normal"), the null hypothesis is rejected, if T_{ij} > \sqrt{2} t_{\alpha,\infty} = \sqrt{2} z_\alpha. Although Nashimoto and Wright (2005) originally did not use any p-adjustment, any method as available by p.adjust.methods can be selected for the adjustment of p-values estimated from the standard normal distribution.

Value

Either a list of class "osrt" if dist = "h" or a list of class "PMCMR" if dist = "Normal".

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated statistic(s)

crit.value

critical values for \alpha = 0.05.

alternative

a character string describing the alternative hypothesis.

parameter

the parameter(s) of the test distribution.

dist

a string that denotes the test distribution.

There are print and summary methods available.

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

The function will give a warning for the unbalanced case and returns the critical value h_{k-1,\alpha,\infty} / \sqrt{2} if applicable.

References

Hayter, A. J.(1990) A One-Sided Studentised Range Test for Testing Against a Simple Ordered Alternative, J Amer Stat Assoc 85, 778–785.

Nashimoto, K., Wright, F.T. (2007) Nonparametric Multiple-Comparison Methods for Simply Ordered Medians. Comput Stat Data Anal 51, 5068–5076.

See Also

Normal, chackoTest, NPMTest

Examples

## Example from Shirley (1977)
## Reaction times of mice to stimuli to their tails.
x <- c(2.4, 3, 3, 2.2, 2.2, 2.2, 2.2, 2.8, 2, 3,
 2.8, 2.2, 3.8, 9.4, 8.4, 3, 3.2, 4.4, 3.2, 7.4, 9.8, 3.2, 5.8,
 7.8, 2.6, 2.2, 6.2, 9.4, 7.8, 3.4, 7, 9.8, 9.4, 8.8, 8.8, 3.4,
 9, 8.4, 2.4, 7.8)
g <- gl(4, 10)

## Shirley's test
## one-sided test using look-up table
shirleyWilliamsTest(x ~ g, alternative = "greater")

## Chacko's global hypothesis test for 'greater'
chackoTest(x , g)

## post-hoc test, default is standard normal distribution (NPT'-test)
summary(chaAllPairsNashimotoTest(x, g, p.adjust.method = "none"))

## same but h-distribution (NPY'-test)
chaAllPairsNashimotoTest(x, g, dist = "h")

## NPM-test
NPMTest(x, g)

## Hayter-Stone test
hayterStoneTest(x, g)

## all-pairs comparisons
hsAllPairsTest(x, g)

Testing against Ordered Alternatives (Chacko's Test)

Description

Performs Chacko's test for testing against ordered alternatives.

Usage

chackoTest(x, ...)

## Default S3 method:
chackoTest(x, g, alternative = c("greater", "less"), ...)

## S3 method for class 'formula'
chackoTest(formula, data, subset, na.action, alternative = alternative, ...)

Arguments

Details

The null hypothesis, H_0: \theta_1 = \theta_2 = \ldots = \theta_k is tested against a simple order hypothesis, H_\mathrm{A}: \theta_1 \le \theta_2 \le \ldots \le \theta_k,~\theta_1 < \theta_k.

Let R_{ij} be the rank of X_{ij}, where X_{ij} is jointly ranked from \left\{1, 2, \ldots, N \right\}, ~~ N = \sum_{i=1}^k n_i, then the test statistic is calculated as

H = \frac{1}{\sigma_R^2} \sum_{i=1}^k n_i \left(\bar{R^*}_i - \bar{R}\right),

where \bar{R^*}_i is the isotonic mean of the i-th group and \sigma_R^2 = N \left(N + 1\right) / 12 the expected variance (without ties). H_0 is rejected, if H > \chi^2_{v,\alpha} with v = k -1 degree of freedom. The p-values are estimated from the chi-square distribution.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Source

The source code for the application of the pool adjacent violators theorem to calculate the isotonic means was taken from the file "pava.f", which is included in the package Iso:

Rolf Turner (2015). Iso: Functions to Perform Isotonic Regression. R package version 0.0-17. https://CRAN.R-project.org/package=Iso.

The file "pava.f" is a Ratfor modification of Algorithm AS 206.1:

Bril, G., Dykstra, R., Pillers, C., Robertson, T. (1984) Statistical Algorithms: Algorithm AS 206: Isotonic Regression in Two Independent Variables, Appl Statist 34, 352–357.

The Algorith AS 206 is available from StatLib https://lib.stat.cmu.edu/apstat/. The Royal Statistical Society holds the copyright to these routines, but has given its permission for their distribution provided that no fee is charged.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

The function does neither check nor correct for ties.

References

Chacko, V. J. (1963) Testing homogeneity against ordered alternatives, Ann Math Statist 34, 945–956.

See Also

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

Chen and Jan Many-to-One Comparisons Test

Description

Performs Chen and Jan nonparametric test for contrasting increasing (decreasing) dose levels of a treatment in a randomized block design.

Usage

chenJanTest(y, ...)

## Default S3 method:
chenJanTest(
  y,
  groups,
  blocks,
  alternative = c("greater", "less"),
  p.adjust.method = c("single-step", "SD1", p.adjust.methods),
  ...
)

Arguments

Details

Chen's test is a non-parametric step-down trend test for testing several treatment levels with a zero control. Let there be k groups including the control and let the zero dose level be indicated with i = 0 and the highest dose level with i = m, then the following m = k - 1 hypotheses are tested:

\begin{array}{ll} \mathrm{H}_{m}: \theta_0 = \theta_1 = \ldots = \theta_m, & \mathrm{A}_{m} = \theta_0 \le \theta_1 \le \ldots \theta_m, \theta_0 < \theta_m \\ \mathrm{H}_{m-1}: \theta_0 = \theta_1 = \ldots = \theta_{m-1}, & \mathrm{A}_{m-1} = \theta_0 \le \theta_1 \le \ldots \theta_{m-1}, \theta_0 < \theta_{m-1} \\ \vdots & \vdots \\ \mathrm{H}_{1}: \theta_0 = \theta_1, & \mathrm{A}_{1} = \theta_0 < \theta_1\\ \end{array}

Let Y_{ij1}, Y_{ij2}, \ldots, Y_{ijn_{ij}} (i = 1, 2, \dots, b, j = 0, 1, \ldots, k ~ \mathrm{and} ~ n_{ij} \geq 1) be a i.i.d. random variable of at least ordinal scale. Further,the zero dose control is indicated with j = 0.

The Mann-Whittney statistic is

T_{ij} = \sum_{u=0}^{j-1} \sum_{s=1}^{n_{ij}} \sum_{r=1}^{n_{iu}} I(Y_{ijs} - Y_{iur}), \qquad i = 1, 2, \ldots, b, ~ j = 1, 2, \ldots, k,

where where the indicator function returns I(a) = 1, ~ \mathrm{if}~ a > 0, 0.5 ~ \mathrm{if} a = 0 otherwise 0.

Let

N_{ij} = \sum_{s=0}^j n_{is} \qquad i = 1, 2, \ldots, b, ~ j = 1, 2, \ldots, k,

and

T_j = \sum_{i=1}^b T_{ij} \qquad j = 1, 2, \ldots, k.

The mean and variance of T_j are

\mu(T_j) = \sum_{i=1}^b n_{ij} ~ N_{ij-1} / 2 \qquad \mathrm{and}

\sigma(T_j) = \sum_{i=1}^b n_{ij} ~ N_{ij-1} \left[ \left(N_{ij} + 1\right) - \sum_{u=1}^{g_i} \left(t_u^3 - t_u \right) / \left\{N_{ij} \left(N_{ij} - 1\right) \right\} \right]/ 2,

with g_i the number of ties in the ith block and t_u the size of the tied group u.

The test statistic T_j^* is asymptotically multivariate normal distributed.

T_j^* = \frac{T_j - \mu(T_j)}{\sigma(T_j)}

If p.adjust.method = "single-step" than the p-values are calculated with the probability function of the multivariate normal distribution with \Sigma = I_k. Otherwise the standard normal distribution is used to calculate p-values and any method as available by p.adjust or by the step-down procedure as proposed by Chen (1999), if p.adjust.method = "SD1" can be used to account for \alpha-error inflation.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Chen, Y.I., Jan, S.L., 2002. Nonparametric Identification of the Minimum Effective Dose for Randomized Block Designs. Commun Stat-Simul Comput 31, 301–312.

See Also

Normal pmvnorm

Examples

## Example from Chen and Jan (2002, p. 306)
## MED is at dose level 2 (0.5 ppm SO2)
y <- c(0.2, 6.2, 0.3, 0.3, 4.9, 1.8, 3.9, 2, 0.3, 2.5, 5.4, 2.3, 12.7,
-0.2, 2.1, 6, 1.8, 3.9, 1.1, 3.8, 2.5, 1.3, -0.8, 13.1, 1.1,
12.8, 18.2, 3.4, 13.5, 4.4, 6.1, 2.8, 4, 10.6, 9, 4.2, 6.7, 35,
9, 12.9, 2, 7.1, 1.5, 10.6)
groups <- gl(4,11, labels = c("0", "0.25", "0.5", "1.0"))
blocks <- structure(rep(1:11, 4), class = "factor",
levels = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10", "11"))

summary(chenJanTest(y, groups, blocks, alternative = "greater"))
summary(chenJanTest(y, groups, blocks, alternative = "greater", p.adjust = "SD1"))

Chen's Many-to-One Comparisons Test

Description

Performs Chen's nonparametric test for contrasting increasing (decreasing) dose levels of a treatment.

Usage

chenTest(x, ...)

## Default S3 method:
chenTest(
  x,
  g,
  alternative = c("greater", "less"),
  p.adjust.method = c("SD1", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
chenTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  p.adjust.method = c("SD1", p.adjust.methods),
  ...
)

Arguments

Details

Chen's test is a non-parametric step-down trend test for testing several treatment levels with a zero control. Let X_{0j} denote a variable with the j-th realization of the control group (1 \le j \le n_0) and X_{ij} the j-the realization in the i-th treatment group (1 \le i \le k). The variables are i.i.d. of a least ordinal scale with F(x) = F(x_0) = F(x_i), ~ (1 \le i \le k). A total of m = k hypotheses can be tested:

\begin{array}{ll} \mathrm{H}_{m}: \theta_0 = \theta_1 = \ldots = \theta_m, & \mathrm{A}_{m} = \theta_0 \le \theta_1 \le \ldots \theta_m, \theta_0 < \theta_m \\ \mathrm{H}_{m-1}: \theta_0 = \theta_1 = \ldots = \theta_{m-1}, & \mathrm{A}_{m-1} = \theta_0 \le \theta_1 \le \ldots \theta_{m-1}, \theta_0 < \theta_{m-1} \\ \vdots & \vdots \\ \mathrm{H}_{1}: \theta_0 = \theta_1, & \mathrm{A}_{1} = \theta_0 < \theta_1\\ \end{array}

The statistics T_i are based on a Wilcoxon-type ranking:

T_i = \sum_{j=0}^{i=1} \sum_{u=1}^{n_i} \sum_{v=1}^{n_j} I(x_{iu} - x_{jv}), \qquad (1 \leq i \leq k),

where the indicator function returns I(a) = 1, ~ \mathrm{if}~ a > 0, 0.5 ~ \mathrm{if} a = 0 otherwise 0.

The expected ith mean is

\mu(T_i) = n_i N_{i-1} / 2,

with N_j = \sum_{j =0}^i n_j and the ith variance:

\sigma^2(T_i) = n_i N_{i-1} / 12 ~ \left\{N_i + 1 - \sum_{j=1}^g t_j \left(t_j^2 - 1 \right) / \left[N_i \left( N_i - 1 \right)\right]\right\}.

The test statistic T_i^* is asymptotically standard normal

T_i^* = \frac{T_i - \mu(T_i)} {\sqrt{\sigma^2(T_i)}}, \qquad (1 \leq i \leq k).

The p-values are calculated from the standard normal distribution. The p-values can be adjusted with any method as available by p.adjust or by the step-down procedure as proposed by Chen (1999), if p.adjust.method = "SD1".

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Chen, Y.-I., 1999, Nonparametric Identification of the Minimum Effective Dose. Biometrics 55, 1236–1240. doi:10.1111/j.0006-341X.1999.01236.x

See Also

wilcox.test, Normal

Examples

## Chen, 1999, p. 1237,
## Minimum effective dose (MED)
## is at 2nd dose level
df <- data.frame(x = c(23, 22, 14,
27, 23, 21,
28, 37, 35,
41, 37, 43,
28, 21, 30,
16, 19, 13),
g = gl(6, 3))
levels(df$g) <- 0:5
ans <- chenTest(x ~ g, data = df, alternative = "greater",
                p.adjust.method = "SD1")
summary(ans)

Cochran Test

Description

Performs Cochran's test for testing an outlying (or inlying) variance.

Usage

cochranTest(x, ...)

## Default S3 method:
cochranTest(x, g, alternative = c("greater", "less"), ...)

## S3 method for class 'formula'
cochranTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  ...
)

Arguments

Details

For normally distributed data the null hypothesis, H_0: \sigma_1^2 = \sigma_2^2 = \ldots = \sigma_k^2 is tested against the alternative (greater) H_{\mathrm{A}}: \sigma_p > \sigma_i ~~ (i \le k, i \ne p) with at least one inequality being strict.

The p-value is computed with the function pcochran.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Cochran, W.G. (1941) The distribution of the largest of a set of estimated variances as a fraction of their total. Ann. Eugen. 11, 47–52.

Wilrich, P.-T. (2011) Critical values of Mandel's h and k, Grubbs and the Cochran test statistic. Adv. Stat. Anal.. doi:10.1007/s10182-011-0185-y.

See Also

bartlett.test, fligner.test.

Examples

data(Pentosan)
cochranTest(value ~ lab, data = Pentosan, subset = (material == "A"))

Testing against Ordered Alternatives (Cuzick's Test)

Description

Performs Cuzick's test for testing against ordered alternatives.

Usage

cuzickTest(x, ...)

## Default S3 method:
cuzickTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  scores = NULL,
  continuity = FALSE,
  ...
)

## S3 method for class 'formula'
cuzickTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  scores = NULL,
  continuity = FALSE,
  ...
)

Arguments

Details

The null hypothesis, H_0: \theta_1 = \theta_2 = \ldots = \theta_k is tested against a simple order hypothesis, H_\mathrm{A}: \theta_1 \le \theta_2 \le \ldots \le \theta_k,~\theta_1 < \theta_k.

The p-values are estimated from the standard normal distribution.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Cuzick, J. (1995) A Wilcoxon-type test for trend, Statistics in Medicine 4, 87–90.

See Also

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

Grubbs Double Outlier Test

Description

Performs Grubbs double outlier test.

Usage

doubleGrubbsTest(x, alternative = c("two.sided", "greater", "less"), m = 10000)

Arguments

Details

Let X denote an identically and independently distributed continuous variate with realizations x_i ~~ (1 \le i \le k). Further, let the increasingly ordered realizations denote x_{(1)} \le x_{(2)} \le \ldots \le x_{(n)}. Then the following model for testing two maximum outliers can be proposed:

x_{(i)} = \left\{ \begin{array}{lcl} \mu + \epsilon_{(i)}, & \qquad & i = 1, \ldots, n - 2 \\ \mu + \Delta + \epsilon_{(j)} & \qquad & j = n-1, n \\ \end{array} \right.

with \epsilon \approx N(0,\sigma). The null hypothesis, H_0: \Delta = 0 is tested against the alternative, H_{\mathrm{A}}: \Delta > 0.

For testing two minimum outliers, the model can be proposed as

x_{(i)} = \left\{ \begin{array}{lcl} \mu + \Delta + \epsilon_{(j)} & \qquad & j = 1, 2 \\ \mu + \epsilon_{(i)}, & \qquad & i = 3, \ldots, n \\ \end{array} \right.

The null hypothesis is tested against the alternative, H_{\mathrm{A}}: \Delta < 0.

The p-value is computed with the function pdgrubbs.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Grubbs, F. E. (1950) Sample criteria for testing outlying observations. Ann. Math. Stat. 21, 27–58.

Wilrich, P.-T. (2011) Critical values of Mandel's h and k, Grubbs and the Cochran test statistic. Adv. Stat. Anal.. doi:10.1007/s10182-011-0185-y.

Examples

data(Pentosan)
dat <- subset(Pentosan, subset = (material == "A"))
labMeans <- tapply(dat$value, dat$lab, mean)
doubleGrubbsTest(x = labMeans, alternative = "less")

Multiple Comparisons of Mean Rank Sums

Description

Performs the all-pairs comparison test for different factor levels according to Dwass, Steel, Critchlow and Fligner.

Usage

dscfAllPairsTest(x, ...)

## Default S3 method:
dscfAllPairsTest(x, g, ...)

## S3 method for class 'formula'
dscfAllPairsTest(formula, data, subset, na.action, ...)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals the DSCF all-pairs comparison test can be used. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: F_i(x) = F_j(x) is tested in the two-tailed test against the alternative A_{ij}: F_i(x) \ne F_j(x), ~~ i \ne j. As opposed to the all-pairs comparison procedures that depend on Kruskal ranks, the DSCF test is basically an extension of the U-test as re-ranking is conducted for each pairwise test.

The p-values are estimated from the studentized range distriburtion.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Douglas, C. E., Fligner, A. M. (1991) On distribution-free multiple comparisons in the one-way analysis of variance, Communications in Statistics - Theory and Methods 20, 127–139.

Dwass, M. (1960) Some k-sample rank-order tests. In Contributions to Probability and Statistics, Edited by: I. Olkin, Stanford: Stanford University Press.

Steel, R. G. D. (1960) A rank sum test for comparing all pairs of treatments, Technometrics 2, 197–207

See Also

Tukey, pairwise.wilcox.test

Duncan's Multiple Range Test

Description

Performs Duncan's all-pairs comparisons test for normally distributed data with equal group variances.

Usage

duncanTest(x, ...)

## Default S3 method:
duncanTest(x, g, ...)

## S3 method for class 'formula'
duncanTest(formula, data, subset, na.action, ...)

## S3 method for class 'aov'
duncanTest(x, ...)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with normally distributed residuals and equal variances Duncan's multiple range test can be performed. Let X_{ij} denote a continuous random variable with the j-the realization (1 \le j \le n_i) in the i-th group (1 \le i \le k). Furthermore, the total sample size is N = \sum_{i=1}^k n_i. A total of m = k(k-1)/2 hypotheses can be tested: The null hypothesis is H_{ij}: \mu_i = \mu_j ~~ (i \ne j) is tested against the alternative A_{ij}: \mu_i \ne \mu_j (two-tailed). Duncan's all-pairs test statistics are given by

t_{(i)(j)} \frac{\bar{X}_{(i)} - \bar{X}_{(j)}} {s_{\mathrm{in}} \left(r\right)^{1/2}}, ~~ (i < j)

with s^2_{\mathrm{in}} the within-group ANOVA variance, r = k / \sum_{i=1}^k n_i and \bar{X}_{(i)} the increasingly ordered means 1 \le i \le k. The null hypothesis is rejected if

\mathrm{Pr} \left\{ |t_{(i)(j)}| \ge q_{vm'\alpha'} | \mathrm{H} \right\}_{(i)(j)} = \alpha' = \min \left\{1,~ 1 - (1 - \alpha)^{(1 / (m' - 1))} \right\},

with v = N - k degree of freedom, the range m' = 1 + |i - j| and \alpha' the Bonferroni adjusted alpha-error. The p-values are computed from the Tukey distribution.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Duncan, D. B. (1955) Multiple range and multiple F tests, Biometrics 11, 1–42.

See Also

Tukey, TukeyHSD tukeyTest

Examples

fit <- aov(weight ~ feed, chickwts)
shapiro.test(residuals(fit))
bartlett.test(weight ~ feed, chickwts)
anova(fit)

## also works with fitted objects of class aov
res <- duncanTest(fit)
summary(res)
summaryGroup(res)

Dunnett's T3 Test

Description

Performs Dunnett's all-pairs comparison test for normally distributed data with unequal variances.

Usage

dunnettT3Test(x, ...)

## Default S3 method:
dunnettT3Test(x, g, ...)

## S3 method for class 'formula'
dunnettT3Test(formula, data, subset, na.action, ...)

## S3 method for class 'aov'
dunnettT3Test(x, ...)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with normally distributed residuals but unequal groups variances the T3 test of Dunnett can be performed. Let X_{ij} denote a continuous random variable with the j-the realization (1 \le j \le n_i) in the i-th group (1 \le i \le k). Furthermore, the total sample size is N = \sum_{i=1}^k n_i. A total of m = k(k-1)/2 hypotheses can be tested: The null hypothesis is H_{ij}: \mu_i = \mu_j ~~ (i \ne j) is tested against the alternative A_{ij}: \mu_i \ne \mu_j (two-tailed). Dunnett T3 all-pairs test statistics are given by

t_{ij} \frac{\bar{X}_i - \bar{X_j}} {\left( s^2_j / n_j + s^2_i / n_i \right)^{1/2}}, ~~ (i \ne j)

with s^2_i the variance of the i-th group. The null hypothesis is rejected (two-tailed) if

\mathrm{Pr} \left\{ |t_{ij}| \ge T_{v_{ij}\rho_{ij}\alpha'/2} | \mathrm{H} \right\}_{ij} = \alpha,

with Welch's approximate solution for calculating the degree of freedom.

v_{ij} = \frac{\left( s^2_i / n_i + s^2_j / n_j \right)^2} {s^4_i / n^2_i \left(n_i - 1\right) + s^4_j / n^2_j \left(n_j - 1\right)}.

The p-values are computed from the studentized maximum modulus distribution that is the equivalent of the multivariate t distribution with \rho_{ii} = 1, ~ \rho_{ij} = 0 ~ (i \ne j). The function pmvt is used to calculate the p-values.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

C. W. Dunnett (1980) Pair wise multiple comparisons in the unequal variance case, Journal of the American Statistical Association 75, 796–800.

See Also

pmvt

Examples

fit <- aov(weight ~ feed, chickwts)
shapiro.test(residuals(fit))
bartlett.test(weight ~ feed, chickwts)
anova(fit)

## also works with fitted objects of class aov
res <- dunnettT3Test(fit)
summary(res)
summaryGroup(res)

Dunnett's Many-to-One Comparisons Test

Description

Performs Dunnett's multiple comparisons test with one control.

Usage

dunnettTest(x, ...)

## Default S3 method:
dunnettTest(x, g, alternative = c("two.sided", "greater", "less"), ...)

## S3 method for class 'formula'
dunnettTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  ...
)

## S3 method for class 'aov'
dunnettTest(x, alternative = c("two.sided", "greater", "less"), ...)

Arguments

Details

For many-to-one comparisons in an one-factorial layout with normally distributed residuals Dunnett's test can be used. Let X_{0j} denote a continuous random variable with the j-the realization of the control group (1 \le j \le n_0) and X_{ij} the j-the realization in the i-th treatment group (1 \le i \le k). Furthermore, the total sample size is N = n_0 + \sum_{i=1}^k n_i. A total of m = k hypotheses can be tested: The null hypothesis is H_{i}: \mu_i = \mu_0 is tested against the alternative A_{i}: \mu_i \ne \mu_0 (two-tailed). Dunnett's test statistics are given by

t_{i} \frac{\bar{X}_i - \bar{X_0}} {s_{\mathrm{in}} \left(1/n_0 + 1/n_i\right)^{1/2}}, ~~ (1 \le i \le k)

with s^2_{\mathrm{in}} the within-group ANOVA variance. The null hypothesis is rejected if |t_{ij}| > |T_{kv\rho\alpha}| (two-tailed), with v = N - k degree of freedom and rho the correlation:

\rho_{ij} = \sqrt{\frac{n_i n_j} {\left(n_i + n_0\right) \left(n_j+ n_0\right)}} ~~ (i \ne j) .

The p-values are computed with the function pDunnett that is a wrapper to the the multivariate-t distribution as implemented in the function pmvt.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Dunnett, C. W. (1955) A multiple comparison procedure for comparing several treatments with a control. Journal of the American Statistical Association 50, 1096–1121.

OECD (ed. 2006) Current approaches in the statistical analysis of ecotoxicity data: A guidance to application - Annexes. OECD Series on testing and assessment, No. 54.

See Also

pmvt pDunnett

Examples

fit <- aov(Y ~ DOSE, data = trout)
shapiro.test(residuals(fit))
bartlett.test(Y ~ DOSE, data = trout)

## works with fitted object of class aov
summary(dunnettTest(fit, alternative = "less"))

All-Pairs Comparisons Test for Balanced Incomplete Block Designs

Description

Performs Conover-Iman all-pairs comparison test for a balanced incomplete block design (BIBD).

Usage

durbinAllPairsTest(y, ...)

## Default S3 method:
durbinAllPairsTest(y, groups, blocks, p.adjust.method = p.adjust.methods, ...)

Arguments

Details

For all-pairs comparisons in a balanced incomplete block design the proposed test of Conover and Imam can be applied. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: \theta_i = \theta_j is tested in the two-tailed test against the alternative A_{ij}: \theta_i \ne \theta_j, ~~ i \ne j.

The p-values are computed from the t distribution. If no p-value adjustment is performed (p.adjust.method = "none"), than a simple protected test is recommended, i.e. the all-pairs comparisons should only be applied after a significant durbinTest. However, any method as implemented in p.adjust.methods can be selected by the user.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Conover, W. J., Iman, R. L. (1979) On multiple-comparisons procedures, Tech. Rep. LA-7677-MS, Los Alamos Scientific Laboratory.

Conover, W. J. (1999) Practical nonparametric Statistics, 3rd. Edition, Wiley.

See Also

durbinTest

Examples

## Example for an incomplete block design:
## Data from Conover (1999, p. 391).
y <- matrix(c(2,NA,NA,NA,3, NA,  3,  3,  3, NA, NA, NA,  3, NA, NA,
  1,  2, NA, NA, NA,  1,  1, NA,  1,  1,
NA, NA, NA, NA,  2, NA,  2,  1, NA, NA, NA, NA,
 3, NA,  2,  1, NA, NA, NA, NA,  3, NA,  2,  2),
ncol=7, nrow=7, byrow=FALSE, dimnames=list(1:7, LETTERS[1:7]))
durbinAllPairsTest(y)

Durbin Test

Description

Performs Durbin's tests whether k groups (or treatments) in a two-way balanced incomplete block design (BIBD) have identical effects.

Usage

durbinTest(y, ...)

## Default S3 method:
durbinTest(y, groups, blocks, ...)

Arguments

Details

For testing a two factorial layout of a balanced incomplete block design whether the k groups have identical effects, the Durbin test can be performed. The null hypothesis, H_0: \theta_i = \theta_j ~ (1 \le i < j \le k), is tested against the alternative that at least one \theta_i \ne \theta_j.

The p-values are computed from the chi-square distribution.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

The function does not test, whether it is a true BIBD. This function does not test for ties.

References

Conover,W. J. (1999) Practical nonparametric Statistics, 3rd. Edition, Wiley.

Heckert, N. A., Filliben, J. J. (2003) NIST Handbook 148: Dataplot Reference Manual, Volume 2: Let Subcommands and Library Functions. National Institute of Standards and Technology Handbook Series, June 2003.

Examples

## Example for an incomplete block design:
## Data from Conover (1999, p. 391).
y <- matrix(c(
2,NA,NA,NA,3, NA,  3,  3,  3, NA, NA, NA,  3, NA, NA,
  1,  2, NA, NA, NA,  1,  1, NA,  1,  1,
NA, NA, NA, NA,  2, NA,  2,  1, NA, NA, NA, NA,
 3, NA,  2,  1, NA, NA, NA, NA,  3, NA,  2,  2
), ncol=7, nrow=7, byrow=FALSE,
dimnames=list(1:7, LETTERS[1:7]))
durbinTest(y)

Testing Several Treatments With One Control

Description

Performs Fligner-Wolfe non-parametric test for simultaneous testing of several locations of treatment groups against the location of the control group.

Usage

flignerWolfeTest(x, ...)

## Default S3 method:
flignerWolfeTest(
  x,
  g,
  alternative = c("greater", "less"),
  dist = c("Wilcoxon", "Normal"),
  ...
)

## S3 method for class 'formula'
flignerWolfeTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  dist = c("Wilcoxon", "Normal"),
  ...
)

Arguments

Details

For a one-factorial layout with non-normally distributed residuals the Fligner-Wolfe test can be used.

Let there be k-1-treatment groups and one control group, then the null hypothesis, H_0: \theta_i - \theta_c = 0 ~ (1 \le i \le k-1) is tested against the alternative (greater), A_1: \theta_i - \theta_c > 0 ~ (1 \le i \le k-1), with at least one inequality being strict.

Let n_c denote the sample size of the control group, N^t = \sum_{i=1}^{k-1} n_i the sum of all treatment sample sizes and N = N^t + n_c. The test statistic without taken ties into account is

W = \sum_{j=1}^{k-1} \sum_{i=1}^{n_i} r_{ij} - \frac{N^t \left(N^t + 1 \right) }{2}

with r_{ij} the rank of variable x_{ij}. The null hypothesis is rejected, if W > W_{\alpha,m,n} with m = N^t and n = n_c.

In the presence of ties, the statistic is

\hat{z} = \frac{W - n_c N^t / 2}{s_W},

where

s_W = \frac{n_c N^t}{12 N \left(N - 1 \right)} \sum_{j=1}^g t_j \left(t_j^2 - 1\right),

with g the number of tied groups and t_j the number of tied values in the jth group. The null hypothesis is rejected, if \hat{z} > z_\alpha (as cited in EPA 2006).

If dist = Wilcoxon, then the p-values are estimated from the Wilcoxon distribution, else the Normal distribution is used. The latter can be used, if ties are present.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

EPA (2006) Data Quality Assessment: Statistical Methods for Practitioners (Guideline No. EPA QA/G-9S), US-EPA.

Fligner, M.A., Wolfe, D.A. (1982) Distribution-free tests for comparing several treatments with a control. Stat Neerl 36, 119–127.

See Also

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

Conover's All-Pairs Comparisons Test for Unreplicated Blocked Data

Description

Performs Conover's all-pairs comparisons tests of Friedman-type ranked data.

Usage

frdAllPairsConoverTest(y, ...)

## Default S3 method:
frdAllPairsConoverTest(
  y,
  groups,
  blocks,
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

For all-pairs comparisons in a two factorial unreplicated complete block design with non-normally distributed residuals, Conover's test can be performed on Friedman-type ranked data.

A total of m = k ( k -1 )/2 hypotheses can be tested. The null hypothesis, H_{ij}: \theta_i = \theta_j, is tested in the two-tailed case against the alternative, A_{ij}: \theta_i \ne \theta_j, ~~ i \ne j.

If p.adjust.method == "single-step" the p-values are computed from the studentized range distribution. Otherwise, the p-values are computed from the t-distribution using any of the p-adjustment methods as included in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Conover, W. J., Iman, R. L. (1979) On multiple-comparisons procedures, Tech. Rep. LA-7677-MS, Los Alamos Scientific Laboratory.

Conover, W. J. (1999) Practical nonparametric Statistics, 3rd. Edition, Wiley.

See Also

friedmanTest, friedman.test, frdAllPairsExactTest, frdAllPairsMillerTest, frdAllPairsNemenyiTest, frdAllPairsSiegelTest

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ##
 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))
 print(y)
 friedmanTest(y)

 ## Eisinga et al. 2017
 frdAllPairsExactTest(y=y, p.adjust = "bonferroni")

 ## Conover's test
 frdAllPairsConoverTest(y=y, p.adjust = "bonferroni")

 ## Nemenyi's test
 frdAllPairsNemenyiTest(y=y)

 ## Miller et al.
 frdAllPairsMillerTest(y=y)

 ## Siegel-Castellan
 frdAllPairsSiegelTest(y=y, p.adjust = "bonferroni")

 ## Irrelevant of group order?
 x <- as.vector(y)
 g <- rep(colnames(y), each = length(x)/length(colnames(y)))
 b <- rep(rownames(y), times = length(x)/length(rownames(y)))
 xDF <- data.frame(x, g, b) # grouped by colnames

 frdAllPairsNemenyiTest(xDF$x, groups = xDF$g, blocks = xDF$b)
 o <- order(xDF$b) # order per block increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])
 o <- order(xDF$x) # order per value increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])

 ## formula method (only works for Nemenyi)
 frdAllPairsNemenyiTest(x ~ g | b, data = xDF)

Exact All-Pairs Comparisons Test for Unreplicated Blocked Data

Description

Performs exact all-pairs comparisons tests of Friedman-type ranked data according to Eisinga et al. (2017).

Usage

frdAllPairsExactTest(y, ...)

## Default S3 method:
frdAllPairsExactTest(
  y,
  groups,
  blocks,
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For all-pairs comparisons in a two factorial unreplicated complete block design with non-normally distributed residuals, an exact test can be performed on Friedman-type ranked data.

A total of m = k ( k -1 )/2 hypotheses can be tested. The null hypothesis, H_{ij}: \theta_i = \theta_j, is tested in the two-tailed case against the alternative, A_{ij}: \theta_i \ne \theta_j, ~~ i \ne j.

The exact p-values are computed using the code of "pexactfrsd.R" that was a supplement to the publication of Eisinga et al. (2017). Additionally, any of the p-adjustment methods as included in p.adjust can be selected, for p-value adjustment.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Source

The function frdAllPairsExactTest uses the code of the file pexactfrsd.R that was a supplement to:

R. Eisinga, T. Heskes, B. Pelzer, M. Te Grotenhuis (2017), Exact p-values for Pairwise Comparison of Friedman Rank Sums, with Application to Comparing Classifiers, BMC Bioinformatics, 18:68.

References

Eisinga, R., Heskes, T., Pelzer, B., Te Grotenhuis, M. (2017) Exact p-values for Pairwise Comparison of Friedman Rank Sums, with Application to Comparing Classifiers, BMC Bioinformatics, 18:68.

See Also

friedmanTest, friedman.test, frdAllPairsConoverTest, frdAllPairsMillerTest, frdAllPairsNemenyiTest, frdAllPairsSiegelTest

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ##
 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))
 print(y)
 friedmanTest(y)

 ## Eisinga et al. 2017
 frdAllPairsExactTest(y=y, p.adjust = "bonferroni")

 ## Conover's test
 frdAllPairsConoverTest(y=y, p.adjust = "bonferroni")

 ## Nemenyi's test
 frdAllPairsNemenyiTest(y=y)

 ## Miller et al.
 frdAllPairsMillerTest(y=y)

 ## Siegel-Castellan
 frdAllPairsSiegelTest(y=y, p.adjust = "bonferroni")

 ## Irrelevant of group order?
 x <- as.vector(y)
 g <- rep(colnames(y), each = length(x)/length(colnames(y)))
 b <- rep(rownames(y), times = length(x)/length(rownames(y)))
 xDF <- data.frame(x, g, b) # grouped by colnames

 frdAllPairsNemenyiTest(xDF$x, groups = xDF$g, blocks = xDF$b)
 o <- order(xDF$b) # order per block increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])
 o <- order(xDF$x) # order per value increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])

 ## formula method (only works for Nemenyi)
 frdAllPairsNemenyiTest(x ~ g | b, data = xDF)

Millers's All-Pairs Comparisons Test for Unreplicated Blocked Data

Description

Performs Miller's all-pairs comparisons tests of Friedman-type ranked data.

Usage

frdAllPairsMillerTest(y, ...)

## Default S3 method:
frdAllPairsMillerTest(y, groups, blocks, ...)

Arguments

Details

For all-pairs comparisons in a two factorial unreplicated complete block design with non-normally distributed residuals, Miller's test can be performed on Friedman-type ranked data.

A total of m = k ( k -1 )/2 hypotheses can be tested. The null hypothesis, H_{ij}: \theta_i = \theta_j, is tested in the two-tailed case against the alternative, A_{ij}: \theta_i \ne \theta_j, ~~ i \ne j.

The p-values are computed from the chi-square distribution.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Bortz J., Lienert, G. A., Boehnke, K. (1990) Verteilungsfreie Methoden in der Biostatistik. Berlin: Springer.

Miller Jr., R. G. (1996) Simultaneous statistical inference. New York: McGraw-Hill.

Wike, E. L. (2006), Data Analysis. A Statistical Primer for Psychology Students. New Brunswick: Aldine Transaction.

See Also

friedmanTest, friedman.test, frdAllPairsExactTest, frdAllPairsConoverTest, frdAllPairsNemenyiTest, frdAllPairsSiegelTest

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ##
 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))
 print(y)
 friedmanTest(y)

 ## Eisinga et al. 2017
 frdAllPairsExactTest(y=y, p.adjust = "bonferroni")

 ## Conover's test
 frdAllPairsConoverTest(y=y, p.adjust = "bonferroni")

 ## Nemenyi's test
 frdAllPairsNemenyiTest(y=y)

 ## Miller et al.
 frdAllPairsMillerTest(y=y)

 ## Siegel-Castellan
 frdAllPairsSiegelTest(y=y, p.adjust = "bonferroni")

 ## Irrelevant of group order?
 x <- as.vector(y)
 g <- rep(colnames(y), each = length(x)/length(colnames(y)))
 b <- rep(rownames(y), times = length(x)/length(rownames(y)))
 xDF <- data.frame(x, g, b) # grouped by colnames

 frdAllPairsNemenyiTest(xDF$x, groups = xDF$g, blocks = xDF$b)
 o <- order(xDF$b) # order per block increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])
 o <- order(xDF$x) # order per value increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])

 ## formula method (only works for Nemenyi)
 frdAllPairsNemenyiTest(x ~ g | b, data = xDF)

Nemenyi's All-Pairs Comparisons Test for Unreplicated Blocked Data

Description

Performs Nemenyi's all-pairs comparisons tests of Friedman-type ranked data.

Usage

frdAllPairsNemenyiTest(y, ...)

## Default S3 method:
frdAllPairsNemenyiTest(y, groups, blocks, ...)

## S3 method for class 'formula'
frdAllPairsNemenyiTest(formula, data, subset, na.action, ...)

Arguments

Details

For all-pairs comparisons in a two factorial unreplicated complete block design with non-normally distributed residuals, Nemenyi's test can be performed on Friedman-type ranked data.

A total of m = k ( k -1 )/2 hypotheses can be tested. The null hypothesis, H_{ij}: \theta_i = \theta_j, is tested in the two-tailed case against the alternative, A_{ij}: \theta_i \ne \theta_j, ~~ i \ne j.

The p-values are computed from the studentized range distribution.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Demsar, J. (2006) Statistical comparisons of classifiers over multiple data sets, Journal of Machine Learning Research 7, 1–30.

Miller Jr., R. G. (1996) Simultaneous statistical inference. New York: McGraw-Hill.

Nemenyi, P. (1963), Distribution-free Multiple Comparisons. Ph.D. thesis, Princeton University.

Sachs, L. (1997) Angewandte Statistik. Berlin: Springer.

See Also

friedmanTest, friedman.test, frdAllPairsExactTest, frdAllPairsConoverTest, frdAllPairsMillerTest, frdAllPairsSiegelTest

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ##
 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))
 print(y)
 friedmanTest(y)

 ## Eisinga et al. 2017
 frdAllPairsExactTest(y=y, p.adjust = "bonferroni")

 ## Conover's test
 frdAllPairsConoverTest(y=y, p.adjust = "bonferroni")

 ## Nemenyi's test
 frdAllPairsNemenyiTest(y=y)

 ## Miller et al.
 frdAllPairsMillerTest(y=y)

 ## Siegel-Castellan
 frdAllPairsSiegelTest(y=y, p.adjust = "bonferroni")

 ## Irrelevant of group order?
 x <- as.vector(y)
 g <- rep(colnames(y), each = length(x)/length(colnames(y)))
 b <- rep(rownames(y), times = length(x)/length(rownames(y)))
 xDF <- data.frame(x, g, b) # grouped by colnames

 frdAllPairsNemenyiTest(xDF$x, groups = xDF$g, blocks = xDF$b)
 o <- order(xDF$b) # order per block increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])
 o <- order(xDF$x) # order per value increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])

 ## formula method (only works for Nemenyi)
 frdAllPairsNemenyiTest(x ~ g | b, data = xDF)

Siegel and Castellan's All-Pairs Comparisons Test for Unreplicated Blocked Data

Description

Performs Siegel and Castellan's all-pairs comparisons tests of Friedman-type ranked data.

Usage

frdAllPairsSiegelTest(y, ...)

## Default S3 method:
frdAllPairsSiegelTest(
  y,
  groups,
  blocks,
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For all-pairs comparisons in a two factorial unreplicated complete block design with non-normally distributed residuals, Siegel and Castellan's test can be performed on Friedman-type ranked data.

A total of m = k ( k -1 )/2 hypotheses can be tested. The null hypothesis, H_{ij}: \theta_i = \theta_j, is tested in the two-tailed case against the alternative, A_{ij}: \theta_i \ne \theta_j, ~~ i \ne j.

The p-values are computed from the standard normal distribution. Any method as implemented in p.adjust can be used for p-value adjustment.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Siegel, S., Castellan Jr., N. J. (1988) Nonparametric Statistics for the Behavioral Sciences. 2nd ed. New York: McGraw-Hill.

See Also

friedmanTest, friedman.test, frdAllPairsExactTest, frdAllPairsConoverTest, frdAllPairsNemenyiTest, frdAllPairsMillerTest

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ##
 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))
 print(y)
 friedmanTest(y)

 ## Eisinga et al. 2017
 frdAllPairsExactTest(y=y, p.adjust = "bonferroni")

 ## Conover's test
 frdAllPairsConoverTest(y=y, p.adjust = "bonferroni")

 ## Nemenyi's test
 frdAllPairsNemenyiTest(y=y)

 ## Miller et al.
 frdAllPairsMillerTest(y=y)

 ## Siegel-Castellan
 frdAllPairsSiegelTest(y=y, p.adjust = "bonferroni")

 ## Irrelevant of group order?
 x <- as.vector(y)
 g <- rep(colnames(y), each = length(x)/length(colnames(y)))
 b <- rep(rownames(y), times = length(x)/length(rownames(y)))
 xDF <- data.frame(x, g, b) # grouped by colnames

 frdAllPairsNemenyiTest(xDF$x, groups = xDF$g, blocks = xDF$b)
 o <- order(xDF$b) # order per block increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])
 o <- order(xDF$x) # order per value increasingly
 frdAllPairsNemenyiTest(xDF$x[o], groups = xDF$g[o], blocks = xDF$b[o])

 ## formula method (only works for Nemenyi)
 frdAllPairsNemenyiTest(x ~ g | b, data = xDF)

House Test

Description

Performs House nonparametric equivalent of William's test for contrasting increasing dose levels of a treatment in a complete randomized block design.

Usage

frdHouseTest(y, ...)

## Default S3 method:
frdHouseTest(y, groups, blocks, alternative = c("greater", "less"), ...)

Arguments

Details

House test is a non-parametric step-down trend test for testing several treatment levels with a zero control. Let there be k groups including the control and let the zero dose level be indicated with i = 0 and the highest dose level with i = m, then the following m = k - 1 hypotheses are tested:

\begin{array}{ll} \mathrm{H}_{m}: \theta_0 = \theta_1 = \ldots = \theta_m, & \mathrm{A}_{m} = \theta_0 \le \theta_1 \le \ldots \theta_m, \theta_0 < \theta_m \\ \mathrm{H}_{m-1}: \theta_0 = \theta_1 = \ldots = \theta_{m-1}, & \mathrm{A}_{m-1} = \theta_0 \le \theta_1 \le \ldots \theta_{m-1}, \theta_0 < \theta_{m-1} \\ \vdots & \vdots \\ \mathrm{H}_{1}: \theta_0 = \theta_1, & \mathrm{A}_{1} = \theta_0 < \theta_1\\ \end{array}

Let Y_{ij} ~ (1 \leq i \leq n, 0 \leq j \leq k) be a i.i.d. random variable of at least ordinal scale. Further, let \bar{R}_0,~\bar{R}_1, \ldots,~\bar{R}_k be Friedman's average ranks and set \bar{R}_0^*, \leq \ldots \leq \bar{R}_k^* to be its isotonic regression estimators under the order restriction \theta_0 \leq \ldots \leq \theta_k.

The statistics is

T_j = \left(\bar{R}_j^* - \bar{R}_0 \right)~ \left[ \left(V_j - H_j \right) \left(2 / n \right) \right]^{-1/2} \qquad (1 \leq j \leq k),

with

V_j = \left(j + 1\right) ~ \left(j + 2 \right) / 12

and

H_j = \left(t^3 - t \right) / \left(12 j n \right),

where t is the number of tied ranks.

The critical t'_{i,v,\alpha}-values as given in the tables of Williams (1972) for \alpha = 0.05 (one-sided) are looked up according to the degree of freedoms (v = \infty) and the order number of the dose level (j).

For the comparison of the first dose level (j = 1) with the control, the critical z-value from the standard normal distribution is used (Normal).

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Chen, Y.-I., 1999. Rank-Based Tests for Dose Finding in Nonmonotonic Dose–Response Settings. Biometrics 55, 1258–1262. doi:10.1111/j.0006-341X.1999.01258.x

House, D.E., 1986. A Nonparametric Version of Williams’ Test for Randomized Block Design. Biometrics 42, 187–190.

See Also

friedmanTest, friedman.test, frdManyOneExactTest, frdManyOneDemsarTest

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ## Assume A is the control.

 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))

 ## Global Friedman test
 friedmanTest(y)

 ## Demsar's many-one test
 summary(frdManyOneDemsarTest(y=y, p.adjust = "bonferroni",
                      alternative = "greater"))

 ## Exact many-one test
 summary(frdManyOneExactTest(y=y, p.adjust = "bonferroni",
                     alternative = "greater"))

 ## Nemenyi's many-one test
 summary(frdManyOneNemenyiTest(y=y, alternative = "greater"))

 ## House test
 frdHouseTest(y, alternative = "greater")

Demsar's Many-to-One Test for Unreplicated Blocked Data

Description

Performs Demsar's non-parametric many-to-one comparison test for Friedman-type ranked data.

Usage

frdManyOneDemsarTest(y, ...)

## Default S3 method:
frdManyOneDemsarTest(
  y,
  groups,
  blocks,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in a two factorial unreplicated complete block design with non-normally distributed residuals, Demsar's test can be performed on Friedman-type ranked data.

Let there be k groups including the control, then the number of treatment levels is m = k - 1. A total of m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: \theta_0 = \theta_i is tested in the two-tailed case against A_i: \theta_0 \ne \theta_i, ~~ (1 \le i \le m).

The p-values are computed from the standard normal distribution. Any of the p-adjustment methods as included in p.adjust can be used for the adjustment of p-values.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Demsar, J. (2006) Statistical comparisons of classifiers over multiple data sets, Journal of Machine Learning Research 7, 1–30.

See Also

friedmanTest, friedman.test, frdManyOneExactTest, frdManyOneNemenyiTest.

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ## Assume A is the control.

 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))

 ## Global Friedman test
 friedmanTest(y)

 ## Demsar's many-one test
 summary(frdManyOneDemsarTest(y=y, p.adjust = "bonferroni",
                      alternative = "greater"))

 ## Exact many-one test
 summary(frdManyOneExactTest(y=y, p.adjust = "bonferroni",
                     alternative = "greater"))

 ## Nemenyi's many-one test
 summary(frdManyOneNemenyiTest(y=y, alternative = "greater"))

 ## House test
 frdHouseTest(y, alternative = "greater")

Exact Many-to-One Test for Unreplicated Blocked Data

Description

Performs an exact non-parametric many-to-one comparison test for Friedman-type ranked data according to Eisinga et al. (2017).

Usage

frdManyOneExactTest(y, ...)

## Default S3 method:
frdManyOneExactTest(y, groups, blocks, p.adjust.method = p.adjust.methods, ...)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in a two factorial unreplicated complete block design with non-normally distributed residuals, an exact test can be performed on Friedman-type ranked data.

Let there be k groups including the control, then the number of treatment levels is m = k - 1. A total of m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: \theta_0 = \theta_i is tested in the two-tailed case against A_i: \theta_0 \ne \theta_i, ~~ (1 \le i \le m).

The exact p-values are computed using the code of "pexactfrsd.R" that was a supplement to the publication of Eisinga et al. (2017). Additionally, any of the p-adjustment methods as included in p.adjust can be selected, for p-value adjustment.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Eisinga, R., Heskes, T., Pelzer, B., Te Grotenhuis, M. (2017) Exact p-values for Pairwise Comparison of Friedman Rank Sums, with Application to Comparing Classifiers, BMC Bioinformatics, 18:68.

See Also

friedmanTest, friedman.test, frdManyOneDemsarTest, frdManyOneNemenyiTest.

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ## Assume A is the control.

 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))

 ## Global Friedman test
 friedmanTest(y)

 ## Demsar's many-one test
 summary(frdManyOneDemsarTest(y=y, p.adjust = "bonferroni",
                      alternative = "greater"))

 ## Exact many-one test
 summary(frdManyOneExactTest(y=y, p.adjust = "bonferroni",
                     alternative = "greater"))

 ## Nemenyi's many-one test
 summary(frdManyOneNemenyiTest(y=y, alternative = "greater"))

 ## House test
 frdHouseTest(y, alternative = "greater")

Nemenyi's Many-to-One Test for Unreplicated Blocked Data

Description

Performs Nemenyi's non-parametric many-to-one comparison test for Friedman-type ranked data.

Usage

frdManyOneNemenyiTest(y, ...)

## Default S3 method:
frdManyOneNemenyiTest(
  y,
  groups,
  blocks,
  alternative = c("two.sided", "greater", "less"),
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in a two factorial unreplicated complete block design with non-normally distributed residuals, Nemenyi's test can be performed on Friedman-type ranked data.

Let there be k groups including the control, then the number of treatment levels is m = k - 1. A total of m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: \theta_0 = \theta_i is tested in the two-tailed case against A_i: \theta_0 \ne \theta_i, ~~ (1 \le i \le m).

The p-values are computed from the multivariate normal distribution. As pmvnorm applies a numerical method, the estimated p-values are seet depended.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Hollander, M., Wolfe, D. A., Chicken, E. (2014), Nonparametric Statistical Methods. 3rd ed. New York: Wiley. 2014.

Miller Jr., R. G. (1996), Simultaneous Statistical Inference. New York: McGraw-Hill.

Nemenyi, P. (1963), Distribution-free Multiple Comparisons. Ph.D. thesis, Princeton University.

Siegel, S., Castellan Jr., N. J. (1988), Nonparametric Statistics for the Behavioral Sciences. 2nd ed. New York: McGraw-Hill.

Zarr, J. H. (1999), Biostatistical Analysis. 4th ed. Upper Saddle River: Prentice-Hall.

See Also

friedmanTest, friedman.test, frdManyOneExactTest, frdManyOneDemsarTest pmvnorm, set.seed

Examples

 ## Sachs, 1997, p. 675
 ## Six persons (block) received six different diuretics
 ## (A to F, treatment).
 ## The responses are the Na-concentration (mval)
 ## in the urine measured 2 hours after each treatment.
 ## Assume A is the control.

 y <- matrix(c(
 3.88, 5.64, 5.76, 4.25, 5.91, 4.33, 30.58, 30.14, 16.92,
 23.19, 26.74, 10.91, 25.24, 33.52, 25.45, 18.85, 20.45,
 26.67, 4.44, 7.94, 4.04, 4.4, 4.23, 4.36, 29.41, 30.72,
 32.92, 28.23, 23.35, 12, 38.87, 33.12, 39.15, 28.06, 38.23,
 26.65),nrow=6, ncol=6,
 dimnames=list(1:6, LETTERS[1:6]))

 ## Global Friedman test
 friedmanTest(y)

 ## Demsar's many-one test
 summary(frdManyOneDemsarTest(y=y, p.adjust = "bonferroni",
                      alternative = "greater"))

 ## Exact many-one test
 summary(frdManyOneExactTest(y=y, p.adjust = "bonferroni",
                     alternative = "greater"))

 ## Nemenyi's many-one test
 summary(frdManyOneNemenyiTest(y=y, alternative = "greater"))

 ## House test
 frdHouseTest(y, alternative = "greater")

Friedman Rank Sum Test

Description

Performs a Friedman rank sum test. The null hypothesis H_0: \theta_i = \theta_j~~(i \ne j) is tested against the alternative H_{\mathrm{A}}: \theta_i \ne \theta_j, with at least one inequality beeing strict.

Usage

friedmanTest(y, ...)

## Default S3 method:
friedmanTest(y, groups, blocks, dist = c("Chisquare", "FDist"), ...)

Arguments

Details

The function has implemented Friedman's test as well as the extension of Conover anf Iman (1981). Friedman's test statistic is assymptotically chi-squared distributed. Consequently, the default test distribution is dist = "Chisquare".

If dist = "FDist" is selected, than the approach of Conover and Imam (1981) is performed. The Friedman Test using the F-distribution leads to the same results as doing an two-way Analysis of Variance without interaction on rank transformed data.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Conover, W.J., Iman, R.L. (1981) Rank Transformations as a Bridge Between Parametric and Nonparametric Statistics. Am Stat 35, 124–129.

Sachs, L. (1997) Angewandte Statistik. Berlin: Springer.

See Also

friedman.test

Examples

## Hollander & Wolfe (1973), p. 140ff.
## Comparison of three methods ("round out", "narrow angle", and
##  "wide angle") for rounding first base.  For each of 18 players
##  and the three method, the average time of two runs from a point on
##  the first base line 35ft from home plate to a point 15ft short of
##  second base is recorded.
RoundingTimes <-
matrix(c(5.40, 5.50, 5.55,
        5.85, 5.70, 5.75,
        5.20, 5.60, 5.50,
        5.55, 5.50, 5.40,
        5.90, 5.85, 5.70,
        5.45, 5.55, 5.60,
        5.40, 5.40, 5.35,
        5.45, 5.50, 5.35,
        5.25, 5.15, 5.00,
        5.85, 5.80, 5.70,
        5.25, 5.20, 5.10,
        5.65, 5.55, 5.45,
        5.60, 5.35, 5.45,
        5.05, 5.00, 4.95,
        5.50, 5.50, 5.40,
        5.45, 5.55, 5.50,
        5.55, 5.55, 5.35,
        5.45, 5.50, 5.55,
        5.50, 5.45, 5.25,
        5.65, 5.60, 5.40,
        5.70, 5.65, 5.55,
        6.30, 6.30, 6.25),
      nrow = 22,
      byrow = TRUE,
      dimnames = list(1 : 22,
                      c("Round Out", "Narrow Angle", "Wide Angle")))

## Chisquare distribution
friedmanTest(RoundingTimes)

## check with friedman.test from R stats
friedman.test(RoundingTimes)

## F-distribution
friedmanTest(RoundingTimes, dist = "FDist")

## Check with One-way repeated measure ANOVA
rmat <- RoundingTimes
for (i in 1:length(RoundingTimes[,1])) rmat[i,] <- rank(rmat[i,])
dataf <- data.frame(
    y = y <- as.vector(rmat),
    g = g <- factor(c(col(RoundingTimes))),
    b = b <- factor(c(row(RoundingTimes))))
summary(aov(y ~ g + Error(b), data = dataf))

Games-Howell Test

Description

Performs Games-Howell all-pairs comparison test for normally distributed data with unequal group variances.

Usage

gamesHowellTest(x, ...)

## Default S3 method:
gamesHowellTest(x, g, ...)

## S3 method for class 'formula'
gamesHowellTest(formula, data, subset, na.action, ...)

## S3 method for class 'aov'
gamesHowellTest(x, ...)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with normally distributed residuals but unequal between-groups variances the Games-Howell Test can be performed. Let X_{ij} denote a continuous random variable with the j-the realization (1 \le j \le n_i) in the i-th group (1 \le i \le k). Furthermore, the total sample size is N = \sum_{i=1}^k n_i. A total of m = k(k-1)/2 hypotheses can be tested: The null hypothesis is H_{ij}: \mu_i = \mu_j ~~ (i \ne j) is tested against the alternative A_{ij}: \mu_i \ne \mu_j (two-tailed). Games-Howell Test all-pairs test statistics are given by

t_{ij} \frac{\bar{X}_i - \bar{X_j}} {\left( s^2_j / n_j + s^2_i / n_i \right)^{1/2}}, ~~ (i \ne j)

with s^2_i the variance of the i-th group. The null hypothesis is rejected (two-tailed) if

\mathrm{Pr} \left\{ |t_{ij}| \sqrt{2} \ge q_{m v_{ij} \alpha} | \mathrm{H} \right\}_{ij} = \alpha,

with Welch's approximate solution for calculating the degree of freedom.

v_{ij} = \frac{\left( s^2_i / n_i + s^2_j / n_j \right)^2} {s^4_i / n^2_i \left(n_i - 1\right) + s^4_j / n^2_j \left(n_j - 1\right)}.

The p-values are computed from the Tukey distribution.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

See Also

Tukey

Examples

fit <- aov(weight ~ feed, chickwts)
shapiro.test(residuals(fit))
bartlett.test(weight ~ feed, chickwts) # var1 = varN
anova(fit)

## also works with fitted objects of class aov
res <- gamesHowellTest(fit)
summary(res)
summaryGroup(res)

Generalized Extreme Studentized Deviate Many-Outlier Test

Description

Performs Rosner's generalized extreme studentized deviate procedure to detect up-to maxr outliers in a univariate sample that follows an approximately normal distribution.

Usage

gesdTest(x, maxr)

Arguments

References

Rosner, B. (1983) Percentage Points for a Generalized ESD Many-Outlier Procedure, Technometrics 25, 165–172.

Examples

## Taken from Rosner (1983):
x <- c(-0.25,0.68,0.94,1.15,1.20,1.26,1.26,
1.34,1.38,1.43,1.49,1.49,1.55,1.56,
1.58,1.65,1.69,1.70,1.76,1.77,1.81,
1.91,1.94,1.96,1.99,2.06,2.09,2.10,
2.14,2.15,2.23,2.24,2.26,2.35,2.37,
2.40,2.47,2.54,2.62,2.64,2.90,2.92,
2.92,2.93,3.21,3.26,3.30,3.59,3.68,
4.30,4.64,5.34,5.42,6.01)

out <- gesdTest(x, 10)

## print method
out

## summary method
summary(out)

Gore Test

Description

Performs Gore's test. The null hypothesis H_0: \theta_i = \theta_j~~(i \ne j) is tested against the alternative H_{\mathrm{A}}: \theta_i \ne \theta_j, with at least one inequality beeing strict.

Usage

goreTest(y, groups, blocks)

Arguments

Details

The function has implemented Gore's test for testing main effects in unbalanced CRB designs, i.e. there are one ore more observations per cell. The statistic is assymptotically chi-squared distributed.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Gore, A. P. (1975) Some nonparametric tests and selection procedures for main effects in two-way layouts. Ann. Inst. Stat. Math. 27, 487–500.

See Also

friedmanTest, skillingsMackTest, durbinTest

Examples

## Crop Yield of 3 varieties on two
## soil classes
X <-c("130,A,Light
115,A,Light
123,A,Light
142,A,Light
117,A,Heavy
125,A,Heavy
139,A,Heavy
108,B,Light
114,B,Light
124,B,Light
106,B,Light
91,B,Heavy
111,B,Heavy
110,B,Heavy
155,C,Light
146,C,Light
151,C,Light
165,C,Light
97,C,Heavy
108,C,Heavy")
con <- textConnection(X)
x <- read.table(con, header=FALSE, sep=",")
close(con)
colnames(x) <- c("Yield", "Variety", "SoilType")
goreTest(y = x$Yield, groups = x$Variety, blocks = x$SoilType)

Grubbs Outlier Test

Description

Performs Grubbs single outlier test.

Usage

grubbsTest(x, alternative = c("two.sided", "greater", "less"))

Arguments

Details

Let X denote an identically and independently distributed continuous variate with realizations x_i ~~ (1 \le i \le k). Further, let the increasingly ordered realizations denote x_{(1)} \le x_{(2)} \le \ldots \le x_{(n)}. Then the following model for a single maximum outlier can be proposed:

x_{(i)} = \left\{ \begin{array}{lcl} \mu + \epsilon_{(i)}, & \qquad & i = 1, \ldots, n - 1 \\ \mu + \Delta + \epsilon_{(n)} & & \\ \end{array} \right.

with \epsilon \approx N(0,\sigma). The null hypothesis, H_0: \Delta = 0 is tested against the alternative, H_{\mathrm{A}}: \Delta > 0.

For testing a single minimum outlier, the model can be proposed as

x_{(i)} = \left\{ \begin{array}{lcl} \mu + \Delta + \epsilon_{(1)} & & \\ \mu + \epsilon_{(i)}, & \qquad & i = 2, \ldots, n \\ \end{array} \right.

The null hypothesis is tested against the alternative, H_{\mathrm{A}}: \Delta < 0.

The p-value is computed with the function pgrubbs.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Grubbs, F. E. (1950) Sample criteria for testing outlying observations. Ann. Math. Stat. 21, 27–58.

Wilrich, P.-T. (2011) Critical values of Mandel's h and k, Grubbs and the Cochran test statistic. Adv. Stat. Anal.. doi:10.1007/s10182-011-0185-y.

Examples

data(Pentosan)
dat <- subset(Pentosan, subset = (material == "A"))
labMeans <- tapply(dat$value, dat$lab, mean)
grubbsTest(x = labMeans, alternative = "two.sided")

Hartley's Maximum F-Ratio Test of Homogeneity of Variances

Description

Performs Hartley's maximum F-ratio test of the null that variances in each of the groups (samples) are the same.

Usage

hartleyTest(x, ...)

## Default S3 method:
hartleyTest(x, g, ...)

## S3 method for class 'formula'
hartleyTest(formula, data, subset, na.action, ...)

Arguments

Details

If x is a list, its elements are taken as the samples to be compared for homogeneity of variances. In this case, the elements must all be numeric data vectors, g is ignored, and one can simply use hartleyTest(x) to perform the test. If the samples are not yet contained in a list, use hartleyTest(list(x, ...)).

Otherwise, x must be a numeric data vector, and g must be a vector or factor object of the same length as x giving the group for the corresponding elements of x.

Hartley's parametric test requires normality and a nearly balanced design. The p-value of the test is calculated with the function pmaxFratio of the package SuppDists.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Hartley, H.O. (1950) The maximum F-ratio as a short cut test for heterogeneity of variance, Biometrika 37, 308–312.

See Also

bartlett.test, pmaxFratio

Examples

hartleyTest(count ~ spray, data = InsectSprays)

Hayter-Stone Test

Description

Performs the non-parametric Hayter-Stone procedure to test against an monotonically increasing alternative.

Usage

hayterStoneTest(x, ...)

## Default S3 method:
hayterStoneTest(
  x,
  g,
  alternative = c("greater", "less"),
  method = c("look-up", "boot", "asympt"),
  nperm = 10000,
  ...
)

## S3 method for class 'formula'
hayterStoneTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  method = c("look-up", "boot", "asympt"),
  nperm = 10000,
  ...
)

Arguments

Details

Let X be an identically and idepentendly distributed variable that was n times observed at k increasing treatment levels. Hayter and Stone (1991) proposed a non-parametric procedure to test the null hypothesis, H: \theta_i = \theta_j ~~ (i < j \le k) against a simple order alternative, A: \theta_i < \theta_j, with at least one inequality being strict.

The statistic for a global test is calculated as,

h = \max_{1 \le i < j \le k} \frac{2 \sqrt{6} \left(U_{ij} - n_i n_j / 2 \right)} {\sqrt{n_i n_j \left(n_i + n_j + 1 \right)}},

with the Mann-Whittney counts:

U_{ij} = \sum_{a=1}^{n_i} \sum_{b=1}^{n_j} I\left\{x_{ia} < x_{ja}\right\}.

Under the large sample approximation, the test statistic h is distributed as h_{k,\alpha,v}. Thus, the null hypothesis is rejected, if h > h_{k,\alpha,v}, with v = \infty degree of freedom.

If method = "look-up" the function will not return p-values. Instead the critical h-values as given in the tables of Hayter (1990) for \alpha = 0.05 (one-sided) are looked up according to the number of groups (k) and the degree of freedoms (v = \infty).

If method = "boot" an asymptotic permutation test is conducted and a p-value is returned.

If method = "asympt" is selected the asymptotic p-value is estimated as implemented in the function pHayStonLSA of the package NSM3.

Value

Either a list of class htest or a list with class "osrt" that contains the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated statistic(s)

crit.value

critical values for \alpha = 0.05.

alternative

a character string describing the alternative hypothesis.

parameter

the parameter(s) of the test distribution.

dist

a string that denotes the test distribution.

There are print and summary methods available.

Source

If method = "asympt" is selected, this function calls an internal probability function pHS. The GPL-2 code for this function was taken from pHayStonLSA of the the package NSM3:

Grant Schneider, Eric Chicken and Rachel Becvarik (2020) NSM3: Functions and Datasets to Accompany Hollander, Wolfe, and Chicken - Nonparametric Statistical Methods, Third Edition. R package version 1.15. https://CRAN.R-project.org/package=NSM3

References

Hayter, A. J.(1990) A One-Sided Studentised Range Test for Testing Against a Simple Ordered Alternative, J Amer Stat Assoc 85, 778–785.

Hayter, A.J., Stone, G. (1991) Distribution free multiple comparisons for monotonically ordered treatment effects. Austral J Statist 33, 335–346.

See Also

osrtTest, hsAllPairsTest, sample, pHayStonLSA

Examples

## Example from Shirley (1977)
## Reaction times of mice to stimuli to their tails.
x <- c(2.4, 3, 3, 2.2, 2.2, 2.2, 2.2, 2.8, 2, 3,
 2.8, 2.2, 3.8, 9.4, 8.4, 3, 3.2, 4.4, 3.2, 7.4, 9.8, 3.2, 5.8,
 7.8, 2.6, 2.2, 6.2, 9.4, 7.8, 3.4, 7, 9.8, 9.4, 8.8, 8.8, 3.4,
 9, 8.4, 2.4, 7.8)
g <- gl(4, 10)

## Shirley's test
## one-sided test using look-up table
shirleyWilliamsTest(x ~ g, alternative = "greater")

## Chacko's global hypothesis test for 'greater'
chackoTest(x , g)

## post-hoc test, default is standard normal distribution (NPT'-test)
summary(chaAllPairsNashimotoTest(x, g, p.adjust.method = "none"))

## same but h-distribution (NPY'-test)
chaAllPairsNashimotoTest(x, g, dist = "h")

## NPM-test
NPMTest(x, g)

## Hayter-Stone test
hayterStoneTest(x, g)

## all-pairs comparisons
hsAllPairsTest(x, g)

Hayter-Stone All-Pairs Comparison Test

Description

Performs the non-parametric Hayter-Stone all-pairs procedure to test against monotonically increasing alternatives.

Usage

hsAllPairsTest(x, ...)

## Default S3 method:
hsAllPairsTest(
  x,
  g,
  alternative = c("greater", "less"),
  method = c("look-up", "boot", "asympt"),
  nperm = 10000,
  ...
)

## S3 method for class 'formula'
hsAllPairsTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("greater", "less"),
  method = c("look-up", "boot", "asympt"),
  nperm = 10000,
  ...
)

Arguments

Details

Let X be an identically and idepentendly distributed variable that was n times observed at k increasing treatment levels. Hayter and Stone (1991) proposed a non-parametric procedure to test the null hypothesis, H: \theta_i = \theta_j ~~ (i < j \le k) against a simple order alternative, A: \theta_i < \theta_j.

The statistic for all-pairs comparisons is calculated as,

S_{ij} = \frac{2 \sqrt{6} \left(U_{ij} - n_i n_j / 2 \right)} {\sqrt{n_i n_j \left(n_i + n_j + 1 \right)}},

with the Mann-Whittney counts:

U_{ij} = \sum_{a=1}^{n_i} \sum_{b=1}^{n_j} I\left\{x_{ia} < x_{ja}\right\}.

Under the large sample approximation, the test statistic S_{ij} is distributed as h_{k,\alpha,v}. Thus, the null hypothesis is rejected, if S_{ij} > h_{k,\alpha,v}, with v = \infty degree of freedom.

If method = "look-up" the function will not return p-values. Instead the critical h-values as given in the tables of Hayter (1990) for \alpha = 0.05 (one-sided) are looked up according to the number of groups (k) and the degree of freedoms (v = \infty).

If method = "boot" an asymetric permutation test is conducted and p-values are returned.

If method = "asympt" is selected the asymptotic p-value is estimated as implemented in the function pHayStonLSA of the package NSM3.

Value

Either a list of class "PMCMR" or a list with class "osrt" that contains the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated statistic(s)

crit.value

critical values for \alpha = 0.05.

alternative

a character string describing the alternative hypothesis.

parameter

the parameter(s) of the test distribution.

dist

a string that denotes the test distribution.

There are print and summary methods available.

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Source

If method = "asympt" is selected, this function calls an internal probability function pHS. The GPL-2 code for this function was taken from pHayStonLSA of the the package NSM3:

Grant Schneider, Eric Chicken and Rachel Becvarik (2020) NSM3: Functions and Datasets to Accompany Hollander, Wolfe, and Chicken - Nonparametric Statistical Methods, Third Edition. R package version 1.15. https://CRAN.R-project.org/package=NSM3

References

Hayter, A. J.(1990) A One-Sided Studentised Range Test for Testing Against a Simple Ordered Alternative, Journal of the American Statistical Association 85, 778–785.

Hayter, A.J., Stone, G. (1991) Distribution free multiple comparisons for monotonically ordered treatment effects. Austral J Statist 33, 335–346.

See Also

hayterStoneTest sample

Examples

## Example from Shirley (1977)
## Reaction times of mice to stimuli to their tails.
x <- c(2.4, 3, 3, 2.2, 2.2, 2.2, 2.2, 2.8, 2, 3,
 2.8, 2.2, 3.8, 9.4, 8.4, 3, 3.2, 4.4, 3.2, 7.4, 9.8, 3.2, 5.8,
 7.8, 2.6, 2.2, 6.2, 9.4, 7.8, 3.4, 7, 9.8, 9.4, 8.8, 8.8, 3.4,
 9, 8.4, 2.4, 7.8)
g <- gl(4, 10)

## Shirley's test
## one-sided test using look-up table
shirleyWilliamsTest(x ~ g, alternative = "greater")

## Chacko's global hypothesis test for 'greater'
chackoTest(x , g)

## post-hoc test, default is standard normal distribution (NPT'-test)
summary(chaAllPairsNashimotoTest(x, g, p.adjust.method = "none"))

## same but h-distribution (NPY'-test)
chaAllPairsNashimotoTest(x, g, dist = "h")

## NPM-test
NPMTest(x, g)

## Hayter-Stone test
hayterStoneTest(x, g)

## all-pairs comparisons
hsAllPairsTest(x, g)

Testing against Ordered Alternatives (Johnson-Mehrotra Test)

Description

Performs the Johnson-Mehrotra test for testing against ordered alternatives in a balanced one-factorial sampling design.

Usage

johnsonTest(x, ...)

## Default S3 method:
johnsonTest(x, g, alternative = c("two.sided", "greater", "less"), ...)

## S3 method for class 'formula'
johnsonTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  ...
)

Arguments

Details

The null hypothesis, H_0: \theta_1 = \theta_2 = \ldots = \theta_k is tested against a simple order hypothesis, H_\mathrm{A}: \theta_1 \le \theta_2 \le \ldots \le \theta_k,~\theta_1 < \theta_k.

The p-values are estimated from the standard normal distribution.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Bortz, J. (1993). Statistik für Sozialwissenschaftler (4th ed.). Berlin: Springer.

Johnson, R. A., Mehrotra, K. G. (1972) Some c-sample nonparametric tests for ordered alternatives. Journal of the Indian Statistical Association 9, 8–23.

See Also

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

Testing against Ordered Alternatives (Jonckheere-Terpstra Test)

Description

Performs the Jonckheere-Terpstra test for testing against ordered alternatives.

Usage

jonckheereTest(x, ...)

## Default S3 method:
jonckheereTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  continuity = FALSE,
  ...
)

## S3 method for class 'formula'
jonckheereTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  continuity = FALSE,
  ...
)

Arguments

Details

The null hypothesis, H_0: \theta_1 = \theta_2 = \ldots = \theta_k is tested against a simple order hypothesis, H_\mathrm{A}: \theta_1 \le \theta_2 \le \ldots \le \theta_k,~\theta_1 < \theta_k.

The p-values are estimated from the standard normal distribution.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Source

The code for the computation of the standard deviation for the Jonckheere-Terpstra test in the presence of ties was taken from:

Kloke, J., McKean, J. (2016) npsm: Package for Nonparametric Statistical Methods using R. R package version 0.5. https://CRAN.R-project.org/package=npsm

Note

jonckheereTest(x, g, alternative = "two.sided", continuity = TRUE) is equivalent to

cor.test(x, as.numeric(g), method = "kendall", alternative = "two.sided", continuity = TRUE)

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Jonckheere, A. R. (1954) A distribution-free k-sample test against ordered alternatives. Biometrica 41, 133–145.

Kloke, J., McKean, J. W. (2015) Nonparametric statistical methods using R. Boca Raton, FL: Chapman & Hall/CRC.

See Also

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

Kruskal-Wallis Rank Sum Test

Description

Performs a Kruskal-Wallis rank sum test.

Usage

kruskalTest(x, ...)

## Default S3 method:
kruskalTest(x, g, dist = c("Chisquare", "KruskalWallis", "FDist"), ...)

## S3 method for class 'formula'
kruskalTest(
  formula,
  data,
  subset,
  na.action,
  dist = c("Chisquare", "KruskalWallis", "FDist"),
  ...
)

Arguments

Details

For one-factorial designs with non-normally distributed residuals the Kruskal-Wallis rank sum test can be performed to test the H_0: F_1(x) = F_2(x) = \ldots = F_k(x) against the H_\mathrm{A}: F_i (x) \ne F_j(x)~ (i \ne j) with at least one strict inequality.

Let R_{ij} be the joint rank of X_{ij}, with R_{(1)(1)} = 1, \ldots, R_{(n)(n)} = N, ~~ N = \sum_{i=1}^k n_i, The test statistic is calculated as

H = \sum_{i=1}^k n_i \left(\bar{R}_i - \bar{R}\right) / \sigma_R,

with the mean rank of the i-th group

\bar{R}_i = \sum_{j = 1}^{n_{i}} R_{ij} / n_i,

the expected value

\bar{R} = \left(N +1\right) / 2

and the expected variance as

\sigma_R^2 = N \left(N + 1\right) / 12.

In case of ties the statistic H is divided by \left(1 - \sum_{i=1}^r t_i^3 - t_i \right) / \left(N^3 - N\right)

According to Conover and Imam (1981), the statistic H is related to the F-quantile as

F = \frac{H / \left(k - 1\right)} {\left(N - 1 - H\right) / \left(N - k\right)}

which is equivalent to a one-way ANOVA F-test using rank transformed data (see examples).

The function provides three different dist for p-value estimation:

Chisquare

p-values are computed from the Chisquare distribution with v = k - 1 degree of freedom.

KruskalWallis

p-values are computed from the pKruskalWallis of the package SuppDists.

FDist

p-values are computed from the FDist distribution with v_1 = k-1, ~ v_2 = N -k degree of freedom.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

References

Conover, W.J., Iman, R.L. (1981) Rank Transformations as a Bridge Between Parametric and Nonparametric Statistics. Am Stat 35, 124–129.

Kruskal, W.H., Wallis, W.A. (1952) Use of Ranks in One-Criterion Variance Analysis. J Am Stat Assoc 47, 583–621.

Sachs, L. (1997) Angewandte Statistik. Berlin: Springer.

See Also

kruskal.test, pKruskalWallis, Chisquare, FDist

Examples

## Hollander & Wolfe (1973), 116.
## Mucociliary efficiency from the rate of removal of dust in normal
## subjects, subjects with obstructive airway disease, and subjects
## with asbestosis.
x <- c(2.9, 3.0, 2.5, 2.6, 3.2) # normal subjects
y <- c(3.8, 2.7, 4.0, 2.4)      # with obstructive airway disease
z <- c(2.8, 3.4, 3.7, 2.2, 2.0) # with asbestosis
g <- factor(x = c(rep(1, length(x)),
                   rep(2, length(y)),
                   rep(3, length(z))),
             labels = c("ns", "oad", "a"))
dat <- data.frame(
   g = g,
   x = c(x, y, z))

## AD-Test
adKSampleTest(x ~ g, data = dat)

## BWS-Test
bwsKSampleTest(x ~ g, data = dat)

## Kruskal-Test
## Using incomplete beta approximation
kruskalTest(x ~ g, dat, dist="KruskalWallis")
## Using chisquare distribution
kruskalTest(x ~ g, dat, dist="Chisquare")

## Not run: 
## Check with kruskal.test from R stats
kruskal.test(x ~ g, dat)

## End(Not run)
## Using Conover's F
kruskalTest(x ~ g, dat, dist="FDist")

## Not run: 
## Check with aov on ranks
anova(aov(rank(x) ~ g, dat))
## Check with oneway.test
oneway.test(rank(x) ~ g, dat, var.equal = TRUE)

## End(Not run)

## Median Test asymptotic
medianTest(x ~ g, dat)

## Median Test with simulated p-values
set.seed(112)
medianTest(x ~ g, dat, simulate.p.value = TRUE)

Conover's All-Pairs Rank Comparison Test

Description

Performs Conover's non-parametric all-pairs comparison test for Kruskal-type ranked data.

Usage

kwAllPairsConoverTest(x, ...)

## Default S3 method:
kwAllPairsConoverTest(
  x,
  g,
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
kwAllPairsConoverTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Conover's non-parametric test can be performed. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: \mu_i(x) = \mu_j(x) is tested in the two-tailed test against the alternative A_{ij}: \mu_i(x) \ne \mu_j(x), ~~ i \ne j.

If p.adjust.method == "single-step" the p-values are computed from the studentized range distribution. Otherwise, the p-values are computed from the t-distribution using any of the p-adjustment methods as included in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Conover, W. J, Iman, R. L. (1979) On multiple-comparisons procedures, Tech. Rep. LA-7677-MS, Los Alamos Scientific Laboratory.

See Also

Tukey, TDist, p.adjust, kruskalTest, kwAllPairsDunnTest, kwAllPairsNemenyiTest

Examples

## Data set InsectSprays
## Global test
kruskalTest(count ~ spray, data = InsectSprays)

## Conover's all-pairs comparison test
## single-step means Tukey's p-adjustment
ans <- kwAllPairsConoverTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "single-step")
summary(ans)

## Dunn's all-pairs comparison test
ans <- kwAllPairsDunnTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "bonferroni")
summary(ans)

## Nemenyi's all-pairs comparison test
ans <- kwAllPairsNemenyiTest(count ~ spray, data = InsectSprays)
summary(ans)

## Brown-Mood all-pairs median test
ans <- medianAllPairsTest(count ~ spray, data = InsectSprays)
summary(ans)

Dunn's All-Pairs Rank Comparison Test

Description

Performs Dunn's non-parametric all-pairs comparison test for Kruskal-type ranked data.

Usage

kwAllPairsDunnTest(x, ...)

## Default S3 method:
kwAllPairsDunnTest(x, g, p.adjust.method = p.adjust.methods, ...)

## S3 method for class 'formula'
kwAllPairsDunnTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Dunn's non-parametric test can be performed. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: \mu_i(x) = \mu_j(x) is tested in the two-tailed test against the alternative A_{ij}: \mu_i(x) \ne \mu_j(x), ~~ i \ne j.

The p-values are computed from the standard normal distribution using any of the p-adjustment methods as included in p.adjust. Originally, Dunn (1964) proposed Bonferroni's p-adjustment method.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Dunn, O. J. (1964) Multiple comparisons using rank sums, Technometrics 6, 241–252.

Siegel, S., Castellan Jr., N. J. (1988) Nonparametric Statistics for The Behavioral Sciences. New York: McGraw-Hill.

See Also

Normal, p.adjust, kruskalTest, kwAllPairsConoverTest, kwAllPairsNemenyiTest

Examples

## Data set InsectSprays
## Global test
kruskalTest(count ~ spray, data = InsectSprays)

## Conover's all-pairs comparison test
## single-step means Tukey's p-adjustment
ans <- kwAllPairsConoverTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "single-step")
summary(ans)

## Dunn's all-pairs comparison test
ans <- kwAllPairsDunnTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "bonferroni")
summary(ans)

## Nemenyi's all-pairs comparison test
ans <- kwAllPairsNemenyiTest(count ~ spray, data = InsectSprays)
summary(ans)

## Brown-Mood all-pairs median test
ans <- medianAllPairsTest(count ~ spray, data = InsectSprays)
summary(ans)

Nemenyi's All-Pairs Rank Comparison Test

Description

Performs Nemenyi's non-parametric all-pairs comparison test for Kruskal-type ranked data.

Usage

kwAllPairsNemenyiTest(x, ...)

## Default S3 method:
kwAllPairsNemenyiTest(x, g, dist = c("Tukey", "Chisquare"), ...)

## S3 method for class 'formula'
kwAllPairsNemenyiTest(
  formula,
  data,
  subset,
  na.action,
  dist = c("Tukey", "Chisquare"),
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Nemenyi's non-parametric test can be performed. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: \theta_i(x) = \theta_j(x) is tested in the two-tailed test against the alternative A_{ij}: \theta_i(x) \ne \theta_j(x), ~~ i \ne j.

Let R_{ij} be the rank of X_{ij}, where X_{ij} is jointly ranked from \left\{1, 2, \ldots, N \right\}, ~~ N = \sum_{i=1}^k n_i, then the test statistic under the absence of ties is calculated as

t_{ij} = \frac{\bar{R}_j - \bar{R}_i} {\sigma_R \left(1/n_i + 1/n_j\right)^{1/2}} \qquad \left(i \ne j\right),

with \bar{R}_j, \bar{R}_i the mean rank of the i-th and j-th group and the expected variance as

\sigma_R^2 = N \left(N + 1\right) / 12.

A pairwise difference is significant, if |t_{ij}|/\sqrt{2} > q_{kv}, with k the number of groups and v = \infty the degree of freedom.

Sachs(1997) has given a modified approach for Nemenyi's test in the presence of ties for N > 6, k > 4 provided that the kruskalTest indicates significance: In the presence of ties, the test statistic is corrected according to \hat{t}_{ij} = t_{ij} / C, with

C = 1 - \frac{\sum_{i=1}^r t_i^3 - t_i}{N^3 - N}.

The function provides two different dist for p-value estimation:

Tukey

The p-values are computed from the studentized range distribution (alias Tukey), \mathrm{Pr} \left\{ t_{ij} \sqrt{2} \ge q_{k\infty\alpha} | mathrm{H} \right\} = \alpha.

Chisquare

The p-values are computed from the Chisquare distribution with v = k - 1 degree of freedom.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Nemenyi, P. (1963) Distribution-free Multiple Comparisons. Ph.D. thesis, Princeton University.

Sachs, L. (1997) Angewandte Statistik. Berlin: Springer.

Wilcoxon, F., Wilcox, R. A. (1964) Some rapid approximate statistical procedures. Pearl River: Lederle Laboratories.

See Also

Tukey, Chisquare, p.adjust, kruskalTest, kwAllPairsDunnTest, kwAllPairsConoverTest

Examples

## Data set InsectSprays
## Global test
kruskalTest(count ~ spray, data = InsectSprays)

## Conover's all-pairs comparison test
## single-step means Tukey's p-adjustment
ans <- kwAllPairsConoverTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "single-step")
summary(ans)

## Dunn's all-pairs comparison test
ans <- kwAllPairsDunnTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "bonferroni")
summary(ans)

## Nemenyi's all-pairs comparison test
ans <- kwAllPairsNemenyiTest(count ~ spray, data = InsectSprays)
summary(ans)

## Brown-Mood all-pairs median test
ans <- medianAllPairsTest(count ~ spray, data = InsectSprays)
summary(ans)

Conover's Many-to-One Rank Comparison Test

Description

Performs Conover's non-parametric many-to-one comparison test for Kruskal-type ranked data.

Usage

kwManyOneConoverTest(x, ...)

## Default S3 method:
kwManyOneConoverTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
kwManyOneConoverTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in an one-factorial layout with non-normally distributed residuals Conover's non-parametric test can be performed. Let there be k groups including the control, then the number of treatment levels is m = k - 1. Then m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: \theta_0 = \theta_i is tested in the two-tailed case against A_i: \theta_0 \ne \theta_i, ~~ (1 \le i \le m).

If p.adjust.method == "single-step" is selected, the p-values will be computed from the multivariate t distribution. Otherwise, the p-values are computed from the t-distribution using any of the p-adjustment methods as included in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Conover, W. J, Iman, R. L. (1979) On multiple-comparisons procedures, Tech. Rep. LA-7677-MS, Los Alamos Scientific Laboratory.

See Also

pmvt, TDist, kruskalTest, kwManyOneDunnTest, kwManyOneNdwTest

Examples

## Data set PlantGrowth
## Global test
kruskalTest(weight ~ group, data = PlantGrowth)

## Conover's many-one comparison test
## single-step means p-value from multivariate t distribution
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "single-step")
summary(ans)

## Conover's many-one comparison test
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Dunn's many-one comparison test
ans <- kwManyOneDunnTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Nemenyi's many-one comparison test
ans <- kwManyOneNdwTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Many one U test
ans <- manyOneUTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Chen Test
ans <- chenTest(weight ~ group, data = PlantGrowth,
                    p.adjust.method = "holm")
summary(ans)

Dunn's Many-to-One Rank Comparison Test

Description

Performs Dunn's non-parametric many-to-one comparison test for Kruskal-type ranked data.

Usage

kwManyOneDunnTest(x, ...)

## Default S3 method:
kwManyOneDunnTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
kwManyOneDunnTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in an one-factorial layout with non-normally distributed residuals Dunn's non-parametric test can be performed. Let there be k groups including the control, then the number of treatment levels is m = k - 1. Then m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: \theta_0 = \theta_i is tested in the two-tailed case against A_i: \theta_0 \ne \theta_i, ~~ (1 \le i \le m).

If p.adjust.method == "single-step" is selected, the p-values will be computed from the multivariate normal distribution. Otherwise, the p-values are computed from the standard normal distribution using any of the p-adjustment methods as included in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Dunn, O. J. (1964) Multiple comparisons using rank sums, Technometrics 6, 241–252.

Siegel, S., Castellan Jr., N. J. (1988) Nonparametric Statistics for The Behavioral Sciences. New York: McGraw-Hill.

See Also

pmvnorm, TDist, kruskalTest, kwManyOneConoverTest, kwManyOneNdwTest

Examples

## Data set PlantGrowth
## Global test
kruskalTest(weight ~ group, data = PlantGrowth)

## Conover's many-one comparison test
## single-step means p-value from multivariate t distribution
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "single-step")
summary(ans)

## Conover's many-one comparison test
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Dunn's many-one comparison test
ans <- kwManyOneDunnTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Nemenyi's many-one comparison test
ans <- kwManyOneNdwTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Many one U test
ans <- manyOneUTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Chen Test
ans <- chenTest(weight ~ group, data = PlantGrowth,
                    p.adjust.method = "holm")
summary(ans)

Nemenyi-Damico-Wolfe Many-to-One Rank Comparison Test

Description

Performs Nemenyi-Damico-Wolfe non-parametric many-to-one comparison test for Kruskal-type ranked data.

Usage

kwManyOneNdwTest(x, ...)

## Default S3 method:
kwManyOneNdwTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
kwManyOneNdwTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

For many-to-one comparisons (pairwise comparisons with one control) in an one-factorial layout with non-normally distributed residuals the Nemenyi-Damico-Wolfe non-parametric test can be performed. Let there be k groups including the control, then the number of treatment levels is m = k - 1. Then m pairwise comparisons can be performed between the i-th treatment level and the control. H_i: \theta_0 = \theta_i is tested in the two-tailed case against A_i: \theta_0 \ne \theta_i, ~~ (1 \le i \le m).

If p.adjust.method == "single-step" is selected, the p-values will be computed from the multivariate normal distribution. Otherwise, the p-values are computed from the standard normal distribution using any of the p-adjustment methods as included in p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

This function is essentially the same as kwManyOneDunnTest, but there is no tie correction included. Therefore, the implementation of Dunn's test is superior, when ties are present.

References

Damico, J. A., Wolfe, D. A. (1989) Extended tables of the exact distribution of a rank statistic for treatments versus control multiple comparisons in one-way layout designs, Communications in Statistics - Theory and Methods 18, 3327–3353.

Nemenyi, P. (1963) Distribution-free Multiple Comparisons, Ph.D. thesis, Princeton University.

See Also

pmvt, TDist, kruskalTest, kwManyOneDunnTest, kwManyOneConoverTest

Examples

## Data set PlantGrowth
## Global test
kruskalTest(weight ~ group, data = PlantGrowth)

## Conover's many-one comparison test
## single-step means p-value from multivariate t distribution
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "single-step")
summary(ans)

## Conover's many-one comparison test
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Dunn's many-one comparison test
ans <- kwManyOneDunnTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Nemenyi's many-one comparison test
ans <- kwManyOneNdwTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Many one U test
ans <- manyOneUTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Chen Test
ans <- chenTest(weight ~ group, data = PlantGrowth,
                    p.adjust.method = "holm")
summary(ans)

Testing against Ordered Alternatives (Le's Test)

Description

Performs Le's test for testing against ordered alternatives.

Usage

leTest(x, ...)

## Default S3 method:
leTest(x, g, alternative = c("two.sided", "greater", "less"), ...)

## S3 method for class 'formula'
leTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  ...
)

Arguments

Details

The null hypothesis, H_0: \theta_1 = \theta_2 = \ldots = \theta_k is tested against a simple order hypothesis, H_\mathrm{A}: \theta_1 \le \theta_2 \le \ldots \le \theta_k,~\theta_1 < \theta_k.

The p-values are estimated from the standard normal distribution.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

Le, C. T. (1988) A new rank test against ordered alternatives in k-sample problems, Biometrical Journal 30, 87–92.

See Also

kruskalTest and shirleyWilliamsTest of the package PMCMRplus, kruskal.test of the library stats.

Examples

## Example from Sachs (1997, p. 402)
x <- c(106, 114, 116, 127, 145,
       110, 125, 143, 148, 151,
       136, 139, 149, 160, 174)
g <- gl(3,5)
levels(g) <- c("A", "B", "C")

## Chacko's test
chackoTest(x, g)

## Cuzick's test
cuzickTest(x, g)

## Johnson-Mehrotra test
johnsonTest(x, g)

## Jonckheere-Terpstra test
jonckheereTest(x, g)

## Le's test
leTest(x, g)

## Spearman type test
spearmanTest(x, g)

## Murakami's BWS trend test
bwsTrendTest(x, g)

## Fligner-Wolfe test
flignerWolfeTest(x, g)

## Shan-Young-Kang test
shanTest(x, g)

Least Significant Difference Test

Description

Performs the least significant difference all-pairs comparisons test for normally distributed data with equal group variances.

Usage

lsdTest(x, ...)

## Default S3 method:
lsdTest(x, g, ...)

## S3 method for class 'formula'
lsdTest(formula, data, subset, na.action, ...)

## S3 method for class 'aov'
lsdTest(x, ...)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with normally distributed residuals and equal variances the least signifiant difference test can be performed after a significant ANOVA F-test. Let X_{ij} denote a continuous random variable with the j-the realization (1 \le j \le n_i) in the i-th group (1 \le i \le k). Furthermore, the total sample size is N = \sum_{i=1}^k n_i. A total of m = k(k-1)/2 hypotheses can be tested: The null hypothesis is H_{ij}: \mu_i = \mu_j ~~ (i \ne j) is tested against the alternative A_{ij}: \mu_i \ne \mu_j (two-tailed). Fisher's LSD all-pairs test statistics are given by

t_{ij} \frac{\bar{X}_i - \bar{X_j}} {s_{\mathrm{in}} \left(1/n_j + 1/n_i\right)^{1/2}}, ~~ (i \ne j)

with s^2_{\mathrm{in}} the within-group ANOVA variance. The null hypothesis is rejected if |t_{ij}| > t_{v\alpha/2}, with v = N - k degree of freedom. The p-values (two-tailed) are computed from the TDist distribution.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

As there is no p-value adjustment included, this function is equivalent to Fisher's protected LSD test, provided that the LSD test is only applied after a significant one-way ANOVA F-test. If one is interested in other types of LSD test (i.e. with p-value adustment) see function pairwise.t.test.

References

Sachs, L. (1997) Angewandte Statistik, New York: Springer.

See Also

TDist, pairwise.t.test

Examples

fit <- aov(weight ~ feed, chickwts)
shapiro.test(residuals(fit))
bartlett.test(weight ~ feed, chickwts)
anova(fit)

## also works with fitted objects of class aov
res <- lsdTest(fit)
summary(res)
summaryGroup(res)

Mack-Wolfe Test for Umbrella Alternatives

Description

Performs Mack-Wolfe non-parametric test for umbrella alternatives.

Usage

mackWolfeTest(x, ...)

## Default S3 method:
mackWolfeTest(x, g, p = NULL, nperm = 1000, ...)

## S3 method for class 'formula'
mackWolfeTest(formula, data, subset, na.action, p = NULL, nperm = 1000, ...)

Arguments

Details

In dose-finding studies one may assume an increasing treatment effect with increasing dose level. However, the test subject may actually succumb to toxic effects at high doses, which leads to decresing treatment effects.

The scope of the Mack-Wolfe Test is to test for umbrella alternatives for either a known or unknown point p (i.e. dose-level), where the peak (umbrella point) is present.

H_i: \theta_0 = \theta_i = \ldots = \theta_k is tested against the alternative A_i: \theta_1 \le \ldots \theta_p \ge \theta_k for some p, with at least one strict inequality.

If p = NULL (peak unknown), the upper-tail p-value is computed via an asymptotic bootstrap permutation test.

If an integer value for p is given (peak known), the upper-tail p-value is computed from the standard normal distribution (pnorm).

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

One may increase the number of permutations to e.g. nperm = 10000 in order to get more precise p-values. However, this will be on the expense of computational time.

References

Chen, I. Y. (1991) Notes on the Mack-Wolfe and Chen-Wolfe Tests for Umbrella Alternatives, Biom. J. 33, 281–290.

Mack, G. A., Wolfe, D. A. (1981) K-sample rank tests for umbrella alternatives, J. Amer. Statist. Assoc. 76, 175–181.

See Also

pnorm, sample.

Examples

## Example from Table 6.10 of Hollander and Wolfe (1999).
## Plates with Salmonella bacteria of strain TA98 were exposed to
## various doses of Acid Red 114 (in mu g / ml).
## The data are the numbers of visible revertant colonies on 12 plates.
## Assume a peak at D333 (i.e. p = 3).
x <- c(22, 23, 35, 60, 59, 54, 98, 78, 50, 60, 82, 59, 22, 44,
  33, 23, 21, 25)
g <- as.ordered(rep(c(0, 100, 333, 1000, 3333, 10000), each=3))
plot(x ~ g)
mackWolfeTest(x=x, g=g, p=3)

Mandel's h Test According to E 691 ASTM

Description

The function calculates the consistency statistics h and corresponding p-values for each group (lab) according to Practice E 691 ASTM.

Usage

mandelhTest(x, ...)

## Default S3 method:
mandelhTest(x, g, ...)

## S3 method for class 'formula'
mandelhTest(formula, data, subset, na.action, ...)

Arguments

Value

A list with class "mandel" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

p.value

the p-value for the test.

statistic

the estimated quantiles of Mandel's statistic.

alternative

a character string describing the alternative hypothesis.

grouplev

a character vector describing the levels of the groups.

nrofrepl

the number of replicates for each group.

References

Practice E 691 (2005) Standard Practice for Conducting an Interlaboratory Study to Determine the Precision of a Test Method, ASTM International.

See Also

qmandelh pmandelh

Examples

data(Pentosan)
mandelhTest(value ~ lab, data=Pentosan, subset=(material == "A"))

Mandel's k Test According to E 691 ASTM

Description

The function calculates the consistency statistics k and corresponding p-values for each group (lab) according to Practice E 691 ASTM.

Usage

mandelkTest(x, ...)

## Default S3 method:
mandelkTest(x, g, ...)

## S3 method for class 'formula'
mandelkTest(formula, data, subset, na.action, ...)

Arguments

Value

A list with class "mandel" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

p.value

the p-value for the test.

statistic

the estimated quantiles of Mandel's statistic.

alternative

a character string describing the alternative hypothesis.

grouplev

a character vector describing the levels of the groups.

nrofrepl

the number of replicates for each group.

References

Practice E 691 (2005) Standard Practice for Conducting an Interlaboratory Study to Determine the Precision of a Test Method, ASTM International.

See Also

qmandelk pmandelk

Examples

data(Pentosan)
mandelkTest(value ~ lab, data=Pentosan, subset=(material == "A"))

Multiple Comparisons with One Control (U-test)

Description

Performs pairwise comparisons of multiple group levels with one control.

Usage

manyOneUTest(x, ...)

## Default S3 method:
manyOneUTest(
  x,
  g,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
manyOneUTest(
  formula,
  data,
  subset,
  na.action,
  alternative = c("two.sided", "greater", "less"),
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

This functions performs Wilcoxon, Mann and Whitney's U-test for a one factorial design where each factor level is tested against one control (m = k -1 tests). As the data are re-ranked for each comparison, this test is only suitable for balanced (or almost balanced) experimental designs.

For the two-tailed test and p.adjust.method = "single-step" the multivariate normal distribution is used for controlling Type 1 error and to calculate p-values. Otherwise, the p-values are calculated from the standard normal distribution with any latter p-adjustment as available by p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

Note

Factor labels for g must be assigned in such a way, that they can be increasingly ordered from zero-dose control to the highest dose level, e.g. integers {0, 1, 2, ..., k} or letters {a, b, c, ...}. Otherwise the function may not select the correct values for intended zero-dose control.

It is safer, to i) label the factor levels as given above, and to ii) sort the data according to increasing dose-levels prior to call the function (see order, factor).

References

OECD (ed. 2006) Current approaches in the statistical analysis of ecotoxicity data: A guidance to application, OECD Series on testing and assessment, No. 54.

See Also

wilcox.test, pmvnorm, Normal

Examples

## Data set PlantGrowth
## Global test
kruskalTest(weight ~ group, data = PlantGrowth)

## Conover's many-one comparison test
## single-step means p-value from multivariate t distribution
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "single-step")
summary(ans)

## Conover's many-one comparison test
ans <- kwManyOneConoverTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Dunn's many-one comparison test
ans <- kwManyOneDunnTest(weight ~ group, data = PlantGrowth,
                             p.adjust.method = "holm")
summary(ans)

## Nemenyi's many-one comparison test
ans <- kwManyOneNdwTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Many one U test
ans <- manyOneUTest(weight ~ group, data = PlantGrowth,
                        p.adjust.method = "holm")
summary(ans)

## Chen Test
ans <- chenTest(weight ~ group, data = PlantGrowth,
                    p.adjust.method = "holm")
summary(ans)

Brown-Mood All Pairs Median Test

Description

Performs Brown-Mood All Pairs Median Test.

Usage

medianAllPairsTest(x, ...)

## Default S3 method:
medianAllPairsTest(x, g, p.adjust.method = p.adjust.methods, ...)

## S3 method for class 'formula'
medianAllPairsTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = p.adjust.methods,
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Brown-Mood non-parametric Median test can be performed. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: \mu_i(x) = \mu_j(x) is tested in the two-tailed test against the alternative A_{ij}: \mu_i(x) \ne \mu_j(x), ~~ i \ne j.

In this procedure the joined median is used for classification, but pairwise Pearson Chisquare-Tests are conducted. Any method as given by p.adjust.methods can be used to account for multiplicity.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Brown, G.W., Mood, A.M., 1951, On Median Tests for Linear Hypotheses, in: Proceedings of the Second Berkeley Symposium on Mathematical Statistics and Probability. University of California Press, pp. 159–167.

See Also

chisq.test.

Examples

## Data set InsectSprays
## Global test
kruskalTest(count ~ spray, data = InsectSprays)

## Conover's all-pairs comparison test
## single-step means Tukey's p-adjustment
ans <- kwAllPairsConoverTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "single-step")
summary(ans)

## Dunn's all-pairs comparison test
ans <- kwAllPairsDunnTest(count ~ spray, data = InsectSprays,
                             p.adjust.method = "bonferroni")
summary(ans)

## Nemenyi's all-pairs comparison test
ans <- kwAllPairsNemenyiTest(count ~ spray, data = InsectSprays)
summary(ans)

## Brown-Mood all-pairs median test
ans <- medianAllPairsTest(count ~ spray, data = InsectSprays)
summary(ans)

Brown-Mood Median Test

Description

Performs Brown-Mood Median Test.

Usage

medianTest(x, ...)

## Default S3 method:
medianTest(x, g, simulate.p.value = FALSE, B = 2000, ...)

## S3 method for class 'formula'
medianTest(
  formula,
  data,
  subset,
  na.action,
  simulate.p.value = FALSE,
  B = 2000,
  ...
)

Arguments

Details

The null hypothesis, H_0: \theta_1 = \theta_2 = \ldots = \theta_k is tested against the alternative, H_\mathrm{A}: \theta_i \ne \theta_j ~~(i \ne j), with at least one unequality beeing strict.

Value

A list with class ‘htest’. For details see chisq.test.

References

Brown, G.W., Mood, A.M., 1951, On Median Tests for Linear Hypotheses, in: Proceedings of the Second Berkeley Symposium on Mathematical Statistics and Probability. University of California Press, pp. 159–167.

See Also

chisq.test.

Examples

## Hollander & Wolfe (1973), 116.
## Mucociliary efficiency from the rate of removal of dust in normal
## subjects, subjects with obstructive airway disease, and subjects
## with asbestosis.
x <- c(2.9, 3.0, 2.5, 2.6, 3.2) # normal subjects
y <- c(3.8, 2.7, 4.0, 2.4)      # with obstructive airway disease
z <- c(2.8, 3.4, 3.7, 2.2, 2.0) # with asbestosis
g <- factor(x = c(rep(1, length(x)),
                   rep(2, length(y)),
                   rep(3, length(z))),
             labels = c("ns", "oad", "a"))
dat <- data.frame(
   g = g,
   x = c(x, y, z))

## AD-Test
adKSampleTest(x ~ g, data = dat)

## BWS-Test
bwsKSampleTest(x ~ g, data = dat)

## Kruskal-Test
## Using incomplete beta approximation
kruskalTest(x ~ g, dat, dist="KruskalWallis")
## Using chisquare distribution
kruskalTest(x ~ g, dat, dist="Chisquare")

## Not run: 
## Check with kruskal.test from R stats
kruskal.test(x ~ g, dat)

## End(Not run)
## Using Conover's F
kruskalTest(x ~ g, dat, dist="FDist")

## Not run: 
## Check with aov on ranks
anova(aov(rank(x) ~ g, dat))
## Check with oneway.test
oneway.test(rank(x) ~ g, dat, var.equal = TRUE)

## End(Not run)

## Median Test asymptotic
medianTest(x ~ g, dat)

## Median Test with simulated p-values
set.seed(112)
medianTest(x ~ g, dat, simulate.p.value = TRUE)

Madhava Rao-Raghunath Test for Testing Treatment vs. Control

Description

The function has implemented the nonparametric test of Madhava Rao and Raghunath (2016) for testing paired two-samples for symmetry. The null hypothesis H: F(x,y) = F(y,x) is tested against the alternative A: F(x,y) \ne F(y,x).

Usage

mrrTest(x, ...)

## Default S3 method:
mrrTest(x, y = NULL, m = NULL, ...)

## S3 method for class 'formula'
mrrTest(formula, data, subset, na.action, ...)

Arguments

Details

Let X_i and Y_i, ~ i \le n denote continuous variables that were observed on the same ith test item (e.g. patient) with i = 1, \ldots n. Let

U_i = X_i + Y_i \qquad V_i = X_i - Y_i

Let U_{(i)} be the ith order statistic, U_{(1)} \le U_{(2)} \le \ldots U_{(n)} and k the number of clusters, with the condition:

n = k ~ m.

Further, let the divider denote d_0 = -\infty, d_k = \infty, and else

d_j = \frac{ U_{(jm)} + U_{(jm+1)} }{2}, ~ 1 \le j \le k -1

The two counts are

n_j^{+} = \left\{ \begin{array}{lr} 1 & \mathrm{if}~ d_{j-1} < u_i < d_j, v_i > 0 \\ 0 & \end{array} \right.

and

n_j^{-} = \left\{ \begin{array}{lr} 1 & \mathrm{if}~ d_{j-1} < u_i < d_j, v_i \le 0 \\ 0 & \end{array} \right.

The test statistic is

M = \sum_{j = 1}^k \frac{\left(n_j^{+} - n_j^{-}\right)^2} {m}

The exact p-values for 5 \le n \le 30 are taken from an internal look-up table. The exact p-values were taken from Table 7, Appendix B of Madhava Rao and Raghunath (2016).

If m = NULL the function uses n = m for all prime numbers, otherwise it tries to find an value for m in such a way, that for k = n / m all variables are integer.

Value

A list with class "htest" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

the estimated quantile of the test statistic.

p.value

the p-value for the test.

parameter

the parameters of the test statistic, if any.

alternative

a character string describing the alternative hypothesis.

estimates

the estimates, if any.

null.value

the estimate under the null hypothesis, if any.

Note

The function returns an error code if a value for m is provided that does not lead to an integer of the ratio k = n /m.

The function also returns an error code, if a tabulated value for given n, m and calculated M can not be found in the look-up table.

References

Madhava Rao, K.S., Ragunath, M. (2016) A Simple Nonparametric Test for Testing Treatment Versus Control. J Stat Adv Theory Appl 16, 133–162. doi:10.18642/jsata_7100121717

Examples

## Madhava Rao and Raghunath (2016), p. 151
## Inulin clearance of living donors
## and recipients of their kidneys
x <- c(61.4, 63.3, 63.7, 80.0, 77.3, 84.0, 105.0)
y <- c(70.8, 89.2, 65.8, 67.1, 87.3, 85.1, 88.1)
mrrTest(x, y)

## formula method
## Student's Sleep Data
mrrTest(extra ~ group, data = sleep)

Lu-Smith All-Pairs Comparison Normal Scores Test

Description

Performs Lu-Smith all-pairs comparison normal scores test.

Usage

normalScoresAllPairsTest(x, ...)

## Default S3 method:
normalScoresAllPairsTest(
  x,
  g,
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

## S3 method for class 'formula'
normalScoresAllPairsTest(
  formula,
  data,
  subset,
  na.action,
  p.adjust.method = c("single-step", p.adjust.methods),
  ...
)

Arguments

Details

For all-pairs comparisons in an one-factorial layout with non-normally distributed residuals Lu and Smith's normal scores transformation can be used prior to an all-pairs comparison test. A total of m = k(k-1)/2 hypotheses can be tested. The null hypothesis H_{ij}: F_i(x) = F_j(x) is tested in the two-tailed test against the alternative A_{ij}: F_i(x) \ne F_j(x), ~~ i \ne j. For p.adjust.method = "single-step" the Tukey's studentized range distribution is used to calculate p-values (see Tukey). Otherwise, the t-distribution is used for the calculation of p-values with a latter p-value adjustment as performed by p.adjust.

Value

A list with class "PMCMR" containing the following components:

method

a character string indicating what type of test was performed.

data.name

a character string giving the name(s) of the data.

statistic

lower-triangle matrix of the estimated quantiles of the pairwise test statistics.

p.value

lower-triangle matrix of the p-values for the pairwise tests.

alternative

a character string describing the alternative hypothesis.

p.adjust.method

a character string describing the method for p-value adjustment.

model

a data frame of the input data.

dist

a string that denotes the test distribution.

References

Lu, H., Smith, P. (1979) Distribution of normal scores statistic for nonparametric one-way analysis of variance. Journal of the American Statistical Association 74, 715–722.

See Also

normalScoresTest, normalScoresManyOneTest, normOrder.