Question: Volcano Plot from DEseq2
1
gravatar for krushnach80
22 months ago by
krushnach80580
krushnach80580 wrote:

Im using this code to make based on log2foldchange and padj value ,im getting the plot but i want those value for my reference how do i extract the same .

alpha <- 0.05 # Threshold on the adjusted p-value
cols <- densCols(res$log2FoldChange, -log10(res$pvalue))
plot(res$log2FoldChange, -log10(res$padj), col=cols, panel.first=grid(),
     main="Volcano plot", xlab="Effect size: log2(fold-change)", ylab="-log10(adjusted p-value)",
     pch=20, cex=0.6)
abline(v=0)
abline(v=c(-1,1), col="brown")
abline(h=-log10(alpha), col="brown")

gn.selected <- abs(res$log2FoldChange) > 2.5 & res$padj < alpha 
text(res$log2FoldChange[gn.selected],
     -log10(res$padj)[gn.selected],
     lab=rownames(res)[gn.selected ], cex=0.4)

when i view gn.selected i get only logical value that is true or false

Any help or suggestion would be highly appreciated

Update I'm doing this

> DF <- DF[DF$log2FoldChange > 1.5 & DF$padj < 0.05,]

is that suffice and am i doing it correctly ?

R • 13k views
ADD COMMENTlink modified 9 months ago by africainterpol1010 • written 22 months ago by krushnach80580
12
gravatar for Kevin Blighe
22 months ago by
Kevin Blighe48k
Kevin Blighe48k wrote:

Edit (October 24, 2018):

This is now a Bioconductor package: EnhancedVolcano: Publication-ready volcano plots with enhanced colouring and labeling

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

Your code appears to run fine on my DESeq2 results objects: yours


I normally do these simple volcano plots a different way:

par(mar=c(5,5,5,5), cex=1.0, cex.main=1.4, cex.axis=1.4, cex.lab=1.4)

topT <- as.data.frame(resultsObject)

#Adjusted P values (FDR Q values)
with(topT, plot(log2FoldChange, -log10(padj), pch=20, main="Volcano plot", cex=1.0, xlab=bquote(~Log[2]~fold~change), ylab=bquote(~-log[10]~Q~value)))

with(subset(topT, padj<0.05 & abs(log2FoldChange)>2), points(log2FoldChange, -log10(padj), pch=20, col="red", cex=0.5))

#with(subset(topT, padj<0.05 & abs(log2FoldChange)>2), text(log2FoldChange, -log10(padj), labels=subset(rownames(topT), topT$padj<0.05 & abs(topT$log2FoldChange)>2), cex=0.8, pos=3))

#Add lines for absolute FC>2 and P-value cut-off at FDR Q<0.05
abline(v=0, col="black", lty=3, lwd=1.0)
abline(v=-2, col="black", lty=4, lwd=2.0)
abline(v=2, col="black", lty=4, lwd=2.0)
abline(h=-log10(max(topT$pvalue[topT$padj<0.05], na.rm=TRUE)), col="black", lty=4, lwd=2.0)

volcano


There is an even better solution that I and colleagues developed using ggplot2, which allows you to easily fit labels into your plot using ggrepel().

ADD COMMENTlink modified 11 months ago • written 22 months ago by Kevin Blighe48k
1

@Kevin wow the author himself ...yes i took the code from this site code

ADD REPLYlink written 22 months ago by krushnach80580
1

Glad to have helped. Note the particular use of the bquote() function in order to get super- and sub-script. These are small modifications but can make a plot feel more professional. bquote() also works with ggplot2 (here I've used it with base functions).

ADD REPLYlink written 22 months ago by Kevin Blighe48k

@Kevin i have some issues Im not able to get the following with the code which you have given as im trying to include

so this is my condition

Up-regulated: =/> 0.58 (Green color)
Down-regulated: =/> -0.59 (Red color)
Significant: =/>1.3 (-Log10 p-Value)

its like all the data points are getting overlapped with same color no distinction

ADD REPLYlink written 19 months ago by krushnach80580
1

Can you paste your code so that I can put this in context?

ADD REPLYlink written 19 months ago by Kevin Blighe48k

okay here is the code im using

par(mar=c(5,5,5,5), cex=1.0, cex.main=1.4, cex.axis=1.4, cex.lab=1.4)
res <- read.csv("volcano.txt", header=TRUE,sep = '\t')

topT <- as.data.frame(res)
head(topT)

with(topT, plot(lfc, -log10(padj), pch=20, main="Volcano plot", cex=1.0, xlab=bquote(~Log[2]~fold~change), ylab=bquote(~-log[10]~P~value)))

with(subset(topT, padj<0.05 & abs(lfc)>=0.58), points(lfc, -log10(padj), pch=20, col="red", cex=0.5))

with(subset(topT, padj<0.05 & abs(lfc)>=-0.59), points(lfc, -log10(padj), pch=20, col="green", cex=0.5))
ADD REPLYlink written 19 months ago by krushnach80580

Yes, the last two lines contradict each other, specifically padj<0.05 & abs(lfc)>=0.58) and padj<0.05 & abs(lfc)>=-0.59).

I think that you want to colour positive fold-change genes as red, and negative as green? Try this:

par(mar=c(5,5,5,5), cex=1.0, cex.main=1.4, cex.axis=1.4, cex.lab=1.4)
res <- read.csv("volcano.txt", header=TRUE,sep = '\t')
topT <- as.data.frame(res)
head(topT)
with(topT, plot(lfc, -log10(padj), pch=20, main="Volcano plot", cex=1.0, xlab=bquote(~Log[2]~fold~change), ylab=bquote(~-log[10]~P~value)))
with(subset(topT, padj<0.05 & lfc>=0.58), points(lfc, -log10(padj), pch=20, col="red", cex=0.5))
with(subset(topT, padj<0.05 & lfc<=-0.59), points(lfc, -log10(padj), pch=20, col="green", cex=0.5))
ADD REPLYlink modified 19 months ago • written 19 months ago by Kevin Blighe48k

"I think that you want to colour positive fold-change genes as red, and negative as green? Try this:" yes thats what i need to show

ADD REPLYlink written 19 months ago by krushnach80580
1

Yep, take a look at the code that I pasted (above). You just needed to remove the abs() function call, and switch the 'greater than' sign to a 'less than' sign on the last line

ADD REPLYlink modified 19 months ago • written 19 months ago by Kevin Blighe48k

Hi Kevin,

I am wondering if you have rscript showing the upreguted and down regulated? Can you share them?

ADD REPLYlink written 15 months ago by bioinfool20

Yes, please take a look here: https://github.com/kevinblighe/EnhancedVolcano

In which program did you conduct differential expression?

ADD REPLYlink written 15 months ago by Kevin Blighe48k

Thanks for the quick reply.

I used edgeR in getting my DE. Once I got them, I pre-filtered my DE list in excel. I only extract FDR<1e-5 which is my significant. Then from those list, my up-reg is log2FC>0 while my down is log2FC<0. I wanted to plot my FDR against log2FC.

Can you embed the script here? Thank you very much!

ADD REPLYlink written 15 months ago by bioinfool20
1

Here is the code for data in an EdgeR results table:

EnhancedVolcanoEdgeR <- function(toptable, NominalCutoff, AdjustedCutoff, LabellingCutoff, FCCutoff, main)
{
    toptable$Significance <- "NS"
    toptable$Significance[(abs(toptable$logFC) > FCCutoff)] <- "FC"
    toptable$Significance[(toptable$FDR<AdjustedCutoff)] <- "FDR"
    toptable$Significance[(toptable$FDR<AdjustedCutoff) & (abs(toptable$logFC)>FCCutoff)] <- "FC_FDR"
    table(toptable$Significance)
    toptable$Significance <- factor(toptable$Significance, levels=c("NS", "FC", "FDR", "FC_FDR"))

    plot <- ggplot(toptable, aes(x=logFC, y=-log10(FDR))) +

        #Add points:
        #   Colour based on factors set a few lines up
        #   'alpha' provides gradual shading of colour
        #   Set size of points
        geom_point(aes(color=factor(Significance)), alpha=1/2, size=0.8) +

        #Choose which colours to use; otherwise, ggplot2 choose automatically (order depends on how factors are ordered in toptable$Significance)
        scale_color_manual(values=c(NS="grey30", FC="forestgreen", FDR="royalblue", FC_FDR="red2"), labels=c(NS="NS", FC=paste("LogFC>|", FCCutoff, "|", sep=""), FDR=paste("FDR Q<", AdjustedCutoff, sep=""), FC_FDR=paste("FDR Q<", AdjustedCutoff, " & LogFC>|", FCCutoff, "|", sep=""))) +

        #Set the size of the plotting window
        theme_bw(base_size=24) +

        #Modify various aspects of the plot text and legend
        theme(legend.background=element_rect(),
            plot.title=element_text(angle=0, size=12, face="bold", vjust=1),

            panel.grid.major=element_blank(),   #Remove gridlines
            panel.grid.minor=element_blank(),   #Remove gridlines

            axis.text.x=element_text(angle=0, size=12, vjust=1),
            axis.text.y=element_text(angle=0, size=12, vjust=1),
            axis.title=element_text(size=12),

            #Legend
            legend.position="top",          #Moves the legend to the top of the plot
            legend.key=element_blank(),     #removes the border
            legend.key.size=unit(0.5, "cm"),    #Sets overall area/size of the legend
            legend.text=element_text(size=8),   #Text size
            title=element_text(size=8),     #Title text size
            legend.title=element_blank()) +     #Remove the title

        #Change the size of the icons/symbols in the legend
        guides(colour = guide_legend(override.aes=list(size=2.5))) +

        #Set x- and y-axes labels
        xlab(bquote(~Log[2]~ "fold change")) +
        ylab(bquote(~-Log[10]~adjusted~italic(P))) +

        #Set the axis limits
        #xlim(-6.5, 6.5) +
        #ylim(0, 100) +

        #Set title
        ggtitle(main) +

        #Tidy the text labels for a subset of genes
        geom_text(data=subset(toptable, FDR<LabellingCutoff & abs(logFC)>FCCutoff),
            aes(label=rownames(subset(toptable, FDR<LabellingCutoff & abs(logFC)>FCCutoff))),
            size=2.25,
            #segment.color="black", #This and the next parameter spread out the labels and join them to their points by a line
            #segment.size=0.01,
            check_overlap=TRUE,
            vjust=1.0) +

        #Add a vertical line for fold change cut-offs
        geom_vline(xintercept=c(-FCCutoff, FCCutoff), linetype="longdash", colour="black", size=0.4) +

        #Add a horizontal line for P-value cut-off
        geom_hline(yintercept=-log10(AdjustedCutoff), linetype="longdash", colour="black", size=0.4)

    return(plot)
}

Requires

  • ggplot2

Execution

  • EnhancedVolcanoDESeq2(topTable, NominalCutoff, AdjustedCutoff, LabellingCutoff, FCCutoff, main)
  • EnhancedVolcanoEdgeR(topTable, NominalCutoff, AdjustedCutoff, LabellingCutoff, FCCutoff, main)

Example

  • EnhancedVolcanoDESeq2(topTable=results, AdjustedCutoff=0.05, LabellingCutoff0.05, FCCutoff=2.0, main="DESeq2 results")
  • EnhancedVolcanoEdgeR(topTable=results, AdjustedCutoff=0.05, LabellingCutoff0.05, FCCutoff=2.0, main="EdgeR results")

Parameters

  • toptable, data-frame of test statistics. Requires at least the following
  • gene names as rownames
  • column named 'log2FoldChange' for DESeq2 or 'logFC' for EdgeR
  • column named 'padj' for DESeq2 or 'FDR' for EdgeR
  • NominalCutoff, nominal p-value cut-off for statistical significance (obsolete)
  • AdjustedCutoff, adjusted p-value cut-off for statistical significance
  • LabellingCutoff, adjusted p-value cut-off for statistical significance for labels
  • FCCutoff, absolute log (base 2) fold change cut-off for statistical significance
  • main, title
ADD REPLYlink written 15 months ago by Kevin Blighe48k

@Kevin

I have this data not row data to do DESeq2 myself

Gene Name   q value A+B: OAC vs normal Log2FC   Minus log10 (q-value) 
PRDM1   0.000113142 -1.18   3.95
TRABD   0.000688272 0.61    3.16
WWC3    0.038536606 -0.46   1.41
LYNX1   0.022322605 -1.54   1.65
ALOX15B 3.97666E-06 -2.02   5.40
RAP2C   0.030696226 -0.35   1.51
FBP2    0.001542994 1.80    2.81
GBP6    1.44846E-08 -2.54   7.84

I don't know how to feed this into your tutorial for volcano plot instead of a top table you have from DESeq2

Can you help me please?

Thank you

ADD REPLYlink written 25 days ago by Angel3.5k
1

Hey, you seem to have all columns that you need

EnhancedVolcano(res1,
    lab = df[,'Gene Name'],
    x = 'Log2FC',
    y = 'q value A+B: OAC vs normal')

Using Q values is not the typical approach for volcanos, but it is no major issue.

ADD REPLYlink written 25 days ago by Kevin Blighe48k

Thanks a lot worked well

Sorry @Kevin I am wondering if I want to put a horizontal line for the threshold used for the q-value (0.05) and two vertical lines (0.6 and -0.6 log2 fold change) for the thresholds used for fold change. The points that represent proteins with q<0.05 and log2FC>0.6 (up-regulated proteins) be red whereas the points with q<0.05 and log2FC<-0.6t (down-regulated proteins) be blue. In your tutorial I am seeing vertical line but how would be with horizontal line?

ADD REPLYlink modified 25 days ago • written 25 days ago by Angel3.5k

Oh, sure, you can add any number of lines with the following parameters:

EnhancedVolcano(
  ...
  hline = c(10e-12, 10e-36, 10e-60, 10e-84),
  hlineCol = c('grey0', 'grey25','grey50','grey75'),
  hlineType = 'longdash',
  hlineWidth = 0.8,
  ...)

For vertical lines, they are:

vline
vlineCol
vlineType
vlineWidth

You can add any number of lines.

Then, there are also the main cut-off lines:

cutoffLineType
cutoffLineCol
cutoffLineWidth

...and the standard ggplot2 engine lines:

gridlines.major = FALSE
gridlines.minor = FALSE

[source: https://github.com/kevinblighe/EnhancedVolcano#adjust-cut-off-lines-and-add-extra-threshold-lines]

ADD REPLYlink modified 25 days ago • written 25 days ago by Kevin Blighe48k

Sorry Kevin

With

> EnhancedVolcano(volcano,
+                 lab = volcano[,'Gene.Name'],
+                 x = 'A.B..OAC.vs.normal.Log2FC',
+                 y = 'q.value',
+                 xlim = c(-8, 8),
+                 title = 'Normal oesophagus versus Tumour',
+                 pCutoff = 0.05,
+                 FCcutoff = 0.6,
+                 transcriptPointSize = 1.5,
+                 transcriptLabSize = 3.0,    col=c('blue', 'blue', 'blue', 'red3'),
+                 colAlpha = 1,    cutoffLineType = 'blank',
+                 cutoffLineCol = 'black',
+                 cutoffLineWidth = 0.8,
+                 hline = c(10e-12, 10e-36, 10e-60, 10e-84),
+                 hlineCol = c('grey0', 'grey25','grey50','grey75'),
+                 hlineType = 'longdash',
+                 hlineWidth = 0.8,
+                 gridlines.major = FALSE,
+                 gridlines.minor = FALSE)
Warning message:
Removed 3 rows containing missing values (geom_hline). 
>

I finished with

enter image description here

But this is far from what I need to put a horizontal line for q-value < 0.05 and two vertical lines for 0.6 and -0.6 log2 fold change

I don't know how to do that can you help please?

ADD REPLYlink modified 20 days ago • written 20 days ago by Angel3.5k
1

Hey F, as you are plotting Q values, you will have to change the y-axis label by adding: ylab = bquote(~-Log[10] ~ italic(Q))

Also, the cut-off lines are not appearing because you have cutoffLineType = 'blank'

Take a look (using a different dataset)

cutoffLineType = 'blank'

EnhancedVolcano(volcano,
  lab = volcano[,'Gene.Name'],
  x = 'A.B..OAC.vs.normal.Log2FC',
  y = 'q.value',
  ylab = bquote(~-Log[10] ~ italic(Q)),
  xlim = c(-8, 8),
  title = 'Normal oesophagus versus Tumour',
  pCutoff = 0.05,
  FCcutoff = 0.6,
  transcriptPointSize = 1.5,
  transcriptLabSize = 3.0,    col=c('blue', 'blue', 'blue', 'red3'),
  colAlpha = 1,
  cutoffLineType = 'blank',
  cutoffLineCol = 'black',
  cutoffLineWidth = 0.8,
  hline = c(10e-12, 10e-36, 10e-60, 10e-84),
  hlineCol = c('grey0', 'grey25','grey50','grey75'),
  hlineType = 'longdash',
  hlineWidth = 0.8,
  gridlines.major = FALSE,
  gridlines.minor = FALSE)

gg

cutoffLineType = 'solid'

EnhancedVolcano(volcano,
  lab = volcano[,'Gene.Name'],
  x = 'A.B..OAC.vs.normal.Log2FC',
  y = 'q.value',
  ylab = bquote(~-Log[10] ~ italic(Q)),
  xlim = c(-8, 8),
  title = 'Normal oesophagus versus Tumour',
  pCutoff = 0.05,
  FCcutoff = 0.6,
  transcriptPointSize = 1.5,
  transcriptLabSize = 3.0,    col=c('blue', 'blue', 'blue', 'red3'),
  colAlpha = 1,
  cutoffLineType = 'solid',
  cutoffLineCol = 'black',
  cutoffLineWidth = 0.8,
  hline = c(10e-12, 10e-36, 10e-60, 10e-84),
  hlineCol = c('grey0', 'grey25','grey50','grey75'),
  hlineType = 'longdash',
  hlineWidth = 0.8,
  gridlines.major = FALSE,
  gridlines.minor = FALSE)

2

Another issue in your code is that the values passed to hline are too low. You could try hline = c(0.10, 0.15, 0.2, 0.25), i.e., in this case, extra cut-offs for Q-values 0.1, 0.15, 0.2, and 0.25

To avail of most recent functionality, you should also install the current dev version:

devtools::install_github('kevinblighe/EnhancedVolcano')
ADD REPLYlink written 20 days ago by Kevin Blighe48k

Thank you so much

I finished with this

enter image description here

Again vertical line is far from and also plot seems a bit crowded that would be nice to show a list of genes instead of all

ADD REPLYlink modified 20 days ago • written 20 days ago by Angel3.5k
1

There are a few ways to further manage the gene labels:

  • use drawConnectors = TRUE
  • use selectLab = c('TGM1', 'EPCAM', 'CLDN3')
  • increase labSize to, e.g., labSize = 5.0
ADD REPLYlink written 20 days ago by Kevin Blighe48k

Sorry Kevin I want

For the non-significant (NS)  light grey 
For the Log2FC  darker grey
For the -Log10Q  darker grey
For the -Log10Q & Log2FC>0  red
For the -Log10Q & Log2FC<0 blue

But this

> EnhancedVolcano(res2,
+                 lab = rownames(res2),
+                 x = 'log2FoldChange',
+                 y = 'pvalue',
+                 xlim = c(-8, 8),
+                 title = 'Tumour versus Normal oesophagus',
+                 pCutoff = 0.05,
+                 FCcutoff = 0.6,col=c('grey75', 'grey50', 'grey25', 'red','blue'),
+                 colAlpha = 1)
>

Gives this enter image description here How I can get blue in left side and red in right side while here I am getting red in both side

ADD REPLYlink modified 11 days ago • written 12 days ago by Angel3.5k
1

i see. In the default colour scheme, there can only be 4 colours. If you need to use 5 colours, you will have to customise it. There is more information here: https://github.com/kevinblighe/EnhancedVolcano#over-ride-colouring-scheme-with-custom-key-value-pairs

Do you think that you can follow that?

ADD REPLYlink written 11 days ago by Kevin Blighe48k

Thank you so much

I did like below but after running R never finished with that I think I am doing something wrong

keyvals <- rep('black', nrow(res2))
 names(keyvals) <- rep('Mid', nrow(res2))
 keyvals[which(res2$log2FoldChange > 0.6 & res2$pvalue < 0.05)] <- 'red'
 names(keyvals)[which(res2$log2FoldChange > 0.6 & res2$pvalue < 0.05)] <- 'high'
 keyvals[which(res2$log2FoldChange < -0.6 & res2$pvalue < 0.05)] <- 'blue'
 names(keyvals)[which(res2$log2FoldChange < -0.6 & res2$pvalue < 0.05)] <- 'low'
 unique(names(keyvals))
[1] "low"  "high" "Mid" 
 unique(keyvals)
[1] "blue"  "red"   "black"
 keyvals[which(abs(res2$log2FoldChange) < 0.6 & res2$pvalue  0.05)] <- 'grey75'
 names(keyvals)[which(abs(res2$log2FoldChange) < 0.6 & res2$pvalue  0.05)] <- 'NS'
 unique(keyvals)
[1] "blue"   "red"    "black"  "grey75"
 keyvals[which(abs(res2$log2FoldChange) > 0.6 & res2$pvalue  0.05)] <- 'grey50'
 names(keyvals)[which(abs(res2$log2FoldChange) > 0.6 & res2$pvalue  0.05)] <- 'log2FoldChange'
 keyvals[which(abs(res2$log2FoldChange) < 0.6 & res2$pvalue < 0.05)] <- 'grey25'
 names(keyvals)[which(abs(res2$log2FoldChange)  < 0.6 & res2$pvalue < 0.05)] <- '-Log10Q'
 unique(keyvals)
[1] "blue"   "red"    "grey25" "grey75" "grey50" "black" 
 unique(names(keyvals))
[1] "low"            "high"           "-Log10Q"       
[4] "NS"             "log2FoldChange" "Mid"           
 EnhancedVolcano(res2,
                                   lab = rownames(res2),
                                   x = 'log2FoldChange',
                                   y = 'pvalue',
                                   selectLab = rownames(res2)[which(names(keyvals) %in% c('NS','log2FoldChange','-Log10Q','low','high'))],
                                   xlim = c(-6.5,6.5),
                                   xlab = bquote(~Log[2]~ 'fold change'),
                                   title = 'Custom colour over-ride',
                                   pCutoff = 0.05,
                                   FCcutoff = 0.6,
                                   colCustom = keyvals,
                                   colAlpha = 1,
                                   legendPosition = 'right',
                                   legendLabSize = 15,
                                   legendIconSize = 5.0,
                                   drawConnectors = TRUE,
                                   widthConnectors = 0.5,
                                   colConnectors = 'grey50',
                                   gridlines.major = TRUE,
                                   gridlines.minor = FALSE,
                                   border = 'partial',
                                   borderWidth = 1.5,
                                   borderColour = 'black')
ADD REPLYlink modified 11 days ago • written 11 days ago by Angel3.5k
1

Hey, here is a reproducible example using the data mentioned on the vignette.

Remember, and assuming that you are using the same data as before, you need to change the ylab parameter to ylab = bquote(~-Log[10] ~ italic(Q)), as you are plotting q-values (adjusted p-values). Also, when deriving the key-value pairs, you need to specify the correct columns, so, for example:

keyvals[which(abs(res2$A.B..OAC.vs.normal.Log2FC) > FC & res2$q.value > p)] <- 'grey50'
names(keyvals)[which(abs(res2$A.B..OAC.vs.normal.Log2FC) > FC & res2$q.value > p)] <- 'log2FoldChange'

Please adapt (modify) the code to suit your data.

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

library('DESeq2')
dds <- DESeqDataSet(airway, design = ~ cell + dex)
dds <- DESeq(dds, betaPrior=FALSE)
res2 <- results(dds,
  contrast = c('cell', 'N061011', 'N61311'))
res2 <- lfcShrink(dds,
  contrast = c('cell', 'N061011', 'N61311'), res=res2)


FC <- 0.6
p <- 10e-3

keyvals <- rep('grey75', nrow(res2))
names(keyvals) <- rep('NS', nrow(res2))

keyvals[which(abs(res2$log2FoldChange) > FC & res2$pvalue > p)] <- 'grey50'
names(keyvals)[which(abs(res2$log2FoldChange) > FC & res2$pvalue > p)] <- 'log2FoldChange'

keyvals[which(abs(res2$log2FoldChange) < FC & res2$pvalue < p)] <- 'grey25'
names(keyvals)[which(abs(res2$log2FoldChange)  < FC & res2$pvalue < p)] <- '-Log10Q'

keyvals[which(res2$log2FoldChange < -FC & res2$pvalue < p)] <- 'blue2'
names(keyvals)[which(res2$log2FoldChange  < -FC & res2$pvalue < p)] <- 'Signif. down-regulated'

keyvals[which(res2$log2FoldChange > FC & res2$pvalue < p)] <- 'red2'
names(keyvals)[which(res2$log2FoldChange > FC & res2$pvalue < p)] <- 'Signif. up-regulated'

unique(keyvals)
unique(names(keyvals))

EnhancedVolcano(res2,
  lab = rownames(res2),
  x = 'log2FoldChange',
  y = 'pvalue',
  #selectLab = rownames(res2)[which(names(keyvals) %in% c('NS','log2FoldChange','-Log10Q','low','high'))],
  xlim = c(-6.5,6.5),
  xlab = bquote(~Log[2]~ 'fold change'),
  ylab = bquote(~-Log[10] ~ italic(P)),
  title = 'Custom colour over-ride',
  pCutoff = 10e-3,
  FCcutoff = 0.6,
  pointSize = 2.5,
  labSize = 4.5,
  #shape = c(6, 4, 2, 11, 15),
  colCustom = keyvals,
  colAlpha = 0.75,
  legendPosition = 'right',
  legendLabSize = 15,
  legendIconSize = 5.0,
  drawConnectors = FALSE,
  widthConnectors = 0.5,
  colConnectors = 'grey50',
  gridlines.major = TRUE,
  gridlines.minor = FALSE,
  border = 'partial',
  borderWidth = 1.5,
  borderColour = 'black')

aaa

ADD REPLYlink written 11 days ago by Kevin Blighe48k

Thank you so much to be this much helpful

I used

  > EnhancedVolcano(res2,
+                 lab = rownames(res2),
+                 x = 'log2FoldChange',
+                 y = 'pvalue',
+                 selectLab = c('DSG2','CALML4','MLEC','PKP2','VIL1','EPCAM','KRT8','FUT8','MUC13','TSPAN8','TBXAS1','GBP6','CLIC3','TACC2','TGM1','CRABP2','MAST4','FASN','ALDH3A1','ACPP','ANXA8L2'),
+                 xlim = c(-6.5,6.5),
+                 xlab = bquote(~Log[2]~ 'fold change'),
+                 ylab = bquote(~-Log[10] ~ italic(Q)),
+                 title = 'Tumour versus normal oesophagus',
+                 pCutoff = 0.05,
+                 FCcutoff = 0.6,
+                 #shape = c(6, 4, 2, 11, 15),
+                 colCustom = keyvals,
+                 colAlpha = 4/5,
+                 legend=c('NS','Log (base 2) fold-change','P value',
+                          'P value & Log (base 2) fold-change'),
+                 legendPosition = 'right',
+                 legendLabSize = 20,
+                 legendIconSize = 10,
+                 drawConnectors = TRUE,
+                 widthConnectors = 1.2,
+                 colConnectors = 'black',boxedLabels=TRUE,labSize =10)

Worked well, but when I want a box around the genes I got error

enter image description here

ADD REPLYlink modified 11 days ago • written 11 days ago by Angel3.5k

boxedLabels was only added in a recent version. You can update to the very latest version with:

devtools::install_github('kevinblighe/EnhancedVolcano')
ADD REPLYlink written 11 days ago by Kevin Blighe48k

Sorry I updated that right now but the same error

Error in EnhancedVolcano(res2, lab = rownames(res2), x = "log2FoldChange",  : 
  unused arguments (labFace = "bold", boxedLabels = TRUE)

Also today I am getting error about labSize argument too

Error in EnhancedVolcano(res2, lab = rownames(res2), x = "log2FoldChange",  : 
  unused arguments (pointSize = 3, labSize = 3)

Is this because of colCustom = keyvals?

ADD REPLYlink modified 11 days ago • written 11 days ago by Angel3.5k
1

I am sorry, but you are therefore not using the latest version. I have just confirmed on my computer that boxedLabels and labSize are functioning correctly.

Please restart your R session, install the development version of EnhancedVolcano, and then confirm (via sessionInfo()) that you are using v1.3.4

Please do not post again until you debug these problems some more.

Thank you! :)

ADD REPLYlink written 11 days ago by Kevin Blighe48k

Sorry you were right, problem was the update

Sorry to be too silly

ADD REPLYlink written 11 days ago by Angel3.5k

No problem, Angel.

ADD REPLYlink written 11 days ago by Kevin Blighe48k

Sorry Kevin

By these lines I plotted the volcano

FC <-0.5
p <- 0.05
keyvals <- rep('grey75', nrow(res_responders))
names(keyvals) <- rep('NS', nrow(res_responders))
keyvals[which(abs(res_responders$log2FoldChange) > FC & res_responders$pvalue > p)] <- 'grey50'
names(keyvals)[which(abs(res_responders$log2FoldChange) > FC & res_responders$pvalue > p)] <- 'log2FoldChange'
keyvals[which(abs(res_responders$log2FoldChange) < FC & res_responders$pvalue < p)] <- 'grey25'
names(keyvals)[which(abs(res_responders$log2FoldChange)  < FC & res_responders$pvalue < p)] <- 'pvalue'
keyvals[which(res_responders$log2FoldChange < -FC & res_responders$pvalue < p)] <- 'blue2'
names(keyvals)[which(res_responders$log2FoldChange  < -FC & res_responders$pvalue < p)] <- 'Signif. down-regulated'
keyvals[which(res_responders$log2FoldChange > FC & res_responders$pvalue < p)] <- 'red2'
names(keyvals)[which(res_responders$log2FoldChange > FC & res_responders$pvalue < p)] <- 'Signif. up-regulated'
unique(keyvals)
unique(names(keyvals))



EnhancedVolcano(res_responders,
                lab = rownames(res_responders),
                x = 'log2FoldChange',
                y = 'pvalue',
                selectLab = c("KRT8","KRT18","EPCAM","AGR2","TSPAN8"),
                xlim = c(-8,8),
                xlab = bquote(~Log[2]~ 'fold change'),
                ylab = bquote(~-Log[10] ~ italic(P)),
                title = 'Responders',
                pCutoff = 0.05,
                FCcutoff = 0.5,
                pointSize = 10,
                labSize = 15,
                #shape = c(6, 4, 2, 11, 15),
                colCustom = keyvals,
                colAlpha = 0.75,
                legendPosition = 'right',
                legendLabSize = 40,
                legendIconSize = 30,
                drawConnectors = FALSE,
                widthConnectors = 0.5,
                colConnectors = 'grey50',
                gridlines.major = TRUE,
                gridlines.minor = FALSE,border = 'partial',
                borderWidth = 3,
                borderColour = 'black',boxedLabels=TRUE,titleLabSize = 50,
                subtitleLabSize = 1,
                captionLabSize = 40,axisLabSize = 40)

enter image description here

But I don't know how to reduce the length of x and particularly y axis

I mean I am seeking for a more shrunk plot

Can you help me please?

ADD REPLYlink modified 1 day ago • written 1 day ago by Angel3.5k
6
gravatar for Renesh
9 months ago by
Renesh1.6k
United States
Renesh1.6k wrote:

Easy Volcano plot for gene expression data in Python https://reneshbedre.github.io/blog/volcano.html

ADD COMMENTlink written 9 months ago by Renesh1.6k

looks pretty cool will give it a try

ADD REPLYlink written 9 months ago by krushnach80580
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1760 users visited in the last hour