I am in the middle of making my heatmap, I change several script to suit my data.
Now I am wondering the possibility of changing the order of the facet. Could you change the order of the body_site of the facet?
For example, in your attached figure, they were ordered like gut;left palm;right palm;tongue. What should be added to the script/command to change it in custom order? E.g right palm;left palm; tongue;gut ??
I try the factor() function with levels the “body_site” but error that dollar sign are not accepted in atomic vector appeared.
library(tidyverse)
library(qiime2R)
metadata<-read_q2metadata("~/GitHub/qiime2R/inst/artifacts/2020.2_moving-pictures/sample-metadata.tsv")
SVs<-read_qza("~/GitHub/qiime2R/inst/artifacts/2020.2_moving-pictures/table.qza")$data
taxonomy<-read_qza("~/GitHub/qiime2R/inst/artifacts/2020.2_moving-pictures/taxonomy.qza")$data
SVs<-apply(SVs, 2, function(x) x/sum(x)*100) #convert to percent
SVsToPlot<-
data.frame(MeanAbundance=rowMeans(SVs)) %>% #find the average abundance of a SV
rownames_to_column("Feature.ID") %>%
arrange(desc(MeanAbundance)) %>%
top_n(30, MeanAbundance) %>%
pull(Feature.ID) #extract only the names from the table
SVs %>%
as.data.frame() %>%
rownames_to_column("Feature.ID") %>%
gather(-Feature.ID, key="SampleID", value="Abundance") %>%
mutate(Feature.ID=if_else(Feature.ID %in% SVsToPlot, Feature.ID, "Remainder")) %>% #flag features to be collapsed
group_by(SampleID, Feature.ID) %>%
summarize(Abundance=sum(Abundance)) %>%
left_join(metadata) %>%
mutate(NormAbundance=log10(Abundance+0.01)) %>% # do a log10 transformation after adding a 0.01% pseudocount. Could also add 1 read before transformation to percent
left_join(taxonomy) %>%
mutate(Feature=paste(Feature.ID, Taxon)) %>%
mutate(Feature=gsub("[kpcofgs]__", "", Feature)) %>% # trim out leading text from taxonomy string
mutate(`body-site`=factor(`body-site`, levels=c("right palm","left palm","gut","tongue"))) %>%
ggplot(aes(x=SampleID, y=Feature, fill=NormAbundance)) +
geom_tile() +
facet_grid(~`body-site`, scales="free_x") +
theme_q2r() +
theme(axis.text.x=element_text(angle=45, hjust=1)) +
scale_fill_viridis_c(name="log10(% Abundance)")
ggsave("heatmap.pdf", height=4, width=11, device="pdf") # save a PDF 3 inches by 4 inches
mutate(SVs$body-site=factor(SVs$body-site, levels=c("right palm","left palm","gut","tongue")))
something like that, and I dont know about where to put the code. Thanks