How can I make boxplot for a single gene from one of the contrasts generated by limma?

Question

[bioinfo8]

Hi,

I am interested to make boxplots for top 20 DEGs (20 separate images) from each of the 5 contrasts. My input data is 1-colour Agilent microarray and I have used limma. I could not find any function equivalent to plotCounts from DESeq.

Please guide.

Thanks.

[cpad0112]

1. Identify the gene of interest
2. Extract transformed/normalized expression values for the gene (significant probe for that gene) for each sample.
3. Do a boxplot and facet by group.

In R, you can loop over the gene IDs/symbols to get 20 boxplots in a single page or as 20 images.

[bioinfo8]

I tried for 1 contrast, but I don't know how to proceed further to extract expression values from fit to make boxplot, would you @[cpad0112] please help:

# Identified top 20 genes.
library(limma)
......................
fit <- eBayes(fit)
df1 <- topTable(fit, coef=1, adjust.method="BH", number=20)
genes_20 <-  as.character(df1$geneIDs)           

# Extract transformed/normalized expression values for the gene 
# (significant probe for that gene) for each sample.

# Do a boxplot and facet by treatment.

Answer 1

[Gordon Smyth]

It is so easy to plot data using R that I don't see any need for limma to provide wrapper functions for everything. For example, if y is an EList object in limma, Group is a factor, and you want to plot gene i vs Group, then just

boxplot(y$E[i,] ~ Group)

If y is a different type of object, then just substitute exprs(y)[i,] where exprs(y) is the matrix of log-expression values.

You seem to be assuming that you can plot data from the linear model fit, but you can't. The expresssion data is stored in the data object not in the fit object.

Your question might be partly about how to set i, but you would need to give more information about your code pipeline for people to help you with that.

The reference to plotCounts is a bit confusing since it doesn't make boxplots and doesn't work with contrasts.