add deseq2 analysis

workflows/deseq2/README.md

+# DESeq2
+
+Workflow used for **differential gene expression analysis**
+
+## DESeq2 docs
+
+- https://bioconductor.org/packages/release/bioc/html/DESeq2.html
+
+## Importing kallisto output with tximport
+
+- https://bioconductor.org/packages/release/bioc/vignettes/tximport/inst/doc/tximport.html#kallisto

workflows/deseq2/dependencies.Rmd

+---
+title: "Install dependencies"
+author: "Dominik Brilhaus"
+date: "`r Sys.Date()`"
+output: html_document
+---
+
+
+```{r}
+if (!require("BiocManager", quietly = TRUE))
+    install.packages("BiocManager")
+
+
+BiocManager::install("DESeq2")
+library("DESeq2")
+
+BiocManager::install("tximport")
+library("tximport")
+
+BiocManager::install("rhdf5")
+library("rhdf5")
+
+```
+

workflows/deseq2/deseq2.Rmd

+---
+title: "deseq2"
+author: "Dominik Brilhaus"
+date: "`r Sys.Date()`"
+output: html_document
+---
+
+
+## Libraries
+
+```{r}
+
+library("DESeq2")
+library("tximport")
+library("rhdf5")
+library("ggplot2")
+
+```
+
+
+
+
+## In-and-out
+
+```{r}
+
+arc <- "../../"
+
+inKallistoResults <- file.path(arc, "runs/kallisto/kallisto_results")
+inMetadataFile <- file.path(arc, "runs/merged_isa_metadata/out/merged_isa.tsv")
+inMetadataSample <- "Source.Name"
+inMetadataFactor <- "Factor..Photosynthesis.mode."
+
+
+```
+
+## Import kallisto count data
+
+```{r}
+
+files <- dir(inKallistoResults, recursive = T, full.names = T ,"abundance.h5")
+names(files) <- dir(inKallistoResults)
+
+txi <- tximport(files, type = "kallisto", txOut = TRUE)
+
+head(txi$counts)
+
+```
+
+
+## Read sample metadata
+
+```{r}
+
+
+samples_metadata <- read.table(file = inMetadataFile, sep = "\t")
+
+samples <- samples_metadata[order(samples_metadata[[inMetadataSample]]), c(inMetadataSample, inMetadataFactor)]
+colnames(samples)[1:2] <- c("sampleID", "condition")
+
+rownames(samples) <- samples$sampleID
+
+
+```
+
+## DESeq
+
+```{r}
+dds <- DESeqDataSetFromTximport(txi,
+                                   colData = samples,
+                                   design = ~ condition)
+
+dds <- DESeq(dds)
+
+res <- results(dds)
+res
+
+plotMA(res, ylim=c(-2,2))
+
+
+```
+
+
+```{r}
+vsd <- vst(dds, blind=FALSE)
+
+pcaData <- plotPCA(vsd, intgroup=c("condition"), returnData=TRUE)
+percentVar <- round(100 * attr(pcaData, "percentVar"))
+ggplot(pcaData, aes(PC1, PC2, color=condition)) +
+  geom_point(size=3) +
+  xlab(paste0("PC1: ",percentVar[1],"% variance")) +
+  ylab(paste0("PC2: ",percentVar[2],"% variance")) + 
+  coord_fixed()
+  
+```
+
+
+