Using the following code to look at differential expression in a couple of datasets. Over and over again the # of up and and down regulated genes is approximately equal. I just wonder if I'm missing something, but I can't find anything about this anywhere. It can't be an artifact of normalization. Any ideas?

library("DESeq2")

countData <- as.matrix(read.csv("genecountmatrix.csv", row.names="geneid")) colData <- read.csv("phenodata.txt", sep="\t", row.names=1)

all(rownames(colData) %in% colnames(countData))

countData <- countData[, rownames(colData)] all(rownames(colData) == colnames(countData))

dds <- DESeqDataSetFromMatrix(countData = countData, colData = colData, design = ~ Treatment)

dds <- DESeq(dds)

res <- results(dds)

(resOrdered <- res[order(res$padj), ])

It's not necessarily an artifact of normalization, although the normalization does seek to minimize the fold changes across all samples through a single scaling factor per sample to account for sequencing depth.

E.g. below I spike in 50 genes with a positive LFC. In the output I have 29:1 ratio of positive to negative LFC at FDR 1%. So I demonstrate it does not have to be balanced 1:1.

However, normalization is a concern as the number of DE genes increases. If you have an experiment with very lopsided effects (e.g. a vast number of genes are up-regulated), you should not rely solely on in silico normalization, that would be a bad experimental design.

> dds <- makeExampleDESeqDataSet(n=5000, dispMeanRel=function(x) 4/x + 0.01)
> counts(dds)[1:50,7:12] <- counts(dds)[1:50,7:12] * 3L
> dds <- DESeq(dds, quiet=TRUE)
> res <- results(dds)
> summary(res, alpha=.01)

out of 4990 with nonzero total read count
adjusted p-value < 0.01
LFC > 0 (up)       : 29, 0.58%
LFC < 0 (down)     : 1, 0.02%