Es ist üblich, Sterne auf Balkendiagrammen oder Boxplots zu setzen, um die Signifikanzstufe (p-Wert) einer oder zwischen zwei Gruppen anzuzeigen. Nachfolgend einige Beispiele:
Die Anzahl der Sterne wird durch den p-Wert definiert. Beispielsweise kann man 3 Sterne für den p-Wert <0,001, zwei Sterne für den p-Wert <0,01 usw. setzen (dies ändert sich jedoch von einem Artikel zum anderen).
Und meine Fragen: Wie generiere ich ähnliche Diagramme? Die Methoden, mit denen Sterne automatisch anhand des Signifikanzniveaus gesetzt werden, sind mehr als willkommen.
Bitte finden Sie meinen Versuch unten.
Zuerst erstellte ich einige Dummy-Daten und ein Balkendiagramm, das nach Wunsch geändert werden kann.
windows(4,4)
dat <- data.frame(Group = c("S1", "S1", "S2", "S2"),
Sub = c("A", "B", "A", "B"),
Value = c(3,5,7,8))
## Define base plot
p <-
ggplot(dat, aes(Group, Value)) +
theme_bw() + theme(panel.grid = element_blank()) +
coord_cartesian(ylim = c(0, 15)) +
scale_fill_manual(values = c("grey80", "grey20")) +
geom_bar(aes(fill = Sub), stat="identity", position="dodge", width=.5)
Das Hinzufügen von Sternchen über einer Spalte ist einfach, wie bereits erwähnt. Erstellen Sie einfach einen data.frame
mit den Koordinaten.
label.df <- data.frame(Group = c("S1", "S2"),
Value = c(6, 9))
p + geom_text(data = label.df, label = "***")
Um die Bögen hinzuzufügen, die einen Untergruppenvergleich anzeigen, berechnete ich die parametrischen Koordinaten eines Halbkreises und fügte sie hinzu, die mit geom_line
verbunden waren. Sternchen brauchen auch neue Koordinaten.
label.df <- data.frame(Group = c(1,1,1, 2,2,2),
Value = c(6.5,6.8,7.1, 9.5,9.8,10.1))
# Define arc coordinates
r <- 0.15
t <- seq(0, 180, by = 1) * pi / 180
x <- r * cos(t)
y <- r*5 * sin(t)
arc.df <- data.frame(Group = x, Value = y)
p2 <-
p + geom_text(data = label.df, label = "*") +
geom_line(data = arc.df, aes(Group+1, Value+5.5), lty = 2) +
geom_line(data = arc.df, aes(Group+2, Value+8.5), lty = 2)
Um den Vergleich zwischen den Gruppen anzuzeigen, habe ich einen größeren Kreis gebildet und ihn oben abgeflacht.
r <- .5
x <- r * cos(t)
y <- r*4 * sin(t)
y[20:162] <- y[20] # Flattens the arc
arc.df <- data.frame(Group = x, Value = y)
p2 + geom_line(data = arc.df, aes(Group+1.5, Value+11), lty = 2) +
geom_text(x = 1.5, y = 12, label = "***")
Ich weiß, dass dies eine alte Frage ist und die Antwort von Jens Tierling bietet bereits eine Lösung für das Problem. Ich habe aber kürzlich eine ggplot-Erweiterung erstellt, die das Hinzufügen von Signifikanzleisten vereinfacht: ggsignif
Anstatt langwierig geom_line
und geom_text
zu Ihrem Plot hinzuzufügen, fügen Sie einfach eine einzelne Ebene hinzu geom_signif
:
library(ggplot2)
library(ggsignif)
ggplot(iris, aes(x=Species, y=Sepal.Length)) +
geom_boxplot() +
geom_signif(comparisons = list(c("versicolor", "virginica")),
map_signif_level=TRUE)
Um ein besseres Diagramm zu erstellen, das dem von Jens Tierling gezeigten ähnlich ist, können Sie Folgendes tun:
dat <- data.frame(Group = c("S1", "S1", "S2", "S2"),
Sub = c("A", "B", "A", "B"),
Value = c(3,5,7,8))
ggplot(dat, aes(Group, Value)) +
geom_bar(aes(fill = Sub), stat="identity", position="dodge", width=.5) +
geom_signif(stat="identity",
data=data.frame(x=c(0.875, 1.875), xend=c(1.125, 2.125),
y=c(5.8, 8.5), annotation=c("**", "NS")),
aes(x=x,xend=xend, y=y, yend=y, annotation=annotation)) +
geom_signif(comparisons=list(c("S1", "S2")), annotations="***",
y_position = 9.3, tip_length = 0, vjust=0.4) +
scale_fill_manual(values = c("grey80", "grey20"))
Die vollständige Dokumentation des Pakets finden Sie unter CRAN .
Es gibt auch eine Erweiterung des Pakets ggsignif mit dem Namen ggpubr , die bei Vergleichen von mehreren Gruppen leistungsfähiger ist. Es baut auf ggsignif auf, behandelt aber auch Anova und Kruskal-wallis sowie paarweise Vergleiche mit dem globalen Mittelwert.
Beispiel:
library(ggpubr)
my_comparisons = list( c("0.5", "1"), c("1", "2"), c("0.5", "2") )
ggboxplot(ToothGrowth, x = "dose", y = "len",
color = "dose", palette = "jco")+
stat_compare_means(comparisons = my_comparisons, label.y = c(29, 35, 40))+
stat_compare_means(label.y = 45)
Meine eigene Funktion gemacht:
ts_test <- function(dataL,x,y,method="t.test",idCol=NULL,paired=F,label = "p.signif",p.adjust.method="none",alternative = c("two.sided", "less", "greater"),...) {
options(scipen = 999)
annoList <- list()
setDT(dataL)
if(paired) {
allSubs <- dataL[,.SD,.SDcols=idCol] %>% na.omit %>% unique
dataL <- dataL[,merge(.SD,allSubs,by=idCol,all=T),by=x] #idCol!!!
}
if(method =="t.test") {
dataA <- eval(parse(text=paste0(
"dataL[,.(",as.name(y),"=mean(get(y),na.rm=T),sd=sd(get(y),na.rm=T)),by=x] %>% setDF"
)))
res<-pairwise.t.test(x=dataL[[y]], g=dataL[[x]], p.adjust.method = p.adjust.method,
pool.sd = !paired, paired = paired,
alternative = alternative, ...)
}
if(method =="wilcox.test") {
dataA <- eval(parse(text=paste0(
"dataL[,.(",as.name(y),"=median(get(y),na.rm=T),sd=IQR(get(y),na.rm=T,type=6)),by=x] %>% setDF"
)))
res<-pairwise.wilcox.test(x=dataL[[y]], g=dataL[[x]], p.adjust.method = p.adjust.method,
paired = paired, ...)
}
#Output the groups
res$p.value %>% dimnames %>% {paste(.[[2]],.[[1]],sep="_")} %>% cat("Groups ",.)
#Make annotations ready
annoList[["label"]] <- res$p.value %>% diag %>% round(5)
if(!is.null(label)) {
if(label == "p.signif"){
annoList[["label"]] %<>% cut(.,breaks = c(-0.1, 0.0001, 0.001, 0.01, 0.05, 1),
labels = c("****", "***", "**", "*", "ns")) %>% as.character
}
}
annoList[["x"]] <- dataA[[x]] %>% {diff(.)/2 + .[-length(.)]}
annoList[["y"]] <- {dataA[[y]] + dataA[["sd"]]} %>% {pmax(lag(.), .)} %>% na.omit
#Make plot
coli="#0099ff";sizei=1.3
p <-ggplot(dataA, aes(x=get(x), y=get(y))) +
geom_errorbar(aes(ymin=len-sd, ymax=len+sd),width=.1,color=coli,size=sizei) +
geom_line(color=coli,size=sizei) + geom_point(color=coli,size=sizei) +
scale_color_brewer(palette="Paired") + theme_minimal() +
xlab(x) + ylab(y) + ggtitle("title","subtitle")
#Annotate significances
p <-p + annotate("text", x = annoList[["x"]], y = annoList[["y"]], label = annoList[["label"]])
return(p)
}
library(ggplot2);library(data.table);library(magrittr);
df_long <- rbind(ToothGrowth[,-2],data.frame(len=40:50,dose=3.0))
df_long$ID <- data.table::rowid(df_long$dose)
ts_test(dataL=df_long,x="dose",y="len",idCol="ID",method="wilcox.test",paired=T)
Ich fand dieses ist nützlich.
library(ggplot2)
library(ggpval)
data("PlantGrowth")
plt <- ggplot(PlantGrowth, aes(group, weight)) +
geom_boxplot()
add_pval(plt, pairs = list(c(1, 3)), test='wilcox.test')