Courtesy Friederike jared.andrews07

NB - please post questions on Bioconductor Support: https://support.bioconductor.org/

[ direct link to create a post: https://support.bioconductor.org/p/new/post/?tag_val=EnhancedVolcano ]

Yep, there is a way to use EnhancedVolcano to also produce an MA plot.

MA plots plot fold-change versus mean expression (as opposed to fold-change versus -log10 p-values in a Volcano). MA plots have been used since the days of the cDNA array and, while array data was typically already log2-transformed as part of normalisation, RNA-seq data is count-based; thus, we typically use natural log (normalised count + 1) for RNA-seq data:

1, load airway data

library(airway)
library(magrittr)

data('airway')
airway$dex %<>% relevel('untrt')

2, convert Ensembl genes to gene symbols

ens <- rownames(airway)

library(org.Hs.eg.db)
symbols <- mapIds(org.Hs.eg.db, keys = ens,
  column = c('SYMBOL'), keytype = 'ENSEMBL')
symbols <- symbols[!is.na(symbols)]
symbols <- symbols[match(rownames(airway), names(symbols))]
rownames(airway) <- symbols
keep <- !is.na(rownames(airway))
airway <- airway[keep,]

3, differential expression with DESeq2

library('DESeq2')
dds <- DESeqDataSet(airway, design = ~ cell + dex)
dds <- DESeq(dds, betaPrior=FALSE)
res <- results(dds,
  contrast = c('dex','trt','untrt'))
res <- lfcShrink(dds,
  contrast = c('dex','trt','untrt'), res=res, type = 'normal')

4, Prepare results object

EnhancedVolcano will internally transform the y-axis values by -log10(); here, we counteract this by creating a new column in the data for base mean.

res$baseMeanNew <- 1 / (10^log(res$baseMean + 1))

Test it:

log(5 + 1)
-log10(1 / (10 ^ log(5 + 1)))

5, generate an EnhancedMA

# define different cell-types that will be shaded
celltype1 <- c('VCAM1','WNT2')

library(EnhancedVolcano)

EnhancedVolcano(res,
    lab = rownames(res),
    title = 'MA plot',
    subtitle = 'EnhancedMA',
    x = 'log2FoldChange',
    y = 'baseMeanNew',
    xlab = bquote(~Log[2]~ 'fold change'),
    ylab = bquote(~Log[e]~ 'base mean + 1'),
    ylim = c(0, 12),
    pCutoff = 1 / (10 ^ 7),
    FCcutoff = 2.0,
    pointSize = 3.5,
    labSize = 4.0,
    boxedLabels = TRUE,
    colAlpha = 1,
    legendLabels = c('NS', expression(Log[2]~FC),
      'Mean expression', expression(Mean-expression~and~log[2]~FC)),
    legendPosition = 'bottom',
    legendLabSize = 16,
    legendIconSize = 4.0,
    drawConnectors = TRUE,
    widthConnectors = 0.75,
    # encircle
      encircle = celltype1,
      encircleCol = 'black',
      encircleSize = 2.5,
      encircleFill = 'pink',
      encircleAlpha = 1/2) + coord_flip()