From bf30e495bba3f460f0df29a6b60c464a5e262b1d Mon Sep 17 00:00:00 2001
From: Dominik Brilhaus <brilhaus@nfdi4plants.org>
Date: Mon, 4 Nov 2024 16:41:45 +0100
Subject: [PATCH] add deseq2 analysis

---
 renv.Rmd                          |  0
 workflows/deseq2/README.md        | 11 ++++
 workflows/deseq2/dependencies.Rmd | 24 ++++++++
 workflows/deseq2/deseq2.Rmd       | 98 +++++++++++++++++++++++++++++++
 4 files changed, 133 insertions(+)
 delete mode 100644 renv.Rmd
 create mode 100644 workflows/deseq2/README.md
 create mode 100644 workflows/deseq2/dependencies.Rmd
 create mode 100644 workflows/deseq2/deseq2.Rmd

diff --git a/renv.Rmd b/renv.Rmd
deleted file mode 100644
index e69de29..0000000
diff --git a/workflows/deseq2/README.md b/workflows/deseq2/README.md
new file mode 100644
index 0000000..5c1de5c
--- /dev/null
+++ b/workflows/deseq2/README.md
@@ -0,0 +1,11 @@
+# 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
diff --git a/workflows/deseq2/dependencies.Rmd b/workflows/deseq2/dependencies.Rmd
new file mode 100644
index 0000000..77b8091
--- /dev/null
+++ b/workflows/deseq2/dependencies.Rmd
@@ -0,0 +1,24 @@
+---
+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")
+
+```
+
diff --git a/workflows/deseq2/deseq2.Rmd b/workflows/deseq2/deseq2.Rmd
new file mode 100644
index 0000000..c2adaca
--- /dev/null
+++ b/workflows/deseq2/deseq2.Rmd
@@ -0,0 +1,98 @@
+---
+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()
+  
+```
+
+
+
-- 
GitLab