R语言生态学数据分析

2018-05-14

概述

Ｒ语言作为强大的数据分析工具，在研究领域发挥重要的作用．在生态学研究中，越来越多的研究采用Ｒ语言做为数据分析的工具，在这里整理常用生态学数据分析代码，供参考．

#载入相应的包 #如果未安装，请按下行代码安装 #install.package(c("ggplot2","reshape2") library(ggplot2) library(reshape2) ＃设置工作路径 setwd("/home/cesar/Ecology_R/20180514/") ＃载入数据，本数据是鱼类物种组成的绝对值数据，在excel中保存为csv格式． taxa<-read.csv("taxa_relative.csv",row.names = 1,header = T) ＃将文件中的NA值改为０ taxa[is.na(taxa)]<-0 ＃将数据中的种类根据数量的多少进行排序 taxa<-taxa[order(rowSums(taxa),decreasing = T),] ＃Ｎ值代表选择数量排前10的物种，将剩下的物种合并成其他 N<-10 taxa_list<-rownames(taxa)[1:N] new_x<-rbind(taxa[rownames(taxa)%in% taxa_list,], 其他=rowSums(taxa[!rownames(taxa) %in% taxa_list,])) ＃合并数据 datm<-melt(cbind(new_x,Taxonomy=rownames(new_x)),id.vars = c('Taxonomy')) ＃作图 ggplot(datm,aes(x=variable,y=value,fill=Taxonomy))+ xlab("")+ ylab("")+ geom_bar(position = "fill",stat = 'identity',width = 0.8)

Taxaplot.png

#清理旧数据 rm(list = ls()) #设置工作路径 setwd("/home/cesar/Ecology_R/20180514/") #载入工作包 library(ggplot2) library(ggdendro) library(vegan) #读取数据 taxa<-read.csv("taxa_relative.csv",header = T,row.names = 1) taxa[is.na(taxa)]<-0 #采用Bray Curtis方法，如需要更换其他方法，可在method参数中调整 beta_bray<-vegdist(t(taxa),method="bray") #建树 hc<-hclust(beta_bray) hcd <- as.dendrogram(hc) dend_data <- dendro_data(hcd, type = "rectangle") ＃绘图 ggplot(dend_data$segments) + theme_dendro()+ scale_x_discrete(expand = c(0,1))+ geom_segment(aes(x = x, y = y, xend = xend, yend = yend))+ geom_text(data = dend_data$labels, aes(x, y, label = label), size = 5,check_overlap = T, nudge_y = -0.02)

Cluster_tree.png

rm(list = ls()) setwd("/home/cesar/Ecology_R/20180514/") library(phyloseq) library(ggplot2) #α多样性需要物种的绝对值，同时载入一个关于样品分组的文件 taxa<-read.csv("taxa.csv",header = T,row.names = 1) meta<-read.csv("map.csv",header = T,row.names = 1) #将数据制成phyloseq的格式 TAXA<-otu_table(taxa,taxa_are_rows = T) SAD<-sample_data(meta) phy<-phyloseq(TAXA,SAD) #绘图，measures可以选择不同的方法 plot_richness(phy,x="Group",color="Month",measures=c("Observed", "Chao1", "ACE", "Shannon", "Simpson", "InvSimpson", "Fisher"))

alpha_diversity.png

library(ggbiplot) setwd("/home/cesar/Ecology_R/20180514/") #读取数据 taxa<-read.csv("taxa_relative2.csv",row.names = 1,header = T) meta<-read.csv("map2.csv",header = T,row.names = 1) taxa<-t(taxa) taxa[is.na(taxa)]<-0 #计算PCA值 pca<-prcomp(taxa,scale. = T) #作图 ggbiplot(pca, obs.scale = 2, var.scale = 1,var.axes = F, groups = meta$Group, ellipse = TRUE, circle = TRUE) + scale_color_discrete(name = '') + theme(legend.direction = 'horizontal', legend.position = 'top')

PCA.png