使用R语言绘制九象限图的教程

#===========================================================
#  1. 读入两个差异分析表（分别是转录组和翻译组）并合并；  
##==========================================================
RNA =read.csv("rna.csv")
Ribo=read.csv("lribo.csv")

dim(RNA)
dim(Ribo)

##合并两个表格。选择表头为"id"的列进行合并；
combine= merge(RNA,Ribo,
               by.x="id",
               by.y="id",
               suffixes = c(".RNA",".Ribo") ,
               all.x=F,
               all.y=F)

dim(combine)

#保存合并后的数据；
write.table(combine,"merge.txt",sep="\t",row.names=FALSE)

#====================
#  2. 数据准备 
#====================

dat = read.table("merge.txt",header=T,sep="\t")

# 查看表头;
colnames(dat) 

#取RNA的差异倍数，将作为X轴;
#取Ribo的差异倍数，将作为Y轴;
x=dat$RNA_log2FC   
y=dat$H3K27me3_log2FC 

#计算两个组学差异倍数的相关性，并取2位小数;
cor = round(cor(x,y),2) 
#准备作为图形的标题;
main = paste("correlation=",cor,sep="")  
main

#数据筛选：
#RNAdown:转录组显著下调的行;
RNAdown<-intersect(which(dat$RNA_log2FC <= -2),
                   which(dat$RNA_Pvalue  <= 0.05)) 
#RNAup：转录组显著上调的行;
RNAup<-intersect(which(dat$RNA_log2FC >= 2), 
                 which(dat$RNA_Pvalue   <=0.05)) 
# RNA显著差异的行；
RNA_diff = intersect(which(abs(dat$RNA_log2FC) >= 2),
                     which(dat$RNA_Pvalue <=0.05)) 
# RNA不显著差异的行；
RNA_nodiff = union(which(abs(dat$RNA_log2FC) < 2),
                   which(dat$RNA_Pvalue > 0.05)) 
#Ribo显著下调的行 ；
Ribdown<-intersect(which(dat$Ribo_log2FC <= -2),
                   which(dat$RNA_Pvalue  <= 0.05)) 
#Ribo显著上调的行；
Ribup<-intersect(which(dat$H3K27me3_log2FC >= 2),
                 which(dat$H3K27me3_Pvalue <=0.05)) 
# Ribo显著差异的行；
Ribo_diff = intersect(which(abs(dat$H3K27me3_log2FC) >= 2),
                      which(dat$H3K27me3_Pvalue <=0.05))
# Ribo不显著差异的行；
Ribo_nodiff = union(which(abs(dat$H3K27me3_log2FC) < 2),
                    which(dat$H3K27me3_Pvalue > 0.05))  
# 两个组学都下调；
samedown<-intersect(RNAdown,Ribdown)
# 两个组学都上调；
sameup<-intersect(RNAup,Ribup) 
# RNA下调，Ribo上调；
diff1<-intersect(RNAdown,Ribup) 
# RNA上调，Ribo下调；
diff2<-intersect(RNAup,Ribdown) 

#差异方向相同的基因数量；
homo=length(union(sameup,samedown)) 
#差异方向相反的基因数量；
opp=length(union(diff1,diff2)) 
# 仅仅RNA有差异的基因数量；
RNAonly=length(intersect(Ribo_nodiff,RNA_diff))
# 仅仅Ribo有差异的基因数量；
Riboonly=length(intersect(Ribo_diff,RNA_nodiff)) 

#生成图例所用标签；
RNA_change=paste("Transcription(",RNAonly,")",sep="")
Ribo_change=paste("H3K27me3(",Riboonly,")",sep="")
homo_change=paste("Homodirection(",homo,")",sep="")
opp_change=paste("Opposite(",opp,")",sep="")

#给图例准备颜色
id<-c(RNA_change,Ribo_change,homo_change,opp_change)
co<-c("red","blue","brown","yellow")  

#=====================
#   3. 图形绘制
#=====================
# 绘制所有点散点图，注意可以修改X轴和y轴范围；
# 可利用type = "n",记得改坐标轴范围大小
plot(x,y,pch=16,cex = 1,
     xlim=c(-9,6),ylim=c(-8,8),
     cex.lab=1.2,main=main ,col="grey80",
     xlab="gene expression",
     ylab="ribo expression")

#利用次级函数绘制不同象限的图，并加入图例；

points(x[RNA_diff],y[RNA_diff],pch=16,cex=0.3,col="red")

points(x[Ribo_diff],y[Ribo_diff],pch=16,cex=0.3,col="blue")

points(x[c(sameup,samedown)],y[c(sameup,samedown)],pch=16,cex=0.3,col="#brown")

points(x[c(diff1,diff2)],y[c(diff1,diff2)],pch=16,cex=0.3,col="yellow")

legend("topright",legend=id,col=co,pch=16,bty = "n",cex = 0.8,text.font = 2)

#dev.off()