diff --git a/renv.Rmd b/renv.Rmd deleted file mode 100644 index e69de29bb2d1d6434b8b29ae775ad8c2e48c5391..0000000000000000000000000000000000000000 diff --git a/workflows/deseq2/README.md b/workflows/deseq2/README.md new file mode 100644 index 0000000000000000000000000000000000000000..5c1de5c392472eea78490002ea2f55e3116e2176 --- /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 0000000000000000000000000000000000000000..77b80914474c31b2713ed085c34d2ef160cb7660 --- /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 0000000000000000000000000000000000000000..c2adacadd543bc78e8ad34c3287f189e496bed15 --- /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() + +``` + + +