Question: Plot normalized expression or count values for a specific gene in limma
1
3.5 years ago by
Yury Bukhman20
Yury Bukhman20 wrote:

Hi,

I am using limma to analyze a barcode sequencing dataset. The raw data are counts, similar to RNA-seq. I am using edgeR::calcNormFactors followed by voom to normalize and transform the count data. I am interested in plotting normalized data for a specific single gene. I have looked around and found plotCounts function from DESeq2. Is there a similar function in limma? If not, is there another package that provides such a function for EList or ExpressionSet objects?

Thanks.

Yury

limma • 641 views
modified 3.5 years ago by Aaron Lun25k • written 3.5 years ago by Yury Bukhman20
Answer: Plot normalized expression or count values for a specific gene in limma
3
3.5 years ago by
Aaron Lun25k
Cambridge, United Kingdom
Aaron Lun25k wrote:

Well, assuming your DGEList or count matrix is y, grouping is a factor indicating which group each library belongs to, and gene is an integer/string specifying the gene of interest, you could do:

norm <- edgeR::cpm(y) # to get normalized expression values.
plot(jitter(as.integer(grouping)), norm[gene,], ylab="CPM", xaxt="n",
xlab="", xlim=c(0.5, nlevels(grouping)+0.5))
axis(side=1, at=seq_len(nlevels(grouping)), label=levels(grouping))


Admittedly, it's not a nice easy function, but you can customize this to your heart's content. I usually end up writing all my own plotting code anyway to get publication-quality figures I'm satisfied with.

ADD COMMENTlink modified 3.5 years ago • written 3.5 years ago by Aaron Lun25k

Thanks, Aaron! I can certainly use this. If I want to visualize the data that are actually used by limma's lmFit function, do I extract values from the E slot of my EList object and multiply them by the edgeR scaling factors? Is this the same as the norm values in your script? I suppose I could also use voom-computed weights to make symbols of different sizes...

1

The expression values in the E slot already incorporate edgeR's normalization factors (assuming you got them from voom) so no extra work is required. Then you can just replace norm with blah$E in the code above, assuming blah is the output from voom. If you want to use weight-based sizes, you could add something like cex=blah$weights*C to the plot call (where C is some constant to scale everything to a nice size).

ADD REPLYlink modified 3.5 years ago • written 3.5 years ago by Aaron Lun25k