Skip to content

Commit 977f763

Browse files
authored
Merge pull request #14 from biosustain/zuvale_branch
🎨 minor changes to downstream analysis-tutorial R-markdowns
2 parents fe975d5 + ece2057 commit 977f763

5 files changed

Lines changed: 65 additions & 43 deletions

tutorial/01_Alpha-diversity_tutorial.Rmd

Lines changed: 16 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -112,19 +112,22 @@ library(FSA)
112112
**Load species abundance, taxonomical annotation and metadata file**
113113

114114
```{r, warning=FALSE}
115-
abundance <- readRDS(file = "../data/MetaphlanAbundance_Species.rds")
115+
TRAINING_DIR <- "/workspaces/dsp_metagenomics_training/"
116+
```
117+
118+
```{r, warning=FALSE}
119+
abundance <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAbundance_Species.rds"))
116120
rownames(abundance) <- gsub("_SRR_db1.metaphlan", "", rownames(abundance))
117121
118122
# Taxonomical annotation
119-
annotation <- readRDS(file = "../data/MetaphlanAnnotations_Species.rds")
123+
annotation <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAnnotations_Species.rds"))
120124
121125
# Metadata
122-
metadata <- read.table(file = "../data/metadata.tsv", header = TRUE,
126+
metadata <- read.table(file = paste0(TRAINING_DIR, "data/metadata.tsv"), header = TRUE,
123127
sep = "\t",
124128
quote = "",
125129
row.names = NULL)
126130
rownames(metadata) <- metadata$Sample
127-
128131
```
129132

130133

@@ -147,7 +150,7 @@ abundance <- abundance[rownames(metadata), ]
147150
## Creating a directory for the related results
148151

149152
```{r, warning=FALSE}
150-
results_dir <- "../results/report/01_Alpha-diversity"
153+
results_dir <- paste0(TRAINING_DIR, "results/report/01_Alpha-diversity/")
151154
152155
dir.create(results_dir,
153156
recursive = TRUE,
@@ -196,7 +199,7 @@ colnames(plot_data) <-
196199
c("Richness", "Shannon", "Area", "Sampling_Date", "Area_Sampling")
197200
198201
write.csv(plot_data,
199-
file = "../results/report/01_Alpha-diversity/01_Alpha_diversity.csv",
202+
file = paste0(results_dir, "02_Alpha_diversity.csv"),
200203
row.names = TRUE,
201204
quote = FALSE)
202205
```
@@ -214,7 +217,8 @@ richness_boxplot <- ggplot(plot_data,
214217
y = Richness)) +
215218
216219
geom_boxplot(aes(fill = Area_Sampling),
217-
outliers = FALSE,
220+
#outliers = FALSE,
221+
outlier.shape = NA,
218222
alpha = 0.7) +
219223
220224
geom_point(aes(fill = Area_Sampling),
@@ -245,7 +249,8 @@ shannon_boxplot <- ggplot(plot_data,
245249
y = Shannon)) +
246250
247251
geom_boxplot(aes(fill = Area_Sampling),
248-
outliers = FALSE,
252+
#outliers = FALSE,
253+
outlier.shape = NA,
249254
alpha = 0.7) +
250255
251256
geom_point(aes(fill = Area_Sampling),
@@ -272,7 +277,7 @@ richness_boxplot +
272277
### Save
273278

274279
```{r, warning=FALSE}
275-
ggsave(filename = "../results/report/01_Alpha-diversity/02_alphadiv_boxplots.png",
280+
ggsave(filename = paste0(results_dir, "01_alphadiv_boxplots.png"),
276281
width = 8,
277282
height = 4,
278283
dpi = 300)
@@ -340,7 +345,7 @@ kruskal.test(Shannon ~ Area,
340345

341346
**Almost significant!**
342347

343-
When it is not significant, we don't needed to test between groups, but lets double check:
348+
When it is not significant, we don't need to test between groups, but lets double check:
344349

345350
```{r, warning=FALSE}
346351
dunnTest(Shannon ~ Area,
@@ -399,4 +404,4 @@ sessionInfo()
399404
rmarkdown::render("01_Alpha-diversity_tutorial.Rmd")
400405
# then run this to convert HTML to PDF (if needed)
401406
#pagedown::chrome_print("core-report.html",output="core-report.pdf")
402-
```
407+
```

tutorial/02_Beta-diversity_tutorial.Rmd

Lines changed: 14 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -116,17 +116,21 @@ library(vegan)
116116

117117
**Load species abundance, taxonomical annotation and metadata file**
118118

119+
```{r}
120+
TRAINING_DIR <- "/workspaces/dsp_metagenomics_training/"
121+
```
122+
119123
```{r}
120124
# Abundance table
121-
abundance <- readRDS(file = "../data/MetaphlanAbundance_Species.rds")
125+
abundance <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAbundance_Species.rds"))
122126
rownames(abundance) <- gsub("_SRR_db1.metaphlan", "", rownames(abundance))
123127
124128
# Taxonomical annotation
125-
annotation <- readRDS(file = "../data/MetaphlanAnnotations_Species.rds")
129+
annotation <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAnnotations_Species.rds"))
126130
127131
# Metadata
128-
metadata <- read.table(file = "../data/metadata.tsv",
129-
header = TRUE, sep = "\t",
132+
metadata <- read.table(file = paste0(TRAINING_DIR, "data/metadata.tsv"), header = TRUE,
133+
sep = "\t",
130134
quote = "",
131135
row.names = NULL)
132136
rownames(metadata) <- metadata$Sample
@@ -137,7 +141,7 @@ abundance <- abundance[rownames(metadata), ]
137141

138142
# Creating a directory for the PCoA related results
139143
```{r}
140-
results_dir <- "../results/report/02_Beta-diversity"
144+
results_dir <- paste0(TRAINING_DIR, "results/report/02_Beta-diversity/")
141145
142146
dir.create(results_dir, recursive = TRUE, showWarnings = FALSE)
143147
```
@@ -265,7 +269,7 @@ p <- ggplot(scores,
265269
p
266270
267271
# Save it as a png
268-
ggsave("../results/report/02_Beta-diversity/01_PCoA_plot.png",
272+
ggsave(paste0(results_dir, "01_PCoA_plot.png"),
269273
width = 8, height = 6, dpi = 300)
270274
```
271275

@@ -320,7 +324,7 @@ p <- ggplot(scores,
320324
p
321325
322326
# Save it as a png
323-
ggsave("../results/report/02_Beta-diversity/02_PCoA_with_loadings.png",
327+
ggsave(paste0(results_dir, "02_PCoA_with_loadings.png"),
324328
width = 8, height = 6, dpi = 300)
325329
```
326330

@@ -395,7 +399,7 @@ dispersion_box <- ggplot(dispersion_plot_data, aes(x = group,
395399
dispersion_box
396400
397401
# Save it as a png
398-
ggsave("../results/report/02_Beta-diversity/03_dispersion_plot.png",
402+
ggsave(paste0(results_dir, "03_dispersion_plot.png"),
399403
width = 6, height = 4, dpi = 300)
400404
```
401405

@@ -448,7 +452,7 @@ dispersion_box <- ggplot(dispersion_plot_data,
448452
dispersion_box
449453
450454
# Save it as a png
451-
ggsave("../results/report/02_Beta-diversity/04_dispersion_plot_sampling.png",
455+
ggsave(paste0(results_dir, "04_dispersion_plot_sampling.png"),
452456
width = 5, height = 4, dpi = 300)
453457
```
454458

@@ -481,7 +485,7 @@ adonis_table <- adonis2(form,
481485
by = "margin")
482486
483487
write.csv(adonis_table,
484-
file = "../results/report/02_Beta-diversity/05_Adonis_table.csv",
488+
file = paste0(results_dir, "05_Adonis_table.csv"),
485489
row.names = TRUE, quote = FALSE)
486490
```
487491

tutorial/03_Abundance-overview_tutorial.Rmd

Lines changed: 12 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -71,17 +71,21 @@ library(RColorBrewer)
7171

7272
**Load species abundance, taxonomical annotation and metadata file**
7373

74+
```{r}
75+
TRAINING_DIR <- "/workspaces/dsp_metagenomics_training/"
76+
```
77+
7478
```{r}
7579
# Abundance table
76-
abundance <- readRDS(file = "../data/MetaphlanAbundance_Species.rds")
80+
abundance <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAbundance_Species.rds"))
7781
rownames(abundance) <- gsub("_SRR_db1.metaphlan", "", rownames(abundance))
7882
7983
# Taxonomical annotation
80-
annotation <- readRDS(file = "../data/MetaphlanAnnotations_Species.rds")
84+
annotation <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAnnotations_Species.rds"))
8185
8286
# Metadata
83-
metadata <- read.table(file = "../data/metadata.tsv",
84-
header = TRUE, sep = "\t",
87+
metadata <- read.table(file = paste0(TRAINING_DIR, "data/metadata.tsv"), header = TRUE,
88+
sep = "\t",
8589
quote = "",
8690
row.names = NULL)
8791
rownames(metadata) <- metadata$Sample
@@ -99,7 +103,7 @@ abundance <- abundance[rownames(abundance) %in%
99103
## Creating a directory for the related results
100104

101105
```{r}
102-
results_dir <- "../results/report/03_Abundance-overview"
106+
results_dir <- paste0(TRAINING_DIR, "results/report/03_Abundance-overview/")
103107
104108
dir.create(results_dir,
105109
recursive = TRUE,
@@ -195,7 +199,7 @@ pbar <- ggplot(abundance_long, aes(x = Sample,
195199
pbar
196200
197201
# Save it as a png
198-
ggsave("../results/report/03_Abundance-overview/01_top10_species.png",
202+
ggsave(paste0(results_dir, "01_top10_species.png"),
199203
width = 8, height = 6, dpi = 300)
200204
```
201205

@@ -253,7 +257,7 @@ pbar_sum <- ggplot(abundance_summary,
253257
pbar_sum
254258
255259
# Save it as a png
256-
ggsave("../results/report/03_Abundance-overview/02_top10_species_per_area_sampling.png",
260+
ggsave(paste0(results_dir, "02_top10_species_per_area_sampling.png"),
257261
width = 8, height = 6, dpi = 300)
258262
259263
```
@@ -278,4 +282,4 @@ sessionInfo()
278282
rmarkdown::render("03_Abundance-overview_tutorial.Rmd")
279283
# then run this to convert HTML to PDF (if needed)
280284
#pagedown::chrome_print("core-report.html",output="core-report.pdf")
281-
```
285+
```

tutorial/04_Group-comparison_tutorial.Rmd

Lines changed: 12 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -73,17 +73,21 @@ library(FSA)
7373

7474
**Load species abundance, taxonomical annotation and metadata file**
7575

76+
```{r}
77+
TRAINING_DIR <- "/workspaces/dsp_metagenomics_training/"
78+
```
79+
7680
```{r}
7781
# Abundance table
78-
abundance <- readRDS(file = "../data/MetaphlanAbundance_Species.rds")
82+
abundance <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAbundance_Species.rds"))
7983
rownames(abundance) <- gsub("_SRR_db1.metaphlan", "", rownames(abundance))
8084
8185
# Taxonomical annotation
82-
annotation <- readRDS(file = "../data/MetaphlanAnnotations_Species.rds")
86+
annotation <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAnnotations_Species.rds"))
8387
8488
# Metadata
85-
metadata <- read.table(file = "../data/metadata.tsv",
86-
header = TRUE, sep = "\t",
89+
metadata <- read.table(file = paste0(TRAINING_DIR, "data/metadata.tsv"), header = TRUE,
90+
sep = "\t",
8791
quote = "",
8892
row.names = NULL)
8993
rownames(metadata) <- metadata$Sample
@@ -95,7 +99,7 @@ abundance <- abundance[rownames(metadata),]
9599
## Creating a directory for the related results
96100

97101
```{r}
98-
results_dir <- "../results/report/04_Group-comparison"
102+
results_dir <- paste0(TRAINING_DIR, "results/report/04_Group-comparison")
99103
100104
dir.create(results_dir, recursive = TRUE, showWarnings = FALSE)
101105
```
@@ -107,7 +111,7 @@ As there are three areas, we will use a non-parametric test
107111
(not assuming any distribution) for three groups: Kruskal-Wallis test.
108112
We need to run this test for all species comparing among Areas (1, 2, 3),
109113
therefore we need to code a loop to run the test for all species and
110-
later we need to correct for multiple testing (Bonferroni-Hochberg).
114+
later we need to correct for multiple testing (Benjamini-Hochberg).
111115
For the significant ones we will run a pot-hoc dunnTest to figure out
112116
between which groups the differences exist.
113117
:::
@@ -186,7 +190,7 @@ rownames(kruskal_df) <- kruskal_df$Taxon
186190

187191
## Multiple testing correction
188192

189-
Adjusting for multiple testing by Bonferroni and Hochberg
193+
Adjusting for multiple testing by Benjamini and Hochberg
190194
method (False Discovery rate)
191195

192196

@@ -371,4 +375,4 @@ sessionInfo()
371375
rmarkdown::render("04_Group-comparison_tutorial.Rmd")
372376
# then run this to convert HTML to PDF (if needed)
373377
#pagedown::chrome_print("core-report.html",output="core-report.pdf")
374-
```
378+
```

tutorial/05_Microb-networks_tutorial.Rmd

Lines changed: 11 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -74,19 +74,24 @@ library(dplyr)
7474
```
7575

7676
## Load data
77+
78+
```{r}
79+
TRAINING_DIR <- "/workspaces/dsp_metagenomics_training/"
80+
```
81+
7782
```{r, warning=FALSE}
7883
# Abundance table
79-
abundance <- readRDS(file = "../data/MetaphlanAbundance_Species.rds")
84+
abundance <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAbundance_Species.rds"))
8085
rownames(abundance) <- gsub("_SRR_db1.metaphlan", "", rownames(abundance))
8186
dim(abundance)
8287
8388
# Taxonomical annotation
84-
annotation <- readRDS(file = "../data/MetaphlanAnnotations_Species.rds")
89+
annotation <- readRDS(file = paste0(TRAINING_DIR, "data/MetaphlanAnnotations_Species.rds"))
8590
dim(annotation)
8691
8792
# Metadata
88-
metadata <- read.table(file = "../data/metadata.tsv",
89-
header = TRUE, sep = "\t",
93+
metadata <- read.table(file = paste0(TRAINING_DIR, "data/metadata.tsv"), header = TRUE,
94+
sep = "\t",
9095
quote = "",
9196
row.names = NULL)
9297
rownames(metadata) <- metadata$Sample
@@ -100,7 +105,7 @@ dim(abundance)
100105
## Create a directory to store the results
101106

102107
```{r, warning=FALSE}
103-
results_dir <- "../results/report/05_Microbial-association-networks"
108+
results_dir <- paste0(TRAINING_DIR, "results/report/05_Microbial-association-networks/")
104109
105110
dir.create(results_dir, recursive = TRUE, showWarnings = FALSE)
106111
```
@@ -375,4 +380,4 @@ to demonstrate how to use this package using the same dataset from this course.
375380
rmarkdown::render("05_Microbial-association-networks.Rmd")
376381
# then run this to convert HTML to PDF (if needed)
377382
#pagedown::chrome_print("core-report.html",output="core-report.pdf")
378-
```
383+
```

0 commit comments

Comments
 (0)