Samples cross QC (Reads above 2000)

data = read.csv("/pita/users/rotem/eyes/human/Feaces/16s_map.tsv",sep = "\t")

Total number of families that passed QC

data %>% 
  filter(reads_number > 2000) %>% 
  select(Family) %>% 
  unique() %>% 
  nrow()
## [1] 26

PCoA not included the families from old DB3

pca = read_core_pcoa("/pita/users/rotem/eyes/human/Feaces/res/core-metrics-results/unweighted_unifrac_pcoa_results.qza")
names(pca) = paste0("PCoA_", 1:length(names(pca)))

faith = read_faith_qzv("res/core-metrics-results/faith-pd-correlation.qzv")
names(faith)[1] =  "sampleid"
faith$sampleid = faith$sampleid %>% as.character()
sum(data$reads_number > 2000)
## [1] 135
data = 
data %>% 
  left_join(pca %>% rownames_to_column("sampleid")) 
## Joining, by = "sampleid"
## Warning: Column `sampleid` joining factor and character vector, coercing into
## character vector
data = data %>% left_join(faith %>% select(sampleid, faith_pd), by = names(faith)[1])
create_dt(data)
pl = 
ggplot(data, aes(x = PCoA_1, y = PCoA_2, color = Age, name = sample_ID)) + 
  geom_point() + 
  theme_rh
plotly::ggplotly(pl)
pl = 
data %>% 
  ggplot(aes(x = Age, y = faith_pd, fill = Family, color = Family)) +
  geom_point(show.legend = F) + 
  theme_rh + theme(axis.text.x = element_text(angle = -45, hjust = 0))
plotly::ggplotly(pl)
# data %>% 
#   ggplot(aes(fill = Retinal_degeneration_type, x = Age)) + 
#   geom_histogram(stat = 'count') + theme_rh
library(forcats)
pl =
data %>%
  filter(grepl("FR.*", Family)) %>% 
  arrange(Date_of_diagnosis) %>% 
  mutate(sample_ID = factor(sample_ID, levels = sample_ID)) %>% 
  mutate(Date_of_diagnosis = replace_na(Date_of_diagnosis, "Healthy")) %>%
  group_by(Family) %>%
  mutate(  index = min(Date_of_diagnosis, na.rm = T)
       , family_disease = paste(Family, na.omit(Retinal_degeneration_type))) %>% 
  ggplot(aes(
      x = sample_ID
    , y = Age
    , sample = sample_ID)) +
  geom_point() +
  geom_point(aes(  y = Age_on_diagnosis
                 , fill = Retinal_degeneration_type
                 , color = Retinal_degeneration_type)) +
  geom_linerange(aes(
                   x = sample_ID
                   , ymin = Age_on_diagnosis
                   , ymax = Age
                   , color = Retinal_degeneration_type)) +
  # geom_point(data = data %>% group_by(Family) %>%
  #   slice(which.min(Date_of_diagnosis)), aes(size = 1/Date_of_diagnosis), color = "red") +
  # facet_wrap(~Date_of_diagnosis)+
  
  facet_wrap(~family_disease, nrow = 4, scales = 'free_x') +
  theme_rh + 
  theme(  strip.background = element_blank()
        , strip.text = element_text(color = "red")
        , axis.title.x = element_blank()
        , axis.text.x = element_blank())
# pl
plotly::ggplotly(pl, width = 1000, height = 400)
pl = 
ggplot(data, aes(x = PCoA_1, y = PCoA_2, color = Family, name = sample_ID)) + 
  geom_point() + 
  theme_rh

plotly::ggplotly(pl)
pl = 
ggplot(data, aes(x = PCoA_1, y = PCoA_2, color = Gender, name = sample_ID)) + 
  geom_point() + 
  theme_rh
plotly::ggplotly(pl)
# faith = read_faith_qzv("res/core-metrics-results/faith-pd-correlation.qzv")
names(faith)[1] =  "sampleid"
faith$sampleid = faith$sampleid %>% as.character()

data = data %>% left_join(faith %>% select(sampleid, contains("faith_pd")), by = names(faith)[1])

pl = 
data %>% 
  ggplot(aes(x = Family, y = faith_pd.x, fill = Family)) +
  geom_boxplot(show.legend = F) + 
  theme_rh + theme(axis.text.x = element_text(angle = -45, hjust = 0))
pl
## Warning: Removed 1 rows containing non-finite values (stat_boxplot).

library(ggpubr)
## Loading required package: magrittr
## 
## Attaching package: 'magrittr'
## The following object is masked from 'package:purrr':
## 
##     set_names
## The following object is masked from 'package:tidyr':
## 
##     extract
families_with_healthy_and_ill = 
data %>% 
    filter(Age > 20) %>% 
    filter(grepl("FR.*", Family)) %>% 
    group_by(Family) %>% 
    summarise(n = n()) %>% 
    filter(n > 1) %>% pull(Family)

pl = 
data %>% 
    filter(Family %in% families_with_healthy_and_ill) %>% 
    filter(Age > 20) %>% 
  ggplot(aes(x = Retinal_degeneration_type, fill = Retinal_degeneration_type
             , y = faith_pd.x)) + 
  geom_boxplot() + 
  geom_jitter(alpha = .5)+
  facet_wrap(~Family, scales="free") + 
  # stat_compare_means(method = "t.test", label = "p.signif") +
  theme_rh +
  ggtitle("Alpha Diversity within familis above 20") +
  theme(  strip.background = element_blank())
# pl

plotly::ggplotly(pl, width = 1000)
my_comparisons = list(c("healthy","stargardt"), c("Healthy", "RP"))
my_comparisons = list(c(1,2), c(1,3))
data %>% 
    filter(Family %in% families_with_healthy_and_ill) %>% 
    filter (Age > 20) %>% 
    mutate(Retinal_degeneration_type = replace_na(Retinal_degeneration_type %>% as.character(), "Healthy")) %>% 
    filter(Retinal_degeneration_type %in% c('Healthy', 'RP', 'stargardt')) %>%
  ggplot(aes(x = Retinal_degeneration_type, fill = Retinal_degeneration_type
             , y = faith_pd.x)) + 
  geom_boxplot() + 
  geom_jitter(alpha = .5)+
  stat_compare_means(comparisons = my_comparisons) +
  theme_rh +
  ggtitle("Alpha Diversity within diseases, age above 20") +
  theme(  strip.background = element_blank())