相信大家做微生物项目的时候，可能会碰到CCA以及RDA分析，来计算物种与环境因子的相关性什么的。楼主之前要画图，写了一个代码，贴出来分享给大家（本文数据纯属编造）。
1
画图准备
vegan包以及ggplot2包。
一般的16s项目等，常常需要研究不同样本-微生物-环境因子之间的相关关系，CCA或者RDA分析是经典的分析方法。如果使用R语言进行分析的话，需要的程序包为vegan。此函数需要的数据通常为2个表格，一份为样本-物种表格，一份为样本-环境因子表格。具体格式见下图：
2
导入vegan包及两个表格的数据
   library(vegan)
sp <-read.table(file=file.choose(),sep="\t",header=T,row.names=1)
sp
se <-read.table(file=file.choose(),sep="\t",header=T,row.names=1)
se
3
分析策略选择
CCA和RDA像是双胞胎，做出的图也差不多，那么选择用RDA还是CCA分析呢？
其实两种分析方式的区别在于对应的数据模型不一样，我们可以用“样本-物种”文件做DCA分析（decorana(sp)），根据分析结果中Axis  Lengths的第一轴的大小，如果大于4.0，就应选CCA（基于单峰模型，典范对应分析）；如果在3.0-4.0之间，选RDA和CCA均可；如果小于3.0，RDA的结果会更合理（基于线性模型，冗余分析）。
如下图所示，根据DCA分析结果显示，我们使用RDA分析更合理。
4
绘制RDA图
sp0 <- rda(sp ~ 1, se)  
sp0
plot(sp0)
sp1 <- rda(sp ~ ., se)  
sp1
plot(sp1)
出来的图片也不算难看呀。但是想象一下，假如你的数据中有20个左右的微生物，再加上非常长的拉丁名字，RDA图恐怕就是春运现场了。
5
提取并转换ggplot2所需要的数据
在美化RDA图之前，我们需要先提取vegan包产出的RDA图的各类元素的数据，并在保持原来数据结构的基础上（很重要，我们是美化图形，不是修改数据），转换为ggplot2所需要的数据类型。
new<-sp1$CCA
new
samples<-data.frame(sample=row.names(new$u),RDA1=new$u[,1],RDA2=new$u[,2])
samples
species<-data.frame(spece=row.names(new$v),RDA1=new$v[,1],RDA2=new$v[,2])
species
envi<-data.frame(en=row.names(new$biplot),RDA1=new$biplot[,1],RDA2=new$biplot[,2])
envi
line_x =c(0,envi[1,2],0,envi[2,2],0,envi[3,2],0,envi[4,2],0,envi[5,2],0,envi[6,2])
line_x
line_y =c(0,envi[1,3],0,envi[2,3],0,envi[3,3],0,envi[4,3],0,envi[5,3],0,envi[6,3])
line_y
line_g =c("pH","pH","T","T","S2","S2","NH4","NH4","NO2","NO2","Fe2","Fe2")
line_g
line_data =data.frame(x=line_x,y=line_y,group=line_g)
line_data
6
重绘美化RDA图
好了，接下来就是属于ggplot2的时刻了：
library(ggplot2)

ggplot(data=samples,aes(RDA1,RDA2)) +geom_point(aes(color=sample),size=2) +
geom_point(data=species,aes(shape=spece),size=2) +
geom_text(data=envi,aes(label=en),color="blue") +
geom_hline(yintercept=0) + geom_vline(xintercept=0)+
geom_line(data=line_data,aes(x=x,y=y,group=group),color="green")+
theme_bw() + theme(panel.grid=element_blank())
ggsave("RDA2.PDF")
噔噔噔噔，比之前的图好看一丢丢。（数据多了你就看出来哪个美了，手头有此类数据的可以试试，没有的就编一份数据试试o-o）。
最后，代码及数据都在论坛贴附件里（点击“阅读原文”可以跳转下载）。欢迎讨论，并祝大家好好利用R语言，提高自己的姿势水平。
感谢网友“风陵渡口”的分享，如果你有其他好方法也欢迎分享出来哟~今天就酱紫了……
继续阅读
阅读原文