vlambda博客
学习文章列表

1.1 线性相关分析

（1) 基于相关系数的相关分析(Pearson相关，秩相关) 。

（2) 基于模型的相关分析

     ✳线性回归分析（ Pearson相关是线性回归的一种特殊情况）

     ✳偏相关 

     ✳基于线性模型和降维手段的RDA分析 

1.2 非线性相关分析

大多基于模型，logit回归模型，生存分析模型，基于峰值模型和降维手段的CCA等。

1.3 常见的相关系数

（1）Pearson相关：更适合分析两个服从正态分布的变量的相关。

（2）秩相关：如Spearman秩相关、kendal秩相关。它们只考虑排序的相关程度，不考虑具体的值，应用相对广泛，且不要求正态分布。

1.4 相关系数的可视化

（1）相关网络：用有无线条、线条粗细、线条颜色等体现相关系数和显著性。

（2）相关性热图：用颜色深浅，符号等体现相关系数和显著性。

2.1 分析包的安装与加载

1install.packages(“pheatmap”)
2install.packages(“psych”)
3install.packages(“stringr”)
4library(pheatmap)
5library(psych)
6library(stringr)

2.2 导入数据

1#首先需要设置自己的工作目录
2setwd("E:/lesson3")
3#导入otu属水平相对丰度表
4species_table <- read.table("otu_table.Genus.relative.txt",comment.char="",check.names=F,stringsAsFactors=F, header = TRUE, sep = "\t")
5#读取map文件
6#na.strings=“”设置缺失值的字符
7map<-read.table("mapping_file.txt",sep="\t",na.strings="",header = T,row.names=1,comment.char = "",check.names = F,stringsAsFactors = F)

2.3 清洗数据

 1#提取我们需要分析的环境因子所在列数据
 2map<-map[,10:14]
 3colname_of_map<-colnames(map)
 4#区分重复的属名（由于GreenGene数据库的关系，某些不同科分类下的属，会出现同名不同类的情况）
 5species_table$Taxonomy[duplicated(species_table$Taxonomy)]<-paste(species_table$Taxonomy[duplicated(species_table$Taxonomy)],"_1",sep = "")
 6#剔除最后一行（其为无效分类unclassified）
 7species_table<-species_table[-nrow(species_table),]
 8#将行名设置为属分类名
 9rownames(species_table)<-species_table$Taxonomy
10#去除表格里不需要的注释信息，删除第一列（行名）和最后一列（物种注释信息）
11species_table<-species_table[,-c(1,ncol(species_table))]

2.4 合并数据

1#转置
2species_table<-t(species_table)
3sum_of_species<-colSums(species_table)
4#重排species_table行的顺序，使其和map行名一致
5species_table<-species_table[match(rownames(map),rownames(species_table)),]
6#将物种丰度和环境因子信息合并
7merged_table<-data.frame(map,species_table,check.names = F,check.rows = T)

2.5 计算相关系数和显著性

1#计算spearman秩相关
2correlation_results<-corr.test(merged_table,method ="spearman",adjust="fdr")
3#计算pearson相关（此处我们使用spearman秩相关进行后续的分析，需要计算Pearson相关系数的可以参考此代码，不需要的请忽略）
4#correlation_results<-corr.test(merged_table,method ="pearson",adjust="fdr")
5#method ="spearman"指明用秩相关的方法
6#adjust="fdr"校正错误发现率
7#提取相关矩阵和p值矩阵
8r<-correlation_results$r
9p<-correlation_results$p

2.6 筛选相关系数和显著性

1#剔除环境因子之间、微生物之间的相关系数，因为此处我们只需要环境因子和微生物的相关系数）
2r<-r[-c(1:5),-c(6:192)]
3p<-p[-c(1:5),-c(6:192)]
4#选择相关系数显著次数最多的20个物种
5selected_position_of_species<-head(order(colSums(t(p)<0.05),decreasing =T),20)
6#得到筛选后的相关系数，p值矩阵
7r<-r[selected_position_of_species,]
8p<-p[selected_position_of_species,]

2.7 显著性标记

1#自定义显著性标记函数，此处显著性标记使用“*”，你也可以根据自己的习惯设置其他的显著性标记。
2sig_label<-function(x){ifelse(x<0.001,"***",ifelse(x<0.01,"**",ifelse(x<0.05,"*","")))}
3#得到显著性标记矩阵
4sig_matrix<-sig_label(p)

2.8 相关系数可视化

 1pheatmap(r,fontsize=15,border_color = "black",
 2         display_numbers = sig_matrix,fontsize_row =15,fontsize_col = 15,
 3         fontsize_number = 22,
 4         #显著性标记的符号大小
 5         cluster_rows=T,clustering_distance_rows="correlation",
 6         #指明行聚类，聚类依据的距离
 7         cluster_cols=T,clustering_distance_cols="euclidean",
 8         #指明列聚类，聚类依据的距离
 9         clustering_method="centroid")
10         #聚类方法

可视化结果如下图所示：

本文所分享的代码是相关性分析中利用spearman秩相关系数进行菌群和环境因子之间的相关性分析，并通过热图形式进行可视化，相关性热图在微生物组学数据分析中十分常见，是很好的关联分析可视化工具，对这方面感兴趣的小伙伴赶紧操练起来吧！