Volcano Plot in R
2
1
Entering edit mode
Rashi ▴ 10
@77197b50
Last seen 13 months ago
India

Hi, this code is not running to generate a volcano plot. So what should I do?

  GSE51674_series_matrix_txt[,c("ID_REF","GSM1581506","GSM1581507","GSM1581508","GSM1581509")]
View(Control)
Control_mean<-rowMeans(Filter(is.numeric, Control))
Control$Control_mean<-Control_mean

DN<-
  GSE51674_series_matrix_txt[,c("ID_REF","GSM1581510", "GSM1581511", "GSM1581512", "GSM1581513", "GSM1581514", "GSM1581517")]
View(DN)
DN_mean<-rowMeans(Filter(is.numeric, DN))
DN$DN_mean<-DN_mean

#Filtering out data for Control(0) and DN(1) 
gsmsA<-
  "0000111111XXXXXX"
smlA<-strsplit(gsmsA, split="")[[1]]
View(smlA)
selA<- which(smlA !="X")
smlA<- smlA[selA]
gseA<- gse[, selA]
View(gseA)

#DEGs Identification for Group A (control vs DL)
gsA <- factor(smlA)
groupsA<- make.names(c("Control","DN"))
levels(gsA)<- groupsA
gseA$group<- gsA
designA<- model.matrix(~group + 0,gseA)
View(designA)
colnames(designA)<- levels(gsA)
View(designA)
fitA<- lmFit(gseA, designA)
View(fitA)
ctsA<- paste(groupsA[1], groupsA[2], sep="-")
cont.matrixA<- makeContrasts(contrasts = ctsA, levels=designA)
fit2A<- contrasts.fit(fitA, cont.matrixA)
View(cont.matrixA)
View(fit2A)
fit2A<- limma::eBayes(fit2A)
limma::topTable(fit2A)
resultsA<- limma::decideTests(fit2A)
summary(resultsA)
limma::topTable(fit2A)
full_resultsA<- limma::topTable(fit2A, number=Inf)
full_resultsA<- tibble::rownames_to_column(full_resultsA,"ID1")
summary(full_resultsA)
full_resultsA$expression = ifelse(full_resultsA$P.Value <= 0.05 & abs(full_resultsA$logFC)>=1, ifelse(full_resultsA$logFC>1,'Up','Down'), 'Insignificant')

#Volcano plot generation for group A
mycolors<- c("blue","grey","red")
names(mycolors)<- c("Down", "Insignificant", "Up")
vol_plotA<- 
  ggplot(data=full_resultsA, aes(x=logFC, y=-log10(P.Value), color=expression, label=GENE_SYMBOL)) + 
  geom_point(alpha=0.4, size=3.5) + geom_text_repel(data=full_resultsA %>% filter(expression == "Down"), hjust=2, min.segment.length=0, max.overlaps=Inf, size=5) + 
  geom_text_repel(data=full_resultsA %>% filter(expression == "Up"), hjust=-2, min.segment.length=0, max.overlaps=Inf, size=5) + 
  scale_color_manual(values=mycolors) + 
  xlim(c(-4.5,4.5)) + 
  geom_vline(xintercept=c(-1,1),lty=4,col="black",lwd=0.8) + 
  geom_hline(yintercept=1.301, lty=4, col="black", lwd=0.8) + 
  labs(x="log2(Fold change)", y="-log10(P-value)", title="Differential expression") + 
  theme_bw() + 
  theme(plot.title=element_text(hjust=0.5), legend.position="right", legend.title=element_blank())
vol_plotA                                                    
separationA<-as_tibble(full_resultsA)
data_up<-separationA %>% filter(full_resultsA['expression']=="Up")
View(data_upA)
data_downA<- separationA %>% filter(full_resultsA['expression']=="Down")
View(data_downA)

Error in geom_point(): ! Problem while computing aesthetics. ℹ Error occurred in the 1st layer. Caused by error in FUN(): ! object 'GENE_SYMBOL' not found Run rlang::last_error() to see where the error occurred.

cellxgenedp • 812 views
ADD COMMENT
1
Entering edit mode
@james-w-macdonald-5106
Last seen 12 hours ago
United States

This support site is intended to provide help for people using Bioconductor packages. While it's true that you are using Bioconductor packages, you are essentially asking people to debug your ggplot2 code, which is what is failing. That is off-topic for this site. Perhaps you could try biostars.org or stackoverflow.com.

An alternative would be to re-read the error message you got and figure out what it is telling you. It seems pretty clear, and learning to read, understand, and act upon error messages is a good skill to have if you plan to use R much.

ADD COMMENT
0
Entering edit mode
@gordon-smyth
Last seen 3 hours ago
WEHI, Melbourne, Australia

Alternatively you could try

volcanoplot(fit2A)

or

plotMD(fit2a, status=resultsA)
ADD COMMENT

Login before adding your answer.

Traffic: 492 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

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

Powered by the version 2.3.6