1#  File src/library/stats/R/friedman.test.R
2#  Part of the R package, https://www.R-project.org
3#
4#  Copyright (C) 1995-2015 The R Core Team
5#
6#  This program is free software; you can redistribute it and/or modify
7#  it under the terms of the GNU General Public License as published by
8#  the Free Software Foundation; either version 2 of the License, or
9#  (at your option) any later version.
10#
11#  This program is distributed in the hope that it will be useful,
12#  but WITHOUT ANY WARRANTY; without even the implied warranty of
13#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14#  GNU General Public License for more details.
15#
16#  A copy of the GNU General Public License is available at
17#  https://www.R-project.org/Licenses/
18
19friedman.test <- function(y, ...) UseMethod("friedman.test")
20
21friedman.test.default <-
22function(y, groups, blocks, ...)
23{
24    DNAME <- deparse1(substitute(y))
25    if (is.matrix(y)) {
26        groups <- factor(c(col(y)))
27        blocks <- factor(c(row(y)))
28    }
29    else {
30        if (anyNA(groups) || anyNA(blocks))
31            stop("NA's are not allowed in 'groups' or 'blocks'")
32        if (any(diff(c(length(y), length(groups), length(blocks))) != 0L))
33            stop("'y', 'groups' and 'blocks' must have the same length")
34        DNAME <- paste0(DNAME, ", ", deparse1(substitute(groups)),
35                        " and ", deparse1(substitute(blocks)))
36        if (any(table(groups, blocks) != 1))
37            stop("not an unreplicated complete block design")
38        groups <- factor(groups)
39        blocks <- factor(blocks)
40        ## Need to ensure consistent order of observations within
41        ## blocks.
42        o <- order(groups, blocks)
43        y <- y[o]
44        groups <- groups[o]
45        blocks <- blocks[o]
46    }
47
48    k <- nlevels(groups)
49    ## <FIXME split.matrix>
50    y <- matrix(unlist(split(c(y), blocks)), ncol = k, byrow = TRUE)
51    y <- y[complete.cases(y), ]
52    n <- nrow(y)
53    r <- t(apply(y, 1L, rank))
54    ## <FIXME split.matrix>
55    TIES <- tapply(c(r), row(r), table)
56    STATISTIC <- 12 * sum((colSums(r) - n * (k + 1) / 2)^2) /
57        (n * k * (k + 1) - sum(unlist(lapply(TIES, function(u) u^3 - u))) / (k-1))
58    PARAMETER <- k - 1
59    PVAL <- pchisq(STATISTIC, PARAMETER, lower.tail = FALSE)
60    names(STATISTIC) <- "Friedman chi-squared"
61    names(PARAMETER) <- "df"
62
63    structure(list(statistic = STATISTIC,
64                   parameter = PARAMETER,
65                   p.value = PVAL,
66                   method = "Friedman rank sum test",
67                   data.name = DNAME),
68              class = "htest")
69}
70
71friedman.test.formula <-
72function(formula, data, subset, na.action, ...)
73{
74    if(missing(formula))
75        stop("formula missing")
76    ## <FIXME>
77    ## Maybe put this into an internal rewriteTwoWayFormula() when
78    ## adding support for strata()
79    if((length(formula) != 3L)
80       || (length(formula[[3L]]) != 3L)
81       || (formula[[3L]][[1L]] != as.name("|"))
82       || (length(formula[[3L]][[2L]]) != 1L)
83       || (length(formula[[3L]][[3L]]) != 1L))
84        stop("incorrect specification for 'formula'")
85    formula[[3L]][[1L]] <- as.name("+")
86    ## </FIXME>
87    m <- match.call(expand.dots = FALSE)
88    m$formula <- formula
89    if(is.matrix(eval(m$data, parent.frame())))
90        m$data <- as.data.frame(data)
91    ## need stats:: for non-standard evaluation
92    m[[1L]] <- quote(stats::model.frame)
93    mf <- eval(m, parent.frame())
94    DNAME <- paste(names(mf), collapse = " and ")
95    names(mf) <- NULL
96    y <- do.call("friedman.test", as.list(mf))
97    y$data.name <- DNAME
98    y
99}
100