{"type":"doc","content":[{"type":"extension","attrs":{"layout":"default","extensionType":"com.atlassian.confluence.macro.core","extensionKey":"toc","parameters":{"macroParams":{"style":{"value":"none"}},"macroMetadata":{"macroId":{"value":"c7230970-1314-40b7-86f5-7f4e0eb6238c"},"schemaVersion":{"value":"1"},"title":"Table of Contents"}},"localId":"cebddf3a-629e-4040-8582-dd55853da3d0"}},{"type":"paragraph"},{"type":"heading","attrs":{"level":1},"content":[{"text":"Overview","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"This section examines the community structure - i.e. the proportion, variance and abundance of taxonomic groups found in each sample and treatment group.","type":"text"}]},{"type":"paragraph","content":[{"text":"Amplicon sequence variants","type":"text","marks":[{"type":"strong"}]}]},{"type":"paragraph","content":[{"text":"Taxonomic groups in this study are based on Amplicon sequence variants (ASV), inferred using the ","type":"text"},{"text":"DADA2 software package","type":"text","marks":[{"type":"strong"},{"type":"link","attrs":{"href":"https://www.nature.com/articles/nmeth.3869"}}]},{"text":" by matching the sample sequences to the ","type":"text"},{"text":"Silva ribosomal RNA sequence database","type":"text","marks":[{"type":"strong"},{"type":"link","attrs":{"href":"https://www.arb-silva.de/documentation/"}}]},{"text":".","type":"text"}]},{"type":"paragraph","content":[{"text":"DADA2","type":"text","marks":[{"type":"strong"}]},{"text":" infers sample sequences exactly and resolves differences of as little as 1 nucleotide.","type":"text"}]},{"type":"paragraph","content":[{"text":"SILVA","type":"text","marks":[{"type":"strong"}]},{"text":" provides comprehensive, quality checked and regularly updated datasets of aligned small (16S/18S, SSU) and large subunit (23S/28S, LSU) ribosomal RNA (rRNA) sequences for all three domains of life (Bacteria, Archaea and Eukarya).","type":"text"}]},{"type":"paragraph"},{"type":"heading","attrs":{"level":1},"content":[{"text":"6a. Results summary","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"You can see a summary of your data by typing the name of your ampvis2 object:","type":"text"}]},{"type":"codeBlock","content":[{"text":"ampvisdata","type":"text"}]},{"type":"paragraph","content":[{"text":"This will provide you a summary of sequence reads per sample (minimum, maximum, average, etc) and the number and percentage of ASVs resolved to each taxonomic level.","type":"text"}]},{"type":"paragraph"},{"type":"heading","attrs":{"level":1},"content":[{"text":"6b. Taxonomic abundance bar plots","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"This section generates stacked bar plots showing the proportion of taxa per sample, for each taxonomic level.","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"First, you can output the raw ASV abundance data as a csv file. You can add this as supplemental data to a manuscript or use it in other plotting or statistical programs:","type":"text"}]},{"type":"codeBlock","content":[{"text":"amp_export_otutable(ampvisdata, filename = \"ASVtable\", normalise = T, sep = \",\")","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Now create the bar plots.","type":"text"}]},{"type":"paragraph","content":[{"text":"Choose the taxonomic level you want to plot (Choose from Phylum, Class, Order, Family, Genus, Species)","type":"text"}]},{"type":"codeBlock","content":[{"text":"taxgroup <- \"Class\"","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Generate the data for plotting. > 1% relative abundance only","type":"text"}]},{"type":"codeBlock","content":[{"text":"# Run amp_boxplot to generate aggregated data for the chosen taxonomic level, but output as a data object rather than plotting\nx <- amp_boxplot(ampvisdata, tax_aggregate = taxgroup, tax_empty = \"remove\")\n# Filter out any taxa with > 3% relative abundance\n# First calculate mean on aggregated taxa\nx2 <- aggregate(x$data$Abundance, list(x$data$Display), mean)\n# Select taxa only with abundance > 1%\nx3 <- x2[x2$x > 1,]\n# Pull out just those taxa in x$data\nx4 <- x$data[x$data$Display %in% x3$Group.1, ]\n# Remove duplicate rows (means already calculated, duplicates only cause bias)\nx5 <- x4[!duplicated(x4), ]","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Generate a stacked bar plot.","type":"text"}]},{"type":"paragraph","content":[{"text":"As with previous plots, you can modify colours using ","type":"text"},{"text":"scale_fill_manual()","type":"text","marks":[{"type":"code"}]},{"text":".","type":"text"}]},{"type":"paragraph","content":[{"text":"You can also change axis labels, text size and angle, ","type":"text"},{"text":"default theme","type":"text","marks":[{"type":"link","attrs":{"href":"https://ggplot2.tidyverse.org/reference/ggtheme.html"}}]},{"text":", etc.","type":"text"}]},{"type":"paragraph","content":[{"text":"Feel free to add additional modifications to this or any other plot in this Notebook. Here is a good guide for doing this: ","type":"text"},{"type":"inlineCard","attrs":{"url":"http://r-statistics.co/Complete-Ggplot2-Tutorial-Part2-Customizing-Theme-With-R-Code.html"}},{"text":" ","type":"text"}]},{"type":"codeBlock","content":[{"text":"p <- ggplot(x5, aes(x = Sample, y = Abundance, fill = Display))\np <- p + geom_bar(position=\"fill\", stat=\"identity\") +\n labs(y = \"Abundance (%)\") +\n theme_classic() +\n scale_y_continuous(label = label_percent()) + \n scale_fill_manual(values = c(\"red\", \"blue\", \"green\")) + \n theme(text = element_text(size = 18), axis.text.x = element_text(angle = 90, size=12), axis.text.y = element_text(size=14)) + \n labs(fill = taxgroup)\np","type":"text"}]},{"type":"paragraph","content":[{"text":"To plot a different taxonomic level, You can change ","type":"text"},{"text":"taxgroup <- \"....\"","type":"text","marks":[{"type":"code"}]},{"text":" above to another taxonomic level, then re-run the code from that point.","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Save your plot as a 300dpi (i.e. publication quality) tiff or pdf file","type":"text"}]},{"type":"codeBlock","content":[{"text":"tiff_exp <- paste0(\"community_tax_abund_bar_plot_\", taxgroup, \"_samples_.tiff\")\nggsave(file = tiff_exp, dpi = 300, compression = \"lzw\", device = \"tiff\", plot = p, width = 20, height = 20, units = \"cm\")\npdf_exp <- paste0(\"community_tax_abund_bar_plot_\", taxgroup, \"_samples_.pdf\")\nggsave(file = pdf_exp, device = \"pdf\", plot = p, width = 20, height = 20, units = \"cm\")","type":"text"}]},{"type":"paragraph"},{"type":"heading","attrs":{"level":1},"content":[{"text":"6c. Taxonomic abundance box and whisker plots","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"This section contains box and whisker plots for each taxonomic level (phylum … species).","type":"text"}]},{"type":"paragraph","content":[{"text":"As with the bar plots in the previous section, these B&W plots show the relative abundance (in %) of taxa. But whereas heatmaps show the mean %, box and whisker plots show the interquartile range between samples (","type":"text"},{"text":"max, min, median, quartiles","type":"text","marks":[{"type":"link","attrs":{"href":"https://courses.lumenlearning.com/fscj-introstats1/chapter/box-plots/"}}]},{"text":"), thus providing additional statistical information.","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Choose the taxonomic level you want to plot (Choose from Phylum, Class, Order, Family, Genus, Species)","type":"text"}]},{"type":"codeBlock","content":[{"text":"taxgroup <- \"Order\"","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Generate the plot (and modify the attributes as you like)","type":"text"}]},{"type":"paragraph","content":[{"text":"Note the ","type":"text"},{"text":"tax_show = 10","type":"text","marks":[{"type":"code"}]},{"text":" attribute. This says to show the top 10 taxa. Change this as desired.","type":"text"}]},{"type":"codeBlock","content":[{"text":"p <- amp_boxplot(ampvisdata, tax_aggregate = taxgroup, tax_empty = \"remove\", tax_show = 10) +\ntheme_bw() +\ntheme(text = element_text(size = 18), axis.text.x = element_text(size=16), axis.text.y = element_text(size=14)) +\nlabs(x = taxgroup)\np","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Export the plots as tiff and pdf:","type":"text"}]},{"type":"codeBlock","content":[{"text":"tiff_exp <- paste0(\"community_tax_abund_BW_plot_\", taxgroup, \"_samples_.tiff\")\nggsave(file = tiff_exp, dpi = 300, compression = \"lzw\", device = \"tiff\", plot = p, width = 20, height = 20, units = \"cm\")\npdf_exp <- paste0(\"community_tax_abund_BW_plot_\", taxgroup, \"_samples_.pdf\")\nggsave(file = pdf_exp, device = \"pdf\", plot = p, width = 20, height = 20, units = \"cm\")","type":"text"}]},{"type":"paragraph"},{"type":"heading","attrs":{"level":1},"content":[{"text":"6d. Choosing a variable to analyse","type":"text"}]},{"type":"paragraph","content":[{"text":" ","type":"text"}]},{"type":"paragraph","content":[{"text":"In the Alpha Diversity section, you already imported your nfcore/ampliseq results and converted them to an ampvis2 object. To have another quick look at your ampvis object:","type":"text"}]},{"type":"codeBlock","content":[{"text":"ampvisdata","type":"text"}]},{"type":"paragraph","content":[{"text":" As in the Alpha Diversity section, you need to choose a variable to work with.","type":"text"}]},{"type":"paragraph","content":[{"text":"You can view your variables as column names in your samples_table:","type":"text"}]},{"type":"codeBlock","content":[{"text":"colnames(samples_table)","type":"text"}]},{"type":"paragraph","content":[{"text":" Now enter the variable you want to analyse:","type":"text"}]},{"type":"codeBlock","content":[{"text":"group <- \"Nose_size\"","type":"text"}]},{"type":"paragraph","content":[{"text":" ","type":"text"}]},{"type":"paragraph","content":[{"text":"And once again you’ll need to order the groups in your variable. See the ‘Ordering your variable’ section in the Alpha diversity section for more details.","type":"text"}]},{"type":"paragraph","content":[{"text":"Choose how you want to order your groups:","type":"text"}]},{"type":"codeBlock","content":[{"text":"lev <- c(\"Small\", \"Medium\", \"Big\")","type":"text"}]},{"type":"paragraph","content":[{"text":"Then run the following to apply the levels to your data:","type":"text"}]},{"type":"codeBlock","content":[{"text":"ampvisdata$metadata[[group]] <- factor(ampvisdata$metadata[[group]], levels = lev)","type":"text"}]},{"type":"paragraph","content":[{"text":" ","type":"text"}]},{"type":"heading","attrs":{"level":1},"content":[{"text":"6e. Community structure by variable","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"In this section you can generate heatmaps and box and whisker plots for the variable you selected above, and for each taxonomic level.","type":"text"}]},{"type":"paragraph","content":[{"text":"If you want to analyse another variable here, go back to the previous section, choose another variable and re-run the code from that point.","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Choose the taxonomic level you want to plot (Choose from Phylum, Class, Order, Family, Genus, Species). This applies to both the heatmaps and B&W plots below.","type":"text"}]},{"type":"codeBlock","content":[{"text":"taxgroup <- \"Order\"","type":"text"}]},{"type":"paragraph"},{"type":"heading","attrs":{"level":2},"content":[{"text":"Heatmaps","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"This section contains heatmaps for each taxonomic level (class - species) for your selected variable. The number in each cell is the relative proportion (in %) of taxa per group.","type":"text"}]},{"type":"paragraph","content":[{"text":"There are a variety of attributes you can modify:","type":"text"}]},{"type":"paragraph","content":[{"text":"Note the ","type":"text"},{"text":"tax_show = 10","type":"text","marks":[{"type":"code"}]},{"text":" attribute. This says to show the top 10 taxa.","type":"text"}]},{"type":"paragraph","content":[{"text":"The ","type":"text"},{"text":"plot_values_size =","type":"text","marks":[{"type":"code"}]},{"text":" defines the size of the text in the heatmap cells.","type":"text"}]},{"type":"paragraph","content":[{"text":"The ","type":"text"},{"text":"tax_add = NULL","type":"text","marks":[{"type":"code"}]},{"text":" adds an additional taxonomic group to the taxa names (e.g. changing this to ","type":"text"},{"text":"tax_add = \"Phylum\"","type":"text","marks":[{"type":"code"}]},{"text":" and plotting Genus will name the taxa as 'phylum:genus' but leaving it as ","type":"text"},{"text":"tax_add = NULL","type":"text","marks":[{"type":"code"}]},{"text":" will just give the genus name.","type":"text"}]},{"type":"paragraph","content":[{"text":"showRemainingTaxa = T","type":"text","marks":[{"type":"code"}]},{"text":" will aggregate, in a single row, the remaining taxa not shown on this heatmap. Change to ","type":"text"},{"text":"F","type":"text","marks":[{"type":"code"}]},{"text":" if you don't want to see this.","type":"text"}]},{"type":"paragraph","content":[{"text":"color_vector = NULL","type":"text","marks":[{"type":"code"}]},{"text":" uses the default colour range (orange -> blue). You can change this by providing your own colour range, e.g. ","type":"text"},{"text":"color_vector = c(\"red\", \"green\")","type":"text","marks":[{"type":"code"}]},{"text":". You can choose from the huge number of R colours here: ","type":"text"},{"text":"http://www.stat.columbia.edu/~tzheng/files/Rcolor.pdf","type":"text","marks":[{"type":"link","attrs":{"href":"http://www.stat.columbia.edu/~tzheng/files/Rcolor.pdf"}}]}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Generate the plot:","type":"text"}]},{"type":"codeBlock","content":[{"text":"p <- amp_heatmap(ampvisdata, tax_show = 10, plot_values_size = 4, tax_add = NULL, showRemainingTaxa = T, color_vector = NULL, group_by = group, tax_aggregate = taxgroup, tax_empty = \"remove\")\np","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Add some additional plot modifications. Change or remove these as desired (","type":"text"},{"text":"or add your own modifications","type":"text","marks":[{"type":"link","attrs":{"href":"http://r-statistics.co/Complete-Ggplot2-Tutorial-Part2-Customizing-Theme-With-R-Code.html"}}]},{"text":")","type":"text"}]},{"type":"codeBlock","content":[{"text":"p <- p +\ntheme_bw() +\ntheme(text = element_text(size = 18), axis.text.x = element_text(angle = 90, size=16), axis.text.y = element_text(size=14)) +\nlabs(y = taxgroup)\np","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Export as pdf and tiff","type":"text"}]},{"type":"codeBlock","content":[{"text":"tiff_exp <- paste0(\"community_tax_abund_heatmap_\", group, \"_\", taxgroup, \"_samples_.tiff\")\nggsave(file = tiff_exp, dpi = 300, compression = \"lzw\", device = \"tiff\", plot = p, width = 20, height = 20, units = \"cm\")\npdf_exp <- paste0(\"community_tax_abund_heatmap_\", group, \"_\", taxgroup, \"_samples_.pdf\")\nggsave(file = pdf_exp, device = \"pdf\", plot = p, width = 20, height = 20, units = \"cm\")","type":"text"}]},{"type":"paragraph"},{"type":"heading","attrs":{"level":2},"content":[{"text":"Box and whisker plots","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"This section contains B&W plots for each taxonomic level (class - species) for your selected variable. This is similar to the heatmaps you generated in section 6, but separates the B&W into groups per variable, thus enabling an examination of how the variable groups differ in terms of taxonomic abundance.","type":"text"}]},{"type":"paragraph","content":[{"text":"Again, there are a variety of attributes you can modify:","type":"text"}]},{"type":"paragraph","content":[{"text":"Note the ","type":"text"},{"text":"tax_show = 7","type":"text","marks":[{"type":"code"}]},{"text":" attribute. This says to show the top 7 taxa (","type":"text"},{"text":"Tip:","type":"text","marks":[{"type":"strong"}]},{"text":" to look at specific taxa you can provide a vector of taxa names, e.g. ","type":"text"},{"text":"tax_show = c(\"Bacteroidales\", \"Fibrobacterales\")","type":"text","marks":[{"type":"code"}]},{"text":").","type":"text"}]},{"type":"paragraph","content":[{"text":"The ","type":"text"},{"text":"tax_add = NULL","type":"text","marks":[{"type":"code"}]},{"text":" adds an additional taxonomic group to the taxa names (e.g. changing this to ","type":"text"},{"text":"tax_add = \"Phylum\"","type":"text","marks":[{"type":"code"}]},{"text":" and plotting Genus will name the taxa as 'phylum:genus' but leaving it as ","type":"text"},{"text":"tax_add = NULL","type":"text","marks":[{"type":"code"}]},{"text":" will just give the genus name.","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Generate the plot:","type":"text"}]},{"type":"codeBlock","content":[{"text":"p <- amp_boxplot(ampvisdata, tax_show = 7, tax_add = NULL, tax_aggregate = taxgroup, tax_empty = \"remove\", group_by = group)\np","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"You can change additional properties of the plot here:","type":"text"}]},{"type":"codeBlock","content":[{"text":"p$mapping$fill <- as.name(\".Group\")\np <- p + theme_bw() +\nscale_color_manual(values=c(\"Red\", \"Green\", \"Blue\")) + scale_fill_manual(values=c(\"Red\", \"Green\", \"Blue\")) +\nlabs(fill = group, x = taxgroup) +\ntheme(text = element_text(size = 18), axis.text.x = element_text(size=18), axis.text.y = element_text(size=16)) +\nguides(color = \"none\")\np","type":"text"}]},{"type":"paragraph"},{"type":"paragraph","content":[{"text":"Export as pdf and tiff","type":"text"}]},{"type":"codeBlock","content":[{"text":"tiff_exp <- paste0(\"community_tax_abund_BW_\", group, \"_\", taxgroup, \"_samples_.tiff\")\nggsave(file = tiff_exp, dpi = 300, compression = \"lzw\", device = \"tiff\", plot = p, width = 20, height = 20, units = \"cm\")\npdf_exp <- paste0(\"community_tax_abund_BW_\", group, \"_\", taxgroup, \"_samples_.pdf\")\nggsave(file = pdf_exp, device = \"pdf\", plot = p, width = 20, height = 20, units = \"cm\")","type":"text"}]},{"type":"paragraph"},{"type":"paragraph"}],"version":1}