DEseq2 usage

DEseq2 installation

source("https://bioconductor.org/biocLite.R")
biocLite("DESeq2")



wd="D:/gene_expression/"
setwd(wd)
# parse gene expression from file
count_table<-read.table(file="gene_expressed_all_samples.csv",header=TRUE)
#count_table like:
GeneID control1 control2 treatment1 treatment2
Gene1                   1              2                 4                   6
Gene2                   3              8                 5                   6
Gene3                   1              2                 8                   9

# Create matrix for Count table of DEseq
cts <- count_table[,c("control1 ","control2 ","treatment1 ","treatment2")]
cts <- as.matrix(cts)
row.names(cts)<-count_none[,1]
cts<-cts[rowSums(cts)>0,]

# Create library's information(condition/design) for DEseq
coldata <-data.frame(
row.names = c("control1","control2","treatment1","treatment2"), condition = c("control","control","treatment","treatment"))

library(DESeq2)
dds <- DESeqDataSetFromMatrix(countData = cts,
                              colData = coldata,
                              design = ~ condition)
featureData <- data.frame(gene=rownames(cts))
mcols(dds) <- DataFrame(mcols(dds), featureData)
dds$condition <- relevel(dds$condition, ref="control")
dds <- DESeq(dds)
res <- results(dds, alpha=0.05)
DEres <- as.matrix(res)
count.nor <- as.matrix( counts(dds,normalized=TRUE) )
resMerge <- cbind(DEres,count.nor)
resMerge<-resMerge[order(resMerge[,colnames(resMerge)=="padj"],decreasing = F),]
write.csv(resMerge, file=paste0(wd,"condition_treated_results_count.csv",collapse=TRUE))