The solutions given below are just one way of obtaining the desired plots! There are probably several different ways you could code to get the same plots.
For the exercises in the first two lab sessions: R introduction
and the ggplot basics 1
, all you need to do is to replace the filename raw_counts.txt
with each of the different files to look at the differences between different normalization methods.
Task Plot 1:
gc_raw <- read.table(file = "data/counts_raw.txt", sep = "\t", header = T)
gc_filt <- read.table(file = "data/counts_filtered.txt", sep = "\t", header = T)
gc_vst <- read.table(file = "data/counts_vst.txt", sep = "\t", header = T)
gc_deseq <- read.table(file = "data/counts_deseq2.txt", sep = "\t", header = T)
md <- read.table("data/metadata.csv", header = T, sep = ";")
gene_counts_all <-
gc_raw %>% gather(Sample_ID, Raw, -Gene) %>%
full_join(gc_filt %>% gather(Sample_ID, Filtered, -Gene), by = c("Gene", "Sample_ID")) %>%
full_join(gc_vst %>% gather(Sample_ID, VST, -Gene), by = c("Gene", "Sample_ID")) %>%
full_join(gc_deseq %>% gather(Sample_ID, DESeq2, -Gene), by = c("Gene", "Sample_ID")) %>%
gather(Method, count, Raw:DESeq2) %>%
filter(!is.na(count)) %>%
full_join(md, by = "Sample_ID")
gene_counts_all$Time <- factor(gene_counts_all$Time, levels = c("t0","t2","t6","t24"))
gene_counts_all$Replicate <- factor(gene_counts_all$Replicate, levels = c("A","B","C"))
gene_counts_all$Method <- factor(gene_counts_all$Method, levels = c("Raw","Filtered","DESeq2","VST"))
gene_counts_all %>%
group_by(Time, Replicate, Method) %>%
summarise(mean=mean(log10(count +1)),se=se(log10(count +1))) %>%
ggplot(aes(x= Time, y= mean, fill = Replicate)) +
geom_bar(position = position_dodge2(), stat = "identity") +
geom_errorbar(aes(ymin=mean-se, ymax=mean+se), position = position_dodge2(.9, padding = .6)) +
facet_wrap(~Method, scales = "free")
Task Plot 2:
gene_counts_all %>%
group_by(Time, Replicate, Method) %>%
ggplot() +
geom_boxplot(mapping = aes(x = Sample_Name, y = log10(count + 1), fill = Time)) +
facet_wrap(~Method*Replicate, ncol = 3, scales = "free")
gc_long %>%
group_by(Time, Replicate) %>%
summarise(mean=mean(log10(count +1)),se=se(log10(count +1))) %>%
ggplot(aes(x=Time, y=mean, color = Replicate)) +
facet_wrap(~Replicate) +
geom_line(aes(group=1), stat= "identity", size = 2) +
scale_x_discrete(limits= c("t0", "t2", "t24")) +
scale_y_continuous(limits = c(0.4,0.8), breaks = seq(0.4,0.8,0.05)) +
guides(color="none") +
ylab(label = "mean(log10(count + 1))") +
theme_light() +
theme(axis.text = element_text(face="bold", size=12),
axis.title = element_text(face="bold", color = "#C84DF9", size=14),
axis.ticks = element_blank())
p4 <- ggplot(data=iris,mapping=aes(x=Sepal.Length, y = Sepal.Width, color = Species)) +
geom_point(size = 3, alpha = 0.6) +
theme_classic(base_size = 12) +
border()
d1 <- ggplot(data=iris,mapping=aes(Sepal.Length, fill = Species)) +
geom_density(alpha = 0.6) +
theme_classic() +
clean_theme() +
theme(legend.position = "none")
d2 <- ggplot(data=iris,mapping=aes(Sepal.Width, fill = Species)) +
geom_density(alpha = 0.6) +
theme_classic() +
clean_theme() +
theme(legend.position = "none") +
rotate()
ggarrange(d1, NULL, p4, d2,
ncol = 2, nrow = 2, align = "hv",
widths = c(3, 1), heights = c(1, 3),
common.legend = TRUE)
Eigenvalues <- gc_mds$eig
Variance <- Eigenvalues / sum(Eigenvalues)
Variance1 <- 100 * signif(Variance[1], 3)
Variance2 <- 100 * signif(Variance[2], 3)
gc_mds_long <- gc_mds$points %>%
as.data.frame() %>%
rownames_to_column("Sample_ID") %>%
full_join(md, by = "Sample_ID")
gc_mds_long$Sample_Name <- factor(gc_mds_long$Sample_Name, levels = c("t0_A","t0_B","t0_C","t2_A","t2_B","t2_C","t6_A","t6_B","t6_C","t24_A","t24_B","t24_C"))
gc_mds_long$Time <- factor(gc_mds_long$Time, levels = c("t0","t2","t6","t24"))
gc_mds_long$Replicate <- factor(gc_mds_long$Replicate, levels = c("A","B","C"))
ggplot(gc_mds_long, aes(x=V1, y=V2, color = Time)) +
geom_point(size = 3, aes(shape = Replicate)) +
xlab(paste("PCO1: ", Variance1, "%")) +
ylab(paste("PCO2: ", Variance2, "%")) +
geom_vline(xintercept = 0, linetype=2) +
geom_hline(yintercept = 0, linetype=2) +
theme_bw() +
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank())
sessionInfo()
## R version 4.1.3 (2022-03-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.6 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/openblas/libblas.so.3
## LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so
##
## locale:
## [1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C LC_TIME=C.UTF-8
## [4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
## [7] LC_PAPER=C.UTF-8 LC_NAME=C LC_ADDRESS=C
## [10] LC_TELEPHONE=C LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] grid stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] shiny_1.7.1 ggimage_0.3.1 treeio_1.18.1
## [4] ggtree_3.2.1 pheatmap_1.0.12 swemaps_1.0
## [7] mapdata_2.3.0 maps_3.4.0 gridExtra_2.3
## [10] jpeg_0.1-9 ggpubr_0.4.0 cowplot_1.1.1
## [13] ggthemes_4.2.4 scales_1.2.0 ggrepel_0.9.1
## [16] wesanderson_0.3.6 forcats_0.5.1 stringr_1.4.0
## [19] purrr_0.3.4 readr_2.1.2 tidyr_1.2.0
## [22] tibble_3.1.6 tidyverse_1.3.1 reshape2_1.4.4
## [25] ggplot2_3.3.5 formattable_0.2.1 kableExtra_1.3.4
## [28] dplyr_1.0.8 lubridate_1.8.0 leaflet_2.1.1
## [31] yaml_2.3.5 fontawesome_0.2.2.9000 captioner_2.2.3
## [34] bookdown_0.26 knitr_1.39
##
## loaded via a namespace (and not attached):
## [1] colorspace_2.0-3 ggsignif_0.6.3 ellipsis_0.3.2
## [4] fs_1.5.2 aplot_0.1.3 rstudioapi_0.13
## [7] xaringan_0.24 farver_2.1.0 fansi_1.0.3
## [10] xml2_1.3.3 splines_4.1.3 cachem_1.0.6
## [13] jsonlite_1.8.0 broom_0.8.0 dbplyr_2.1.1
## [16] compiler_4.1.3 httr_1.4.2 backports_1.4.1
## [19] assertthat_0.2.1 Matrix_1.4-0 fastmap_1.1.0
## [22] lazyeval_0.2.2 cli_3.3.0 later_1.3.0
## [25] leaflet.providers_1.9.0 htmltools_0.5.2 tools_4.1.3
## [28] gtable_0.3.0 glue_1.6.2 Rcpp_1.0.8.3
## [31] carData_3.0-5 cellranger_1.1.0 jquerylib_0.1.4
## [34] vctrs_0.4.1 ape_5.6-2 svglite_2.1.0
## [37] nlme_3.1-155 crosstalk_1.2.0 xfun_0.30
## [40] ps_1.7.0 rvest_1.0.2 mime_0.12
## [43] lifecycle_1.0.1 rstatix_0.7.0 promises_1.2.0.1
## [46] hms_1.1.1 parallel_4.1.3 RColorBrewer_1.1-3
## [49] curl_4.3.2 ggfun_0.0.6 yulab.utils_0.0.4
## [52] sass_0.4.1 stringi_1.7.6 highr_0.9
## [55] tidytree_0.3.9 rlang_1.0.2 pkgconfig_2.0.3
## [58] systemfonts_1.0.4 evaluate_0.15 lattice_0.20-45
## [61] patchwork_1.1.1 htmlwidgets_1.5.4 labeling_0.4.2
## [64] processx_3.5.3 tidyselect_1.1.2 plyr_1.8.7
## [67] magrittr_2.0.3 R6_2.5.1 magick_2.7.3
## [70] generics_0.1.2 DBI_1.1.2 pillar_1.7.0
## [73] haven_2.5.0 withr_2.5.0 mgcv_1.8-39
## [76] abind_1.4-5 modelr_0.1.8 crayon_1.5.1
## [79] car_3.0-12 utf8_1.2.2 tzdb_0.3.0
## [82] rmarkdown_2.14 readxl_1.4.0 callr_3.7.0
## [85] reprex_2.0.1 digest_0.6.29 webshot_0.5.3
## [88] xtable_1.8-4 httpuv_1.6.5 gridGraphics_0.5-1
## [91] munsell_0.5.0 viridisLite_0.4.0 ggplotify_0.1.0
## [94] bslib_0.3.1