How to label custom multiple gene lists object into volcano plot?

Question

[choijamtsm]

Hello everyone!

I am new here and I tried to make a basic volcano plot with R.

#import data
res <- read.table("results.txt", header=TRUE)
genes <- read.table("genes.txt", header=TRUE)

# Make a basic volcano plot 
with(res, plot(log2foldchange, -log10(pvalue), pch=20, cex.lab=1.6, main="Volcano Plot", xlim=c(-8,8), col="grey"))

# Add colored points:
with(subset(res, pvalue<.05 & abs(log2foldchange)>1.5), points(log2foldchange, -log10(pvalue), pch=20, col="cyan"))

# Label points with the textxy function from the calibrate plot
library(calibrate)

#Label only 1 gene
with(subset(res, Gene == "WNT7A"), textxy(log2FoldChange, -log10(pvalue), labs=Gene, cex=.8))

#How to insert "genes" into volcano plot (overlap res)?

My question is how to add "genes" which contains multiple gene lists into res? (which should label only genes list)

Thank you so much!

[zx8754]

Provide example input: res and genes objects.

Note, there is a dedicated package for volcanoes - EnhancedVolcano

[choijamtsm]

Thanks for your reply

• Here is the input:

res <- read.table("results.txt", header=TRUE)

gene      log2foldchange    pvalue      padj
Cx3cr1      -8.039238787    1.29E-11     2.78E-116
Trem1       -5.258501816    3.01E-44       1.77E-42
Serpina3f  3.202817674  2.12E-09       2.40E-08
Fcrl1        -4.56714509    2.99E-16       5.91866E-15
AU021092 1.886580391    8.00E-03       0.031416887 and so on..

genes <- read.table("genes.txt", header=TRUE)

gene
Fcrl1
Trem1
Serpina3f and so on..

I am struggling with inject genes$gene into res volcano plot

Thank you

[Kevin Blighe]

Run this:

if (!requireNamespace('BiocManager', quietly = TRUE))
    install.packages('BiocManager')

BiocManager::install('EnhancedVolcano')
require('EnhancedVolcano')

EnhancedVolcano(res,
    lab = res$gene,
    x = 'log2foldchange',
    y = 'pvalue')

More information in the vignette that I wrote: https://bioconductor.org/packages/release/bioc/vignettes/EnhancedVolcano/inst/doc/EnhancedVolcano.html

[choijamtsm]

Thank you so much for your reply. I downloaded 'EnhancedVolcano' package then run above command.

EnhancedVolcano(res,
lab = res$gene,
x = 'log2foldchange',
y = 'pvalue')

But All of the res$gene labeled on volcano plot.

Volcano plot Enhanced Volcano res$gene

Actually my goal is to label custom genes. So I've tried to run following command:

EnhancedVolcano(res,
            lab = res$gene,
            x = 'log2foldchange',
            y = 'pvalue',
            selectLab = gene_list)

But it gives me the following plot! which are only row number's shown? how can I see the gene name?

Trying to labeling custom genes

Thank you

[Kevin Blighe]

But All of the res$gene labeled on volcano plot.

Volcano plot Enhanced Volcano res$gene

You simply need to alter the cut-offs for p-value and log [base 2] fold-change. Please read the manual or take a look at the vignette:

• 4.1 Modify cut-offs for log2FC and P value; specify title; adjust point and label size

---------------------------------------

But it gives me the following plot! which are only row number's shown? how can I see the gene name?

Trying to labeling custom genes

Please show the outputs of the following:

class(res)
head(res)
head(res$gene)
class(gene_list)
head(gene_list)

[choijamtsm]

Thank you,

1. I altered the cut-offs for p-value and log2 fold change and graph looks nice!

altered p-value and log2 fold change

2. Here are the outputs:

class(res)

[1] "data.frame"

head(res)

    gene log2foldchange    pvalue          padj
1    Cx3cr1      -8.039239 1.29e-118 2.780000e-116
2     Trem1      -5.258502  3.01e-44  1.770000e-42
3 Serpina3f       3.202818  2.12e-09  2.400000e-08
4     Fcrl1      -4.567145  2.99e-16  5.918660e-15
5  AU021092       1.886580  8.00e-03  3.141689e-02
6     Gbp11       2.090275  3.65e-14  6.250000e-13

head(res$gene)

[1] Cx3cr1    Trem1     Serpina3f Fcrl1     AU021092  Gbp11    
1594 Levels: 10-Mar 1190002N15Rik 1500017E21Rik 1700003F12Rik 1700009N14Rik 2610528A11Rik 2900026A02Rik 3425401B19Rik 3632451O06Rik ... Zmynd15

class(gene_list)

[1] "character"

head(gene_list)

[1] "Fam122A" "Ccnd2"   "Usp2"    "Arl2"    "Rab9B"   "Dcaf7"

My goal is to only label gene_list in res. I also tried following code but failed!

EnhancedVolcano(res,
            lab = res$gene,
            x = 'log2foldchange',
            y = 'pvalue',
            pCutoff = 10e-16,
            FCcutoff = 1.5,
            selectLab = res$gene[which(res$gene %in% gene_list)])

Thank you~

[Kevin Blighe]

Oh, you just need to do this before running anything else:

res$gene <- as.character(res$gene)

[choijamtsm]

It works now!! thank you so much

[Ram]

Please use the formatting bar (especially the code option) to present your post better. You can use backticks for inline code (`text` becomes text), or select a chunk of text and use the highlighted button to format it as a code block. I've done it for you this time.