4
$\begingroup$

I am trying to make MA-plot for my bulk-RNA-seq dataset and experiencing issues with it. I just copy-pasted the code that can be found here:

https://f1000research.com/articles/4-1070/v2

which is:

library("DESeq2")
plotMA(resLFC1, ylim=c(-5,5))
topGene <- rownames(resLFC1)[which.min(resLFC1$padj)]
with(resLFC1[topGene, ], {
  points(baseMean, log2FoldChange, col="dodgerblue", cex=2, lwd=2)
  text(baseMean, log2FoldChange, topGene, pos=2, col="dodgerblue")
})

However, I am not getting the right plot with the top gene circled, but the following one:

enter image description here

Maybe it is a bug in the plotMA function, not sure. How can I reproduce the figure from the publication?

$\endgroup$
  • $\begingroup$ Are you sure that resLFC1 is created as stated in the article ? There seems to be a lot more points than in the article. Also, do you get any error or message? $\endgroup$ – llrs Apr 13 '18 at 10:23
  • $\begingroup$ Are you sure the weird double quotes (“ versus ") are not the problem? $\endgroup$ – benn Apr 13 '18 at 15:22
  • $\begingroup$ I corrected them, of course. Rstudio complains about it right away. $\endgroup$ – Nikita Vlasenko Apr 13 '18 at 18:39
  • 1
    $\begingroup$ Your plot shows all genes as DE (red points). That’s suspicious. $\endgroup$ – Konrad Rudolph Apr 14 '18 at 17:44
4
$\begingroup$

It's likely that the fold-change of the topGene is greater than 5. Consequently, the circle is being drawn and the text produced, but they're outside the bounds of the plot (the bounds don't update if you add a circle and text label). Look at the log2FoldChange column of resLFC1[topGene, ] and increase the bounds to encompass that.

| improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.