Plot normalized expression or count values for a specific gene in limma
1
1
Entering edit mode
Yury Bukhman ▴ 20
@yury-bukhman-3186
Last seen 5.1 years ago

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 • 1.1k views
3
Entering edit mode
Aaron Lun ★ 27k
@alun
Last seen 3 hours ago
The city by the bay

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.

0
Entering edit mode

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
Entering edit mode

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).

0
Entering edit mode

I am trying to use the same plot, but norm here is a count table where rows are gene name and column is a sample, how "grouping" factor identify each sample here?

Traffic: 182 users visited in the last hour
FAQ
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.