Сломанный t-тест в facet_wrap

Я пытаюсь одновременно использовать facet_wrap и stat_compare_means, но у меня проблема. Две стороны данных не имеют одинакового количества точек. Следовательно, stat_compare_means не работает ... Посмотрите, например, на изображение:

введите описание изображения здесь

в Type1 B есть три точки, тогда как в Type2 B только одна точка. Это несоответствие приводит к тому, что почти все t-тест не проходит и не отображается на графике. Что мне нужно, так это t-тест для групп, у которых есть совпадающее количество точек (в данном случае все t-тест на Type1 и A против C в Type2). Используемый сюжет следующий:

library(RColorBrewer)
library(ggpubr)
library(BBmisc)


adf=read.csv("test1.txt", sep=" ")

myColors <- brewer.pal(length(unique(adf$ID)) ,"Set1")
names(myColors) <- unique(adf$ID)
colScale <- scale_colour_manual(name = "ID",values = myColors)

my_comparison=as.data.frame(combn(unique(adf$sampletype) ,2))
my_comparison=convertColsToList(my_comparison)

ggplot(adf, aes(x=sampletype, y=value, fill=sampletype ))+
  geom_point(aes(group=ID, colour=ID))+
  geom_line(aes(group=ID, colour=ID))+
  facet_wrap(~response, scale="free")+
  colScale+
  ggtitle("Entropy")+
  theme(text = element_text(size=20))+
  stat_compare_means(comparisons = my_comparison, method = "t.test", paired = TRUE)

Данные (сохранены как test1.txt):

sampletype value ID response
A 8.192 gr_6 Type2
B 13.99 gr_6 Type2
C 9.186 gr_6 Type2
A 5.616 gr_5 Type1
B 15.55 gr_5 Type1
C 7.126 gr_5 Type1
A 5.484 gr_4 Type1
B 12.54 gr_4 Type1
C 4.492 gr_4 Type1
A 9.949 gr_3 Type2
C 6.631 gr_3 Type2
A 2.533 gr_7 Type2
C 12.25 gr_7 Type2
A 2.196 gr_2 Type2
C 6.447 gr_2 Type2
A 11.20 gr_1 Type1
B 16.63 gr_1 Type1
C 6.637 gr_1 Type1

Есть ли обходной путь?


person Fabrizio    schedule 28.09.2020    source источник
comment
Вы пробовали построить участок, используя только ggpubr, а не добавлять его к графику, построенному с помощью ggplot? См. Здесь stackoverflow.com/questions/64051045/   -  person NColl    schedule 28.09.2020
comment
Но с этим ответом мне придется построить график ggpubr для каждого ответа вместо использования facet_wrap. Это может решить мою проблему, но на этом этапе я также могу разделить данные и продолжить работу с ggpplot. Или не?   -  person Fabrizio    schedule 28.09.2020
comment
Вы можете выполнить фасетирование в ggpubr, добавив facet.by= к первому вызову ggpubr или добавив его отдельно, как описано здесь rpkgs.datanovia.com/ggpubr/reference/facet.html   -  person NColl    schedule 28.09.2020
comment
Я попытался запустить код на странице, так как я не очень хорошо знаю ggpubr, я оказался в: Ошибка: Проблема с mutate() вводом data. x _tidyr_melt_dataframe не разрешен из текущего пространства имен (tidyr) ℹ Вход data - map(.data$data, .f, ...). Запустите rlang::last_error(), чтобы узнать, где произошла ошибка.   -  person Fabrizio    schedule 28.09.2020


Ответы (2)


У меня есть способ сделать это, но предоставленные данные не прошли t-тест из-за отсутствия дисперсии, поэтому я изменил их.

library(tidyverse)
library(rstatix)
library(ggpubr)
df <- enframe(c("A 8.192 gr_6 Type2",
"B 13.99 gr_6 Type2",
"C 9.186 gr_6 Type2",
"A 5.616 gr_5 Type1",
"B 15.55 gr_5 Type1",
"C 7.126 gr_5 Type1",
"A 5.484 gr_4 Type1",
"B 12.54 gr_4 Type1",
"C 4.492 gr_4 Type1",
"A 9.949 gr_3 Type2",
"C 6.631 gr_3 Type2",
"A 2.533 gr_7 Type2",
"C 12.25 gr_7 Type2",
"A 2.196 gr_2 Type2",
"C 6.447 gr_2 Type2",
"A 11.20 gr_1 Type1",
"B 16.63 gr_1 Type1",
"C 6.637 gr_1 Type1"))
df <- df %>%
  separate(value, into =c("sampletype", "value", "ID", "response"), 
           sep=" ") %>% select(-name) %>%
  mutate(
    value = sample(1:50, 18)
  )

keep_vars <- df %>%
  group_by(response, sampletype) %>% tally() %>% filter(n>1) %>% 
  pivot_wider(names_from = sampletype, values_from = n) %>%split(.$response) %>%
  map(pivot_longer, cols=-c(response)) %>% map(filter,value>=0) %>% bind_rows(.) %>%
  mutate(UID = paste0(response, name)) %>%pull(UID)
  
df_plot <- df %>% mutate(UID = paste0(response, sampletype))  %>%
  filter(UID %in% keep_vars) %>%
  group_by(response) %>%
  t_test(value~sampletype) %>% add_xy_position(x='sampletype')

ggpubr::ggline(df, 
                  x='sampletype', 
                  y='value',
                  color = 'ID',
                  add='jitter',
                  facet.by = 'response'
                  ) +
  stat_pvalue_manual(df_plot, 
                     label = "p.adj.signif", 
                     tip.length = 0.01)

введите описание изображения здесь

person NColl    schedule 28.09.2020
comment
Значения парные, поэтому при использовании t_test (value ~ sampletype, paired = TRUE) вам не нужно изменять значения. Я думаю, что ваше решение работает, позвольте мне немного попробовать, но спасибо за ваши усилия. - person Fabrizio; 29.09.2020

Решение NColl отличное, я тоже его поддержал, но в моем конкретном случае он не удался. В частности, когда есть одно значение для одного из случаев (например, одно значение в Type2, A, gr_X). Это моя вина, в этом примере нет той общности, которую я хотел. Тем не менее, Ncoll заставляет меня понять, что мне нужно смотреть на переменные, которые я хочу сравнить. Я не хотел больше беспокоиться, поэтому реализовал свое решение. На практике я рисую два ggplot и объединяю их с gridExtra:

adf_r=adf[adf$sampletype %in% ctime  & adf$response=="Responders", ]
adf_nr=adf[adf$sampletype %in% ctime &d adf$response=="Non-responders", ]

adrmyColors <- brewer.pal(length(unique(adf$ID)) ,"Set1")
names(myColors) <- unique(adf$ID)
colScale <- scale_colour_manual(name = "ID",values = myColors)


for(st in unique(adf_r$sampletype)){
  if(sum(adf_r$sampletype==st)>=3){
    keepvar1=c(keepvar1,st)
  }
}

my_comparison=NULL    
if(length(keepvar1)>1){
  my_comparison=as.data.frame(combn(keepvar1,2))
  my_comparison=convertColsToList(my_comparison)
}

p1=ggplot(adf_r, aes(x=sampletype, y=value, fill=sampletype ))+
  geom_point(aes(group=ID, colour=ID))+
  geom_line(aes(group=ID, colour=ID))+
   colScale+
   ggtitle("Entropy")+
   theme(text = element_text(size=20))+
   guides(fill=FALSE)+
   stat_compare_means(comparisons = my_comparison, method = "t.test", paired = TRUE)


keepvar2=NULL
for(st in unique(adf_nr$sampletype)){
  if(sum(adf_nr$sampletype==st)>=3){
    keepvar2=c(keepvar2,st)
  }
}

my_comparison=NULL
if(length(keepvar2)>1){
  my_comparison=as.data.frame(combn(keepvar2,2))
  my_comparison=convertColsToList(my_comparison)
}

p2=ggplot(adf_nr, aes(x=sampletype, y=value, fill=sampletype ))+
  geom_point(aes(group=ID, colour=ID))+
  geom_line(aes(group=ID, colour=ID))+
  colScale+
  ggtitle("Entropy")+
  guides(fill=FALSE)+
  theme(text = element_text(size=20))+
  stat_compare_means(comparisons = my_comparison, method = "t.test", paired = TRUE)

final_p=gridExtra::grid.arrange(p1, p2, ncol = 2)
person Fabrizio    schedule 01.10.2020