Hervorheben signifikanter Ergebnisse aus nicht parametrischen Mehrfachvergleichen auf Boxplots

9

Ich habe Box-Plots von 13 Gruppen, die ich in einem Plot zeige. Die Gruppen haben unausgeglichene Bevölkerungsgruppen und sind nicht normal verteilt. Ich möchte zeigen, welche Paare statistisch ähnlich sind (dh kruskal.test p-Wert <0,05 haben), indem ich a, b, c usw. auf übereinstimmende Kästchen setze. Hier ist ein Pseudocode, um zu zeigen, was ich habe:

A = c(1, 5, 8, 17, 16, 3, 24, 19, 6) 
B = c(2, 16, 5, 7, 4, 7, 3) 
C = c(1, 1, 3, 7, 9, 6, 10, 13) 
D = c(2, 15, 2, 9, 7) 
junk = list(g1=A, g2=B, g3=C, g4=D) 
boxplot(junk) 

Hier ist eine Handlung, die ich gefunden habe und die macht, was ich will (außer ich habe 13 Gruppen in einer Reihe):

Ilik
quelle

Antworten:

7

Der einfachste Code, der mir in den Sinn kommt, ist unten dargestellt. Ich bin mir ziemlich sicher, dass es auf CRAN bereits einige Funktionen gibt, aber ich bin zu faul, um danach zu suchen, selbst bei R-seek .

dd <- data.frame(y=as.vector(unlist(junk)), 
                 g=rep(paste("g", 1:4, sep=""), unlist(lapply(junk, length))))

aov.res <- kruskal.test(y ~ g, data=dd)
alpha.level <- .05/nlevels(dd$g)  # Bonferroni correction, but use 
                                  # whatever you want using p.adjust()

# generate all pairwise comparisons
idx <- combn(nlevels(dd$g), 2)

# compute p-values from Wilcoxon test for all comparisons
pval.res <- numeric(ncol(idx))
for (i in 1:ncol(idx))
  # test all group, pairwise
  pval.res[i] <- with(dd, wilcox.test(y[as.numeric(g)==idx[1,i]], 
                                      y[as.numeric(g)==idx[2,i]]))$p.value

# which groups are significantly different (arranged by column)
signif.pairs <- idx[,which(pval.res<alpha.level)]

boxplot(y ~ g, data=dd, ylim=c(min(dd$y)-1, max(dd$y)+1))
# use offset= to increment space between labels, thanks to vectorization
for (i in 1:ncol(signif.pairs))
    text(signif.pairs[,i], max(dd$y)+1, letters[i], pos=4, offset=i*.8-1)

Hier ist ein Beispiel dafür, was der obige Code erzeugen würde (mit signifikanten Unterschieden zwischen den vier Gruppen):

Geben Sie hier die Bildbeschreibung ein

Anstelle des Wilcoxon-Tests könnte man sich auf die in der kruskalmc()Funktion implementierte Prozedur aus dem pgirmess- Paket verlassen (siehe eine Beschreibung der hier verwendeten Prozedur ).

Lesen Sie auch unbedingt die R-Tipps von Rudolf Cardinal zu R : Grundlegende Diagramme 2 (siehe insbesondere Ein weiteres Balkendiagramm mit Anmerkungen ).

chl
quelle
Danke, chl. Ich bin kein Statistiker, daher muss ich Ihre Antwort genauer lernen, aber es sieht so aus, als würde es das tun, was ich brauche.
Ilik
Wenn etwas unklar ist, zögern Sie nicht zu fragen. Hier gibt es nicht viel statistisches Material: Ich habe Wilcoxon-Tests für paarweise Gruppenvergleiche verwendet und einzelne p-Werte korrigiert, um das Gesamtrisiko für falsch positive Ergebnisse auf 5% zu begrenzen. Das npmc- Paket enthält zusätzliche Funktionen für die Verarbeitung nicht parametrischer Mehrfachvergleiche, aber es gibt auch das Münz- Framework. Der Rest ist reiner R-Code mit Basis-R-Grafiken. es könnte auch mit gitter oder ggplot2 gemacht werden.
Chl
Ich entschuldige mich für meine Unkenntnis in der Statistik. Also habe ich Ihren Code ausprobiert und zuerst bemerkt, dass Sie den kruskal.test berechnen, aber das Ergebnis nicht verwenden (aov.res). Ich verstehe jetzt, dass der wilcox.test ein Sonderfall für kruskal für zwei Proben ist. Aber dann habe ich versucht, die Werte in meinen Gruppen zu ändern, um sie (intuitiv) anders zu machen und zu sehen, was dabei herauskommt. A = c (10, 50, 18, 17, 16, 31, 24, 19, 6) B = c (10, 50, 18, 17, 16, 30, 25, 18, 7) C = c (1, 1, 3, 7, 9, 6, 10, 13) D = c (200, 158, 206, 119, 77). g1 & g2 sind jetzt deutlich unterschiedlich ?? und g3 & g4 nicht? Ich bin mir nicht sicher, warum.
Ilik
@Ilik In meinem Code ist etwas schief gelaufen (die with(dd,Anweisung war an der falschen Stelle, was zu einem seltsamen Test führte!). Danke, dass du das verstanden hast! Ja, ich verwende keine Ergebnisse aus dem KW-Test, aber es ist immer eine gute Idee, diese zuerst zu überprüfen, da sonst mehrere Post-hoc-Tests bedeutungslos wären (oder zumindest auf Daten-Snooping hindeuten). Hinweis: Ich habe den Code korrigiert und nichts hat sich als signifikant herausgestellt, aber ich habe das Originalbild aus Gründen der Klarheit mit signifikanten Ergebnissen belassen.
Chl
Für alle Interessierten habe ich den Code ein wenig geändert, um die Ergebnisse in das Boxplot zu schreiben: MyText = rep ('', nlevels (dd g)), rep (max (dd g)), MyText)g))for(iin1:ncol(signif.pairs))MyText[signif.pairs[,i]]=paste(MyText[signif.pairs[,i]],letters[i],sep=)text(c(1:nlevels(ddy)+1,nlevels(dd
Ilik