@@ -28,20 +28,24 @@ knitr::opts_chunk$set(eval = TRUE,
2828 message = FALSE,
2929 engine = "R", # Chunks will always have R code, unless noted
3030 error = TRUE,
31+ fig.height = 7, fig.width = 7,
3132 fig.path="./Figures/") # Set the figure options
3233```
3334
3435
3536
3637# Load Packages & Functions
37- ``` {r setup}
38+ ``` {r load-packages-funs}
39+ # Load functions and packages that are necessary
3840source("functions.R")
3941library(plotly)
4042
4143# We'd like to plot with pretty colors based on national park posters :)
4244# install.packages("devtools")
4345#devtools::install_github("katiejolly/nationalparkcolors")
46+
4447library(nationalparkcolors)
48+ # Assign the 4 colors for the 4 different depths
4549colors4 <- park_palette("Saguaro", 4)
4650```
4751
@@ -101,7 +105,6 @@ t_mat <- t(as.matrix(mag_table))
101105norm_mat <- t_mat/ordered_checkm$exp_genome_size
102106t_norm_mat <- t(norm_mat)
103107
104-
105108# Combine into a normalized phyloseq object
106109tara_norm_physeq <- phyloseq(meta_physeq, tax_physeq, mag_tree,
107110 otu_table(t_norm_mat, taxa_are_rows = TRUE))
@@ -111,59 +114,80 @@ tara_norm_physeq <- phyloseq(meta_physeq, tax_physeq, mag_tree,
111114## Heatmap
112115
113116``` {r heatmap, fig.height=8, fig.width=10}
117+ # A clustered heatmap
114118heatmap(otu_table(tara_norm_physeq))
115119
116120# Visualize only the top 50 taxa
117121top_50MAGs <- names(sort(taxa_sums(tara_norm_physeq), decreasing = TRUE))[1:50]
118122top_50MAGs_physeq <- prune_taxa(top_50MAGs, tara_norm_physeq)
119123
120124# Subset only the bacterial samples
121- bac_top50MAGs <- subset_samples(top_50MAGs_physeq, size_fraction =="bact ")
125+ girus_top50MAGs <- subset_samples(top_50MAGs_physeq, size_fraction =="girus ")
122126
123- otu_long_50_bact <- psmelt(otu_table(bac_top50MAGs)) %>%
124- mutate(log_abund = log2(Abundance+0.0000001))
127+ # Melt into long format and fix zeros to avoid infinity values
128+ otu_long_50_girus <- psmelt(otu_table(girus_top50MAGs)) %>%
129+ mutate(log_abund = log2(Abundance + 0.0000001))
125130
126- heat_plot <- otu_long_50_bact %>%
131+ # Make a "heatmap" with geom_tile that works with plotly :)
132+ heat_plot <- otu_long_50_girus %>%
127133 ggplot(aes(x=Sample, y=OTU)) +
128134 geom_tile(aes(fill = log_abund)) +
129135 scale_fill_distiller(palette = "YlGnBu") +
130136 labs(title = "Top 50 MAGs") +
131137 theme(axis.text.x = element_text(angle = 30, vjust = 1, hjust = 1))
138+
139+ # Plot the plotly plot!
132140ggplotly(heat_plot)
133141```
134142
135143## Beta diversity plots
136- ``` {r betadiv}
144+ ``` {r betadiv, fig.width = 9}
145+ # Calculate the bray curtis distances for all samples
137146tara_norm.ord <- ordinate(tara_norm_physeq, method = "PCoA", distance = "bray")
138- plot_ordination(tara_norm_physeq, tara_norm.ord, color = "province") +
147+
148+ # Make a faceted plot by size fraction
149+ sizeFrac_PCoA_prov <-
150+ plot_ordination(tara_norm_physeq, tara_norm.ord, color = "province") +
139151 facet_wrap(~size_fraction) +
140152 theme_minimal() +
141153 scale_color_brewer(palette = "Paired")
142154
155+ # Make it a plotly plot
156+ ggplotly(sizeFrac_PCoA_prov)
157+ ```
158+
159+ ``` {r betadiv2, fig.height = 4, fig.width=5}
143160# Subset the bacterial samples and color by depth
144161bact <- subset_samples(tara_norm_physeq, size_fraction == "bact")
162+ # Calculate the bray curtis distances for a PCoA
145163bact.ord <- ordinate(bact, method = "PCoA", distance = "bray")
146- plot_ordination(bact, bact.ord, color = "depth")+
164+ # Make the plot
165+ bact_PCoA_depth <- plot_ordination(bact, bact.ord, color = "depth")+
147166 theme_minimal() +
148167 geom_point(size = 2) +
149168 scale_color_manual(values = colors4)
150169
170+ # Make it a plotly plot
171+ ggplotly(bact_PCoA_depth)
151172```
152173
153174## Abundance plots
154175``` {r abundance plots}
155- top_taxa <- names(sort(taxa_sums(bact),decreasing = TRUE))[1:4]
176+ # Select the top 4 taxa
177+ top_taxa <- names(sort(taxa_sums(bact),decreasing = TRUE))[1:4]
156178top_bact <- prune_taxa(top_taxa, bact)
179+ # Melt the data to the long format
157180top_bact_df <- psmelt(top_bact)
158181
159- abundplot <- ggplot(top_bact_df, aes(x=province, y= Abundance, color = depth))+
182+ # Plot it! :)
183+ abundplot <- ggplot(top_bact_df, aes(x=province, y= Abundance, color = depth)) +
160184 facet_wrap(~OTU, scales = "free_y") +
161185 geom_jitter(size = 2) + theme_minimal()+
162186 labs(y = "Normalized MAG Abundance") +
163187 theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 1),
164188 axis.title.x = element_blank()) +
165189 scale_color_manual(values = colors4)
166190
191+ # Make it interactive
167192ggplotly(abundplot)
168-
169193```
0 commit comments