Test hipotez dla różnicy w medianach pomiędzy więcej niż dwiema próbkami

Pytanie

Wyniki testu trzech grup ludzi są zapisywane jako osobne wektory w R.

set.seed(1)
group1 <- rnorm(100, mean = 75, sd = 10)
group2 <- rnorm(100, mean = 85, sd = 10)
group3 <- rnorm(100, mean = 95, sd = 10)

Chcę wiedzieć, czy istnieje znacząca różnica w medianach między tymi grupami. Wiem, że mogłem przetestować grupę 1 w porównaniu z grupą 2 za pomocą testu Wilcoxona.

wilcox.test(group1, group2)

Jednak porównuje to tylko dwie grupy na raz i chciałbym porównać wszystkie trzy jednocześnie. Chciałbym przeprowadzić test statystyczny, który daje wartość ap na poziomie istotności 0,05. Czy ktoś mógłby pomóc?

Edycja # 1 - Mediana testu Mooda

Po sugerowanej odpowiedzi użytkownika Hibernating spróbowałem testu mediany Mooda.

median.test <- function(x, y){
    z <- c(x, y)
    g <- rep(1:2, c(length(x), length(y)))
    m <- median(z)
    fisher.test(z < m, g)$p.value
}

median.test(group1, group2)

Jednak takie podejście pozwala mi przetestować istotną różnicę między medianami tylko dwóch grup jednocześnie. Nie jestem pewien, jak go użyć, aby porównać mediany wszystkich trzech jednocześnie.

Edycja # 2 - Test Kruskala-Wallisa

Sugerowana odpowiedź użytkownika dmartin wydaje się być mniej więcej tym, czego potrzebuję i pozwala mi przetestować wszystkie trzy grupy jednocześnie.

kruskal.test(list(group1, group2, group3))

Edytuj # 3

Użytkownik Greg Snow z pomocą zauważa w swojej odpowiedzi, że test Kruskala-Wallisa jest odpowiedni, o ile zawiera on surowe założenia, które czynią go również testem środków.

Można również zastosować test Kruskala-Wallisa, ponieważ jest to nieparametryczna ANOVA. Ponadto często uważa się go za silniejszy niż test mediany Mooda . Można go zaimplementować w języku R za pomocą funkcji kruskal.test w pakiecie statystyk w języku R.

Aby odpowiedzieć na Twoją edycję, interpretacja KW jest podobna do jednostronnej ANOVA. Znacząca wartość p odpowiada odrzuconemu zeru, że wszystkie trzy średnie są równe. Musisz użyć testu kontrolnego (ponownie, podobnie jak ANOVA), aby odpowiedzieć na pytania dotyczące określonych grup. Zwykle wynika to z określonych pytań badawczych. Wystarczy spojrzeć na parametry symulacji, wszystkie trzy grupy powinny się znacznie różnić od siebie, jeśli wykonasz test kontrolny (ponieważ wszystkie są 1 SD z wyjątkiem N = 100).

Po pierwsze, test Wilcoxona (lub test Manna-Whitneya) nie jest testem median (chyba że przyjmujesz bardzo surowe założenia, które również czynią z niego test środków). Dla porównania więcej niż 2 grup test Wilcoxona może prowadzić do paradoksalnych wyników (patrz Kości Efrona ).

Ponieważ test Wilcoxona jest tylko specjalnym przypadkiem testu permutacji i jesteś szczególnie zainteresowany medianami, proponuję test permutacji na medianach.

Najpierw wybierz miarę różnicy, coś w rodzaju największej z 3 median minus najmniejsza z 3 (lub wariancji 3 median lub MAD itp.).

Teraz oblicz swoje statystyki dla oryginalnych danych.

połącz wszystkie dane w jeden zestaw, a następnie losowo podziel wartości na 3 grupy

te same rozmiary co oryginał i oblicz te same statystyki.

powtarzaj wiele razy (jak 9998)

Porównaj porównanie statystyki z rzeczywistych danych z rozkładem wszystkich statystyk dla twojego testu.

Mediana testu Mooda to test nieparametryczny, który służy do testowania równości median z dwóch lub więcej populacji. Zobacz tutaj część R swojego pytania. Zobacz także powiązane pytanie tutaj . Również stąd :

Test mediany Mooda jest najłatwiejszy do wykonania ręcznie: oblicz ogólną medianę (wszystkich danych) i policz, ile wartości jest powyżej i poniżej mediany w każdej grupie. Jeśli grupy są prawie takie same, obserwacje powinny wynosić około 50-50 powyżej i poniżej ogólnej mediany w każdej grupie ... Liczby poniżej mediany i powyżej mediany ... tworzą tabelę dwukierunkową, która jest następnie analizowany za pomocą testu chi-kwadrat. Mediana testu Mooda jest bardzo podobna do testu znaku uogólnionego na dwie lub więcej grup.

Edycja: w przypadku trzech grup możesz rozważyć tę prostą generalizację kodu R, do którego linkowałem:

median.test2 <- function(x, y, z) {
  a <- c(x, y, z)
  g <- rep(1:3, c(length(x), length(y), length(z)))
  m <- median(a)
  fisher.test(a < m, g)$p.value
}

Wiem, że jest już późno, ale nie mogłem też znaleźć dobrego pakietu dla testu mediany Mooda, więc wziąłem na siebie zadanie stworzenia funkcji w R, która wydaje się załatwić sprawę.

#Mood's median test for a data frame with one column containing data (d),
#and another containing a factor/grouping variable (f)

moods.median = function(d,f) {

    #make a new matrix data frame
    m = cbind(f,d)
    colnames(m) = c("group", "value")


    #get the names of the factors/groups
    facs = unique(f)

    #count the number of factors/groups
    factorN = length(unique(f))


    #Make a 2 by K table that will be saved to the global environment by using "<<-":
    #2 rows (number of values > overall median & number of values <= overall median)
    #K-many columns for each level of the factor
    MoodsMedianTable <<- matrix(NA, nrow = 2, ncol = factorN)

    rownames(MoodsMedianTable) <<- c("> overall median", "<= overall median")
    colnames(MoodsMedianTable) <<- c(facs[1:factorN])
    colnames(MoodsMedianTable) <<- paste("Factor: ",colnames(MoodsMedianTable))


    #get the overall median
    overallmedian = median(d)



    #put the following into the 2 by K table:
    for(j in 1:factorN){ #for each factor level

        g = facs[j] #assign a temporary "group name"


        #count the number of observations in the factor that are greater than
        #the overall median and save it to the table
        MoodsMedianTable[1,j] <<- sum(m[,2][ which(m[,1]==g)] > overallmedian)


        #count the number of observations in the factor that are less than
        # or equal to the overall median and save it to the table
        MoodsMedianTable[2,j] <<- sum(m[,2][ which(m[,1]==g)] <= overallmedian)

    }


    #percent of cells with expected values less than 5
    percLT5 = ((sum(chisq.test(MoodsMedianTable)$expected < 5)) /
        (length(chisq.test(MoodsMedianTable)$expected)))


    #if >20% of cells have expected values less than 5
    #then give chi-squared stat, df, and Fisher's exact p.value
    if (percLT5 > 0.2) {
        return(list(
            "Chi-squared" = chisq.test(MoodsMedianTable)$statistic,
            "df" = chisq.test(MoodsMedianTable)$parameter,
            "Fisher's exact p.value" = fisher.test(MoodsMedianTable)$p.value))

    }


    #if <= 20% of cells have expected values less than 5
    #then give chi-squared stat, df, and chi-squared p.value
    if (percLT5 <= 0.2) {
        return(list(
            "Chi-squared" = chisq.test(MoodsMedianTable)$statistic,
            "df" = chisq.test(MoodsMedianTable)$parameter,
            "Chi-squared p.value" = chisq.test(MoodsMedianTable)$p.value))

    }

}

W przypadku pytania PO najpierw uruchom to, aby utworzyć nową ramkę danych, w której będą przechowywane wartości z trzech wektorów grupy za pomocą dopasowanej zmiennej „group”.

require(reshape2)
df = cbind(group1, group2, group3)
df = melt(df)
colnames(df) = c("observation", "group", "value")

i uruchom funkcję testu mediany Mooda za pomocą moods.median(df$value, df$group)