genome_counts_filt %>%
mutate_at(vars(-genome),~./sum(.)) %>% #apply TSS nornalisation
pivot_longer(-genome, names_to = "sample", values_to = "count") %>% #reduce to minimum number of columns
left_join(., genome_metadata, by = join_by(genome == genome)) %>% #append genome metadata
left_join(., sample_metadata, by = join_by(sample == sample)) %>% #append sample metadata
filter(!is.na(count)) %>%
ggplot(aes(y=count,x=sample, fill=phylum, group=phylum)) + #grouping enables keeping the same sorting of taxonomic units
geom_bar(stat="identity", colour="white", linewidth=0.1) +
scale_fill_manual(values=phylum_colors) +
facet_grid(~ Bd_status, scale="free", space="free")+
guides(fill = guide_legend(ncol = 1)) +
theme(
axis.title.x = element_blank(),
panel.background = element_blank(),
panel.border = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
strip.background = element_rect(fill = "white"),
strip.text = element_text(size = 8, lineheight = 0.6),
strip.placement = "outside",
axis.text.x = element_blank(), axis.ticks.x = element_blank(),
axis.line = element_line(linewidth = 0.5, linetype = "solid", colour = "black")) +
labs(fill="Phylum",y = "Relative abundance",x="Samples")
Phylum relative abundances
phylum_summary <- genome_counts_filt %>%
mutate_at(vars(-genome),~./sum(.)) %>%
pivot_longer(-genome, names_to = "sample", values_to = "count") %>%
left_join(sample_metadata, by = join_by(sample == sample)) %>%
left_join(genome_metadata, by = join_by(genome == genome)) %>%
group_by(sample,phylum,Bd_status) %>%
summarise(relabun=sum(count))
phylum_summary %>%
group_by(phylum) %>%
summarise(
total_mean = mean(relabun * 100, na.rm = T),
total_sd = sd(relabun * 100, na.rm = T),
positive_mean = mean(relabun[Bd_status == "Bd_positive"] * 100, na.rm = T),
positive_sd = sd(relabun[Bd_status == "Bd_positive"] * 100, na.rm = T),
negative_mean = mean(relabun[Bd_status == "Bd_negative"] * 100, na.rm = T),
negative_sd = sd(relabun[Bd_status == "Bd_negative"] * 100, na.rm = T)
) %>%
mutate(
total = str_c(round(total_mean, 3), "±", round(total_sd, 3)),
Positive = str_c(round(positive_mean, 3), "±", round(positive_sd, 3)),
Negative = str_c(round(negative_mean, 3), "±", round(negative_sd, 3))
) %>%
arrange(-total_mean) %>%
dplyr::select(phylum, total, Positive, Negative) %>%
tt()
| phylum |
total |
Positive |
Negative |
| p__Pseudomonadota |
73.437±23.223 |
70.928±18.007 |
74.552±25.433 |
| p__Actinomycetota |
6.249±13.878 |
3.852±7.991 |
7.314±15.832 |
| p__Bacteroidota |
5.885±7.341 |
9.221±10.595 |
4.402±4.883 |
| p__Verrucomicrobiota |
5.635±16 |
2.96±5.359 |
6.824±18.9 |
| p__Patescibacteria |
5.435±14.284 |
7.971±13.383 |
4.308±14.77 |
| p__Cyanobacteriota |
1.422±5.766 |
2.147±7.203 |
1.1±5.126 |
| p__Desulfobacterota_F |
0.416±1.031 |
0.671±1.538 |
0.303±0.713 |
| p__Campylobacterota |
0.336±1.012 |
0.714±1.617 |
0.168±0.545 |
| p__Bacillota |
0.311±1.118 |
0.573±1.879 |
0.195±0.536 |
| p__Chloroflexota |
0.247±0.875 |
0.196±0.68 |
0.269±0.96 |
| p__Gemmatimonadota |
0.216±0.603 |
0.055±0.19 |
0.288±0.706 |
| p__Thermoproteota |
0.208±0.908 |
0.328±1.138 |
0.155±0.805 |
| p__Nitrospirota |
0.071±0.443 |
0.23±0.798 |
0±0 |
| p__Halobacteriota |
0.046±0.232 |
0±0 |
0.067±0.277 |
| p__Bdellovibrionota |
0.037±0.234 |
0±0 |
0.054±0.281 |
| p__Acidobacteriota |
0.025±0.157 |
0.082±0.284 |
0±0 |
| p__Bacillota_A |
0.022±0.136 |
0.071±0.245 |
0±0 |