Add stats to the plot
1
0
Entering edit mode
6 months ago
G.S ▴ 60

Hello,

I would like to add stats to this plot. Any help would be much appreciated and thanks in advance

Here is my code:

data4.ts%>%  ungroup() %>% t.test(data =.,value ~ Condition) %>%
  adjust_pvalue(method = "bonferroni") 

data2.ts %>%  ungroup() %>% t.test(data =.,value ~ Condition) %>%
  adjust_pvalue(method = "bonferroni") 

a %>% 
  filter(name == "norm_CntTs",
         !is.na(coding)) %>% 
  ggplot(aes(x = Condition, y = value, color = Condition))+
  geom_boxplot(width = 0.5,
               show.legend = FALSE)+
  scale_y_log10()+
  scale_color_manual(
    values = c(
      "HRSV" = "#eb4a40",
      "HRSV_RBV" = "#045275"
    )
  )+
  facet_wrap(.~coding,ncol = 1)
+
  labs(
    title = "Ts",
    y = "Average Ts count" 
  )+
  theme_bw()+
  theme(
    text = element_text(size = 12),
    plot.title = element_text(hjust = 0.5),
    #strip.background = element_rect(color =)
    strip.text = element_text(face = "bold")
  )

enter image description here

R • 1.1k views
ADD COMMENT
0
Entering edit mode

Hi, you can add statistical annotations using ggpubr function stat_pvalue_manual. You need to make sure you calculate the p-values first and then use the stat_pvalue_manual function to add these values to the plot. For example:

data4_test <- data4.ts %>%
  ungroup() %>%
  t.test(value ~ Condition) %>%
  broom::tidy() %>%
  mutate(p.adjusted = p.adjust(p.value, method = "bonferroni"),
         group1 = levels(data4.ts$Condition)[1],
         group2 = levels(data4.ts$Condition)[2])

data2_test <- data2.ts %>%
  ungroup() %>%
  t.test(value ~ Condition) %>%
  broom::tidy() %>%
  mutate(p.adjusted = p.adjust(p.value, method = "bonferroni"),
         group1 = levels(data2.ts$Condition)[1],
         group2 = levels(data2.ts$Condition)[2])

stat_test <- bind_rows(data4_test, data2_test)

a_filtered <- a %>%
  filter(name == "norm_CntTs", !is.na(coding))

plot <- ggplot(a_filtered, aes(x = Condition, y = value, color = Condition)) +
  geom_boxplot(width = 0.5, show.legend = FALSE) +
  scale_y_log10() +
  scale_color_manual(values = c("HRSV" = "#eb4a40", "HRSV_RBV" = "#045275")) +
  facet_wrap(. ~ coding, ncol = 1) +
  labs(title = "Ts", y = "Average Ts count") +
  theme_bw() +
  theme(text = element_text(size = 12),
        plot.title = element_text(hjust = 0.5),
        strip.text = element_text(face = "bold"))

plot_stats <- plot + 
  stat_pvalue_manual(stat_test, 
                     label = "p.adjusted",
                     y.position = max(a_filtered$value) * 1.1, 
                     label.size = 3.5,
                     tip.length = 0.01)

Check the code

ADD REPLY
0
Entering edit mode

Thanks. That helpful.

I am getting this error

data4_test <- data4.ts%>%  ungroup() %>% t.test(data =.,value ~ Condition)%>%
+   broom::tidy() %>%
+   mutate(p.adjusted = p.adjust(p.value, method = "bonferroni"),
+          group1 = levels(data4.ts$Condition)[1],
+          group2 = levels(data4.ts$Condition)[2])
> 
> data2_test <-data2.ts %>%  ungroup() %>% t
.test(data =.,value ~ Condition) %>%
+   broom::tidy() %>%
+   mutate(p.adjusted = p.adjust(p.value, method = "bonferroni"),
+          group1 = levels(data2.ts$Condition)[1],
+          group2 = levels(data2.ts$Condition)[2])
> 
> stat_test <- bind_rows(data4_test, data2_test)
> 
> a_filtered <- a %>%
+   filter(name == "norm_CntTs", !is.na(coding))
> 
> plot <- ggplot(a_filtered, aes(x = Condition, y = value, color = Condition)) +
+   geom_boxplot(width = 0.5, show.legend = FALSE) +
+   scale_y_log10() +
+   scale_color_manual(values = c("HRSV" = "#eb4a40", "HRSV_RBV" = "#045275")) +
+   facet_wrap(. ~ coding, ncol = 1) +
+   labs(title = "Ts", y = "Average Ts count") +
+   theme_bw() +
+   theme(text = element_text(size = 12),
+         plot.title = element_text(hjust = 0.5),
+         strip.text = element_text(face = "bold"))
> 
> plot_stats <- plot + 
+   stat_pvalue_manual(stat_test, 
+                      label = "p.adjusted",
+                      y.position = max(a_filtered$value) * 1.1, 
+                      label.size = 3.5,
+                      tip.length = 0.01)



Error in asserttat_group_columns_exists(data) : 
  data should contain group1 and group2 columns
ADD REPLY
1
Entering edit mode

Please use 101010 to format code so it is represented in monospace font. I have done this for you now.

ADD REPLY
0
Entering edit mode

I think this is what caused the error. we do not have group 1 and 2 in the statistical test results????

enter image description here

ADD REPLY
0
Entering edit mode

The error is likely due to the fact that the structure of your data in data4.ts and data2.ts may not contain the Condition values. The group1 and group2 columns are added based on this. The Condition column should be a factor or a categorical variable that distinguishes between the groups you are comparing. Could you show me data4.ts and data2.ts?

ADD REPLY
0
Entering edit mode

Thanks, thats really helpful

enter image description here

enter image description here

ADD REPLY
0
Entering edit mode

Thank you, it seems to be clearer. Your Condition already has two factorial levels. You can replace these names in groups 1 and 2.

data4_test <- data4.ts %>%
  ungroup() %>%
  t.test(value ~ Condition) %>%
  broom::tidy() %>%
  mutate(p.adjusted = p.adjust(p.value, method = "bonferroni"),
         group1 = "HRSV",
         group2 = "HRSV_RBV")

data2_test <- data2.ts %>%
  ungroup() %>%
  t.test(value ~ Condition) %>%
  broom::tidy() %>%
  mutate(p.adjusted = p.adjust(p.value, method = "bonferroni"),
         group1 = "HRSV",
         group2 = "HRSV_RBV")

stat_test <- bind_rows(data4_test, data2_test)

a_filtered <- a %>%
  filter(name == "norm_CntTs", !is.na(coding))

plot <- ggplot(a_filtered, aes(x = Condition, y = value, color = Condition)) +
  geom_boxplot(width = 0.5, show.legend = FALSE) +
  scale_y_log10() +
  scale_color_manual(values = c("HRSV" = "#eb4a40", "HRSV_RBV" = "#045275")) +
  facet_wrap(. ~ coding, ncol = 1) +
  labs(title = "Ts", y = "Average Ts count") +
  theme_bw() +
  theme(text = element_text(size = 12),
        plot.title = element_text(hjust = 0.5),
        strip.text = element_text(face = "bold"))

plot_stats <- plot + 
  stat_pvalue_manual(stat_test, 
                     label = "p.adjusted",
                     y.position = max(a_filtered$value) * 1.1, 
                     label.size = 3.5,
                     tip.length = 0.01)

Try and let me know

ADD REPLY
0
Entering edit mode

if don't work try binding the stats without group.

stat_data4 <- data4.ts %>%
  ungroup() %>%
  t.test(data = ., value ~ Condition) %>%
  broom::tidy() %>%
  mutate(p.adj = p.adjust(p.value, method = "bonferroni"))

stat_data2 <- data2.ts %>%
  ungroup() %>%
  t.test(data = ., value ~ Condition) %>%
  broom::tidy() %>%
  mutate(p.adj = p.adjust(p.value, method = "bonferroni"))
stat_data <- bind_rows(
  stat_data4 %>% mutate(dataset = "data4"),
  stat_data2 %>% mutate(dataset = "data2")
)

your plot + stat_pvalue_manual(stat_data, label = "p.adj", tip.length = 0.01)

the other parts of the code remain the same. I hope this work, let me know.

ADD REPLY
0
Entering edit mode

Thanks for the help. Unfortunately, did not work :(

error in asserttat_group_columns_exists(data) : data should contain group1 and group2 columns

ADD REPLY
0
Entering edit mode
6 months ago
G.S ▴ 60

Finally, It works. I have added this code

my_comparisons=list(c("HRSV", "HRSV_RBV"))

plot + stat_compare_means(method = "wilcox.test",comparisons = my_comparisons, label.y = 4.5, size= 4, hjust= 0.09,vjust=0.7, alternative = "two.sided",
                          paired = F, exact = F, label = "p.value")

enter image description here

ADD COMMENT

Login before adding your answer.

Traffic: 1375 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6