R2-03 第七次作业-火山图

微思微丝 2017-12-17 23:45:32 阅读: 1784

R语言之ggplot画火山图

任务1：

先画一个与主题无关的柱状图：

library(ggplot2)
a<-table(mtcars$gear,mtcars$cyl)
col<-c("lightblue","lightgreen","pink")
png("E:/PNG/R2/R2-7-火山图/R2-03-13.png")
barplot(a,xlab="cyl",ylab="count",col=col,main='R2-03-13')
legend("topleft",legend=rownames(table(mtcars$gear)),pch=22,cex=1.2,
fill=col,bty="n")
dev.off()

再画一个简略图：

library(ggplot2)
a<-ggplot(mtcars,aes(x=factor(cyl),fill=factor(gear),color=factor(gear))
          main="R2-03-14")
b<-a+geom_bar()
b
ggsave("E:/PNG/R2/R2-7-火山图/R2-03-14.png",width=4,height=4)

然后在画一个主题修改图：

library(ggplot2)
mytheme<-
  theme_bw(base_size = 16,base_family = "Times")+##默认为无衬线的Helvetica,默认大小为12. 
  theme(plot.title = element_text(colour = "orange",size=20,face ="bold.italic",hjust=0,vjust = 1))+
  theme(plot.subtitle=element_text(colour="red",size=13,hjust=0,vjust = 0.9))+
  theme(plot.caption = element_text(face="bold.italic",colour="purple",hjust = 1,vjust=0))+
  theme(legend.position = "none")+
  theme(axis.line = element_line(colour = "red"),panel.border = element_blank())+
  theme(panel.grid.major = element_blank(),panel.grid.minor = element_blank())
ggplot(mtcars,aes(x=factor(cyl),fill=factor(gear),color=factor(gear)))+##color=factor(gear)可要可不要
       labs(x="cyl",y="count",title="R2-03-15",
          subtitle="Volcano Plot",caption="HW20171216")+##可写成+xlab("cyl")+ylab("count")
       geom_bar(width = 0.6)+
mytheme  
ggsave("E:/PNG/R2/R2-7-火山图/R2-03-15.png",width=4,height=4)

任务2：

1）通过“ggplot2”和“Cairo”两包来画火山图：（图形太大，不可直接上传，该图为截图）

library(ggplot2)
library(Cairo)
data<-read.csv(file="E:/PNG/R2/R2-7-火山图/R2-第七期作业-火山图/nature23643-s4.csv",header=T)
data$threshold <- as.factor(ifelse(data$q_value < 0.05 & abs(log2(data$FC))>=1,
                                   ifelse(log2(data$FC)> 1 ,'Up','Down'),'Not'))
Cairo(file="E:/PNG/R2/R2-7-火山图/volcan_PNG-4.png", ##Cairo包是R语言的高质量图形渲染库
      units="in",
      bg="white",
      width=4, 
      height=6, 
      pointsize=12, 
      dpi=300)
ggplot(data,aes(x=log2(FC),y =-log2(q_value),colour=threshold,size=threshold))+
  geom_point(alpha=1,size=3)+
  xlim(-2.4,2.6)+ylim(-0.5,15)+
  scale_color_manual(values=c("blue","grey","red"))+
  labs(title = "CMTM6 sgRNA2 vs Control",
       caption = "created by R2-03",
       subtitle="Volcano Plot-R2-03-16")+    
  labs(x = expression(log[2])) +         
  labs(x = expression(paste(log[2],"(fold change)")),
       y = expression(paste(-log[2],"(Q value)")))+
  geom_vline(xintercept=0,lty=5,col="black",lwd=0.6)+
  geom_vline(xintercept=c(-1,1),lty=4,col="grey",lwd=0.5)+ 
  geom_hline(yintercept = -log2(0.5),lty=4,col="grey",lwd=0.5)+
mytheme
dev.off()

2)通过“ggplot2”和“latex2exp”、“ggrepel”三包来画火山图：

library("ggplot2")
library("latex2exp")
library("ggrepel")
data<-read.csv(file="E:/PNG/R2/R2-7-火山图/R2-第七期作业-火山图/nature23643-s4.csv",header=T)
data$threshold <- as.factor(ifelse(data$q_value < 0.05 & abs(log2(data$FC))>=1,
                                   ifelse(log2(data$FC)>1 ,'Up','Down'),'Not'))
data1<-subset(data,threshold== ifelse(log2(data$FC)>1 ,'Up','Down'))##选出上调和下调的数据，最好看下所选数据
data1
Cairo(file="E:/PNG/R2/R2-7-火山图/volcan_R2-03.png",
      units="in",bg="white", width=8, height=6, 
      pointsize=12, dpi=300)
ggplot(data,aes(x = log2(FC) , y = -log2(q_value)))+ 
  geom_point(col = "grey70",alpha = 1)+  
  xlim(-2.5,2.7) +ylim(-0.8,15)+  
  geom_point(data=data1,aes(x = log2(FC),y = -log2(q_value),col = Gene.ID),size =4,
  ##注意有两个geom_point时需注明具体数据名，下同。
             position = position_jitter(height = 0.1,width = 0)) + 
  geom_text_repel(data=data1,aes(label =Gene.ID ),size =3,col = "black",fontface = "italic")+  
  annotate("segment", x = 1 , xend = 2.5, y =-0.6,yend  = -0.6 ,color = "black",
           size = 1,arrow = arrow(angle = 18,ends = "last",type = "closed")) +
  annotate("segment", x = -1 , xend = -2.5, y = -0.6,yend  = -0.6 ,color = "black",
           size = 1,arrow = arrow(angle = 18,ends = "last",type = "closed")) +
  annotate("text",x = 1.7,y = 0.3,label = "Upregulated",size = 4.5) +
  annotate("text",x = -1.7,y = 0.3,label = "Downregulated",size = 4.5)  +
  geom_vline(xintercept=0,lty=5,col="black",lwd=0.6)+
  geom_vline(xintercept=c(-1,1),lty=4,col="grey",lwd=0.5)+ 
  geom_hline(yintercept = -log2(0.5),lty=4,col="grey",lwd=0.5)+
  labs(title = "CMTM6 sgRNA2 vs Control",
       caption = "created by R2-03",
       subtitle="Volcano Plot-R2-03-17")+    
  labs(x = expression(log[2])) +         
  labs(x = expression(paste(log[2],"(fold change)")),
       y = expression(paste(-log[2],"(Q value)")))+
  mytheme
dev.off()