Hi,

I am analyzing mRNA-Seq dataset using EdgeR package and testing filtering by filterByExpr and rowSums that would keep genes. I have question about interpreting Histogram of average log2 CPM in EdgeR?

I tested filtering in 4 different ways, and would like to know how to interpret the plot? Basically, filterByExpr looks good, however, I am interested in creating model.matrix on other variables too like treatment, severity, etc., for comparisons. How to do I decide the cut-off in perhaps rowsums? What does the negative values in the x-axis signifies? Should the graph look like bell shaped distribution?

Thank you in advance.

Best Regards, Toufiq

dge <- DGEList(counts = Counts,  remove.zeros = TRUE)
dge$samples

# Either;
## filterByExpr
keep <- filterByExpr(dge, design). ## ## Pairing and blocking is essential for comparison as different cells are extracted from same subjects
table(keep.keep)

## (OR)
## Filtering to remove low counts
keep <- rowSums(dge$counts) >= 10

## (OR)
## Filtering to remove low counts
<- rowSums(dge$counts) >= 50

dge <- dge[keep, , keep.lib.sizes=FALSE]
dge$counts
dim(dge$counts)

AveLogCPM <- aveLogCPM(dge)
hist(AveLogCPM)

Thank you,

Toufiq

In the edgeR workflow vignette we discuss three possible filtering methods:

1. filterByExpr(y, design)
2. rowSums(Counts) > 50
3. aveLogCPM(y)

The use of filterByExpr with default arguments is the recommended method, and I don't see a reason not to do that. The exact covariates in the design matrix makes little difference to the filtering. You should just include the main covariates and treatment factors. If you have paired design, then include the treatments and covariates in the design for filtering but leave out the subjects.

If you filter by total gene counts, then a value like 50 needs to be used. Using >= 10 is too small.

The third method is to examine the histogram and filter by AveLogCPM. For your data, you would note that the histogram is bimodal and a value like -1 or -2 would separate the modes. Therefore you could filter by

keep <- (AveLogCPM > -2)

in order to remove the mode of small values on the left.

edgeR makes no assumptions about the distribution of the AveLogCPM values, so the shape of the histogram is irrelevant. The only requirement is that there shouldn't be a left tail of very small values.

Negative AveLogCPM values simply mean that the average CPM is less than 1. There's nothing unusual or worrying about that.