2020-10-16

关联测度与矩阵

1 相关系数概念

1.1 关联测度

关联测度，简言之是对关联的指标进行测量，检测指标间是否存在相关性。一种是R型相关测量,主要是对变量之间的依赖性（dependence）测度，如协方差或相关系数；另一种是Q型相关测量，常指距离（distance)（相异性）或相似性(similarity)测量。如欧氏距离、Jaccard相似系数；

1.2 皮尔逊相关系数

皮尔逊相关系数是比欧几里德距离更加复杂的可以判断人们兴趣的相似度的一种方法。假设有两个变量X、Y，那么两变量间的皮尔逊相关系数可通过以下公式计算：

其中E是数学期望，cov表示协方差，分母是两个变量的标准差。

1.3 斯皮尔曼相关系数

斯皮尔曼相关系数被定义成等级变量之间的皮尔逊相关系数。即样本容量为X和Y的N个样本，将被转换成n个X和Y的等级数据，相关系数ρ为：

其中E是数学期望，cov表示协方差，分母是两个变量的标准差。

1.4 肯德尔相关系数

无参数假设检验，它使用计算而得的相关系数去检验两个随机变量的统计依赖性。肯德尔相关系数的取值范围在-1到1之间，当τ为1时，表示两个随机变量拥有一致的等级相关性；当τ为-1时，表示两个随机变量拥有完全相反的等级相关性；当τ为0时，表示两个随机变量是相互立的。

X、Y（也可以看做两个集合），它们的元素个数均为N，两个随即变量取的第i（1<=i<=N）个值分别用Xi、Yi表示。X与Y中的对应元素组成一个元素对集合XY，其包含的元素为(Xi, Yi)（1<=i<=N）。当集合XY中任意两个元素(Xi, Yi)与(Xj, Yj)的排行相同时（也就是说当出现情况1或2时；情况1：Xi>Xj且Yi>Yj，情况2：Xi<Xj且Yi<Yj），这两个元素就被认为是一致的。当出现情况3或4时（情况3：Xi>Xj且Yi<Yj，情况4：XiYj），这两个元素被认为是不一致的。当出现情况5或6时（情况5：Xi=Xj，情况6：Yi=Yj），这两个元素既不是一致的也不是不一致的。计算肯德尔相关系数的值如下：

其中C表示XY中拥有一致性的元素对数（两个元素为一对）；D表示XY中拥有不一致性的元素对数。

注意：这一公式仅适用于集合X与Y中均不存在相同元素的情况（集合中各个元素唯一）。

1.5 Jaccard相似系数

用于比较有限样本集之间的相似性与差异性。Jaccard系数值越大，样本相似度越高。给定两个集合A,B，Jaccard系数定义为A与B交集的大小与A与B并集的大小的比值，定义如下：

当集合A，B都为空时，J(A,B)定义为1。与Jaccard系数相关的指标叫做Jaccard 距离，用于描述集合之间的不相似度。公式定义如下：

Jaccard距离越大，样本相似度越低。

1.6 Bray-Curits相异系数

Bray-Curtis 相异度（Bray-Curtisdissimilarity）是生态学中用来衡量不同样地物种组成差异的测度。由J. Roger Bray and John T.Curtis提出。其计算基于样本中不同物种组成的数量特征（多度，盖度，重要值等）。计算公式为：

其中k是物种数（样地-物种矩阵中的物种数），y表示物种多度，i和j表示相比较的两个样地的编号。

2 相关系数应用分析

2.1 导入数据

对以下程序包用install.packages（“包名”）进行安装，再进行加载。

2.2 加载函数包

library(ade4)
library(vegan)  # 应该先加载ade4再加载vegan，以避免一些冲突
library(gclus)
library(cluster)
library(FD)

2.3 数据结构查看

# 导入CSV格式的数据
spe <- read.csv("C:\\Users\\Lenovo\\Desktop\\data_3\\DoubsSpe.csv", row.names=1)
str(spe)

## 'data.frame':    30 obs. of  27 variables:
##  $ CHA: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ TRU: int  3 5 5 4 2 3 5 0 0 1 ...
##  $ VAI: int  0 4 5 5 3 4 4 0 1 4 ...
##  $ LOC: int  0 3 5 5 2 5 5 0 3 4 ...
##  $ OMB: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ BLA: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ HOT: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ TOX: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ VAN: int  0 0 0 0 5 1 1 0 0 2 ...
##  $ CHE: int  0 0 0 1 2 2 1 0 5 2 ...
##  $ BAR: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ SPI: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ GOU: int  0 0 0 1 2 1 0 0 0 1 ...
##  $ BRO: int  0 0 1 2 4 1 0 0 0 0 ...
##  $ PER: int  0 0 0 2 4 1 0 0 0 0 ...
##  $ BOU: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ PSO: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ ROT: int  0 0 0 0 2 0 0 0 0 0 ...
##  $ CAR: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ TAN: int  0 0 0 1 3 2 0 0 1 0 ...
##  $ BCO: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ PCH: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ GRE: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ GAR: int  0 0 0 0 5 1 0 0 4 0 ...
##  $ BBO: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ ABL: int  0 0 0 0 0 0 0 0 0 0 ...
##  $ ANG: int  0 0 0 0 0 0 0 0 0 0 ...

env <- read.csv("C:\\Users\\Lenovo\\Desktop\\data_3\\DoubsEnv.csv", row.names=1)
str(env)

## 'data.frame':    30 obs. of  11 variables:
##  $ das: num  0.3 2.2 10.2 18.5 21.5 32.4 36.8 49.1 70.5 99 ...
##  $ alt: int  934 932 914 854 849 846 841 792 752 617 ...
##  $ pen: num  48 3 3.7 3.2 2.3 3.2 6.6 2.5 1.2 9.9 ...
##  $ deb: num  0.84 1 1.8 2.53 2.64 2.86 4 1.3 4.8 10 ...
##  $ pH : num  7.9 8 8.3 8 8.1 7.9 8.1 8.1 8 7.7 ...
##  $ dur: int  45 40 52 72 84 60 88 94 90 82 ...
##  $ pho: num  0.01 0.02 0.05 0.1 0.38 0.2 0.07 0.2 0.3 0.06 ...
##  $ nit: num  0.2 0.2 0.22 0.21 0.52 0.15 0.15 0.41 0.82 0.75 ...
##  $ amm: num  0 0.1 0.05 0 0.2 0 0 0.12 0.12 0.01 ...
##  $ oxy: num  12.2 10.3 10.5 11 8 10.2 11.1 7 7.2 10 ...
##  $ dbo: num  2.7 1.9 3.5 1.3 6.2 5.3 2.2 8.1 5.2 4.3 ...

spa <- read.csv("C:\\Users\\Lenovo\\Desktop\\data_3\\DoubsSpa.csv", row.names=1)
str(spa)

## 'data.frame':    30 obs. of  2 variables:
##  $ x: int  88 94 102 100 106 112 114 110 136 168 ...
##  $ y: int  7 14 18 28 39 51 61 76 100 112 ...

通过str查看数据结构，得值数据结构为int和num类型，可进一步对数据进行矩阵测度和相关性检验。

2.4 Bray-Curtis相异矩阵

Bray-Curtis相异矩阵定量与半定量数据的相异和距离测度常用Bray-Curtis。Bray-Curtis相异度是生态学中用来衡量不同样地物种组成差异的测度，其计算基于样本中不同物种组成的数量特征（多度，盖度，重要值等）。Bray-curtis计算时，不仅考虑样本中物种的有无，而且还考虑不同物种的相对丰度。

2.4.1 Bray-Curtis相异系数

默认使用Bray-Curtis相异系数，得到前6行的相异相异系数值为：0.6000000，0.6842105， 0.7500000 ，0.8918919 ，0.7500000 ，0.6842105。

spe.db <- vegdist(spe)  # Bray-Curtis相异系数（默认）
head(spe.db)

## [1] 0.6000000 0.6842105 0.7500000 0.8918919 0.7500000 0.6842105

2.4.2 对数化Bray-Curtis相异系数

为什么要对数据进行对数化？

优点：1. 缩小数据的绝对数值，方便计算。例如，每个数据项的值都很大，许多这样的值进行计算可能对超过常用数据类型的取值范围，这时取对数，就把数值缩小，便于计算；2. 取对数后，可以将乘法计算转换称加法计算；3.数据的整个值域中的在不同区间的差异带来的影响不同，对数值小的部分差异的敏感程度比数值大的部分的差异敏感程度更高；4. 取对数之后不会改变数据的性质和相关关系，但压缩了变量的尺度，数据更加平稳，也消弱了模型的共线性、异方差性等。

当然，如果数据集中有负数当然就不能取对数了。实践中，取对数的一般是水平量，而不是比例数据，例如变化率等。

对数转化后物种数据的Bray-Curtis相异矩阵，Bray-Curtis相异矩阵前6行的结果为：0.5509095 ，0.6280761 ，0.7446012， 0.8664653， 0.7657643 ，0.6519273。发现取对数和默认状态下的Bray-Curtis相比，系数值减小了，则说明未取对数前数据中存在一些数值的差异。

spe.dbln <- vegdist(log1p(spe))
head(spe.dbln)

## [1] 0.5509095 0.6280761 0.7446012 0.8664653 0.7657643 0.6519273

2.5 弦距离矩阵

弦距离矩阵，简言之，范数标准化后再计算欧氏距离。即将所有的点假想分布在圆上，对圆上各点的距离进行计算。

对各点进行弦距离矩阵计算，前6行结果为 0.7653669 ,0.9235374, 1.0432236 ,1.2786132 ,1.1153793 ,0.8873116.

spe.norm <- decostand(spe, "nor")
spe.dc <- dist(spe.norm)
head(spe.dc)

## [1] 0.7653669 0.9235374 1.0432236 1.2786132 1.1153793 0.8873116

2.6 Hellinger距离矩阵

Hellinger Distance 被用来衡量两个概率分布之间的相似性.在概率论和统计理论中，Hellinger距离被用来度量两个概率分布的相似度。它是f散度的一种（f散度——度量两个概率分布相似度的指标）。Hellinger距离被定义成Hellinger积分的形式，这种形式由Ernst Hellinger在1909年引进。原理：除以样方多度总和再取平均值后计算欧氏距离。

Hellinger距离矩阵前6行结果为0.8420247， 0.9391305， 1.0616631， 1.2308244 ，1.1153793 ，0.9391305。

spe.hel <- decostand(spe, "hel")
spe.dh <- dist(spe.hel)
head(spe.dh)

## [1] 0.8420247 0.9391305 1.0616631 1.2308244 1.1153793 0.9391305

2.7 Jaccard相异矩阵

Jaccard在分析时，只考虑样本中物种的有无，而不考虑丰度。

注意：所有的二元距离函数在计算系数时，均会自动对数据进行二元转化因此这里的数据不需要二元转化（decostand(，“pa”)）。函数dist.binary（）会自动对数据进行二元转化，但函数vegist（）需要设定参数binary=TRUE。

2.7.1 vegdist函数计算Jaccard相异矩阵

vegdist函数计算Jaccard相异矩阵如下，前6行的计算结果为0.6666667， 0.7500000， 0.8750000 ，0.9090909 ，0.9000000，0.8000000。

spe.dj <- vegdist(spe, "jac", binary=TRUE)
head(spe.dj)

## [1] 0.6666667 0.7500000 0.8750000 0.9090909 0.9000000 0.8000000

2.7.2 dist函数计算Jaccard相异矩阵

dist函数计算Jaccard相异矩阵如下，前6行的计算结果为0.6666667， 0.7500000 ，0.8750000 ，0.9090909 ，0.9000000 ，0.8000000。

spe.dj2 <- dist(spe, "binary")
head(spe.dj2)

## [1] 0.6666667 0.7500000 0.8750000 0.9090909 0.9000000 0.8000000

2.7.3 dist.binary计算Jaccard相异矩阵

使用dist.binary的Jaccard相异矩阵，前6行的计算结果为0.8164966， 0.8660254， 0.9354143 ，0.9534626， 0.9486833 ，0.8944272。

spe.dj3 <- dist.binary(spe, method=1)
head(spe.dj3 )

## [1] 0.8164966 0.8660254 0.9354143 0.9534626 0.9486833 0.8944272

2.7.4 dist.binary计算Sorensen相异矩阵

使用dist.binary（）函数计算Sorensen相异矩阵,得到前6行结果为0.7071068， 0.7745967， 0.8819171， 0.9128709 ，0.9045340， 0.8164966

spe.ds <- dist.binary(spe, method=5)
head(spe.ds)

## [1] 0.7071068 0.7745967 0.8819171 0.9128709 0.9045340 0.8164966

###vegdist计算Sorensen相异矩阵

使用vegdist（）函数计算Sorensen相异矩阵,得到前6行的结果为0.5000000， 0.6000000， 0.7777778， 0.8333333， 0.8181818 ，0.6666667。

spe.ds2 <- vegdist(spe, binary=TRUE)
head(spe.ds2)

## [1] 0.5000000 0.6000000 0.7777778 0.8333333 0.8181818 0.6666667

2.8 环境变量之间的相关系数

2.8.1 Pearson相关系数

环境变量之间的Pearson线性相关系数

env.pearson <- cor(env) #默认 method = "pearson“
round(env.pearson, 2)

##       das   alt   pen   deb    pH   dur   pho   nit   amm   oxy   dbo
## das  1.00 -0.94 -0.38  0.95  0.01  0.70  0.48  0.75  0.41 -0.51  0.39
## alt -0.94  1.00  0.44 -0.87 -0.04 -0.74 -0.44 -0.76 -0.38  0.36 -0.34
## pen -0.38  0.44  1.00 -0.34 -0.22 -0.53 -0.19 -0.31 -0.17  0.31 -0.18
## deb  0.95 -0.87 -0.34  1.00  0.02  0.70  0.39  0.61  0.29 -0.36  0.25
## pH   0.01 -0.04 -0.22  0.02  1.00  0.09 -0.08 -0.05 -0.12  0.18 -0.15
## dur  0.70 -0.74 -0.53  0.70  0.09  1.00  0.36  0.51  0.29 -0.38  0.34
## pho  0.48 -0.44 -0.19  0.39 -0.08  0.36  1.00  0.80  0.97 -0.72  0.89
## nit  0.75 -0.76 -0.31  0.61 -0.05  0.51  0.80  1.00  0.80 -0.63  0.64
## amm  0.41 -0.38 -0.17  0.29 -0.12  0.29  0.97  0.80  1.00 -0.72  0.89
## oxy -0.51  0.36  0.31 -0.36  0.18 -0.38 -0.72 -0.63 -0.72  1.00 -0.84
## dbo  0.39 -0.34 -0.18  0.25 -0.15  0.34  0.89  0.64  0.89 -0.84  1.00

2.8.2 Spearman相关系数

环境变量之间的Spearman相关系数

env.spearman <- cor(env, method = "spearman")   
round(env.spearman, 2)

##       das   alt   pen   deb    pH   dur   pho   nit   amm   oxy   dbo
## das  1.00 -1.00 -0.83  0.99  0.03  0.72  0.77  0.85  0.60 -0.51  0.46
## alt -1.00  1.00  0.83 -0.99 -0.03 -0.72 -0.77 -0.85 -0.60  0.51 -0.46
## pen -0.83  0.83  1.00 -0.82 -0.12 -0.60 -0.77 -0.76 -0.63  0.43 -0.38
## deb  0.99 -0.99 -0.82  1.00  0.04  0.70  0.77  0.83  0.56 -0.46  0.44
## pH   0.03 -0.03 -0.12  0.04  1.00  0.19  0.17  0.01 -0.06  0.18 -0.11
## dur  0.72 -0.72 -0.60  0.70  0.19  1.00  0.64  0.57  0.39 -0.39  0.42
## pho  0.77 -0.77 -0.77  0.77  0.17  0.64  1.00  0.78  0.72 -0.64  0.59
## nit  0.85 -0.85 -0.76  0.83  0.01  0.57  0.78  1.00  0.78 -0.53  0.49
## amm  0.60 -0.60 -0.63  0.56 -0.06  0.39  0.72  0.78  1.00 -0.81  0.68
## oxy -0.51  0.51  0.43 -0.46  0.18 -0.39 -0.64 -0.53 -0.81  1.00 -0.86
## dbo  0.46 -0.46 -0.38  0.44 -0.11  0.42  0.59  0.49  0.68 -0.86  1.00

2.8.3 Kendall相关系数

环境变量之间的Kendall相关系数

env.kendall <- cor(env, method = "kendall") 
round(env.kendall, 2)

##       das   alt   pen   deb    pH   dur   pho   nit   amm   oxy   dbo
## das  1.00 -1.00 -0.63  0.95  0.02  0.57  0.57  0.68  0.42 -0.35  0.31
## alt -1.00  1.00  0.63 -0.95 -0.02 -0.57 -0.57 -0.68 -0.42  0.35 -0.31
## pen -0.63  0.63  1.00 -0.62 -0.08 -0.45 -0.59 -0.58 -0.49  0.31 -0.27
## deb  0.95 -0.95 -0.62  1.00  0.03  0.55  0.56  0.64  0.38 -0.30  0.29
## pH   0.02 -0.02 -0.08  0.03  1.00  0.16  0.14  0.00 -0.05  0.15 -0.10
## dur  0.57 -0.57 -0.45  0.55  0.16  1.00  0.47  0.44  0.29 -0.28  0.30
## pho  0.57 -0.57 -0.59  0.56  0.14  0.47  1.00  0.62  0.58 -0.48  0.45
## nit  0.68 -0.68 -0.58  0.64  0.00  0.44  0.62  1.00  0.64 -0.37  0.36
## amm  0.42 -0.42 -0.49  0.38 -0.05  0.29  0.58  0.64  1.00 -0.67  0.54
## oxy -0.35  0.35  0.31 -0.30  0.15 -0.28 -0.48 -0.37 -0.67  1.00 -0.69
## dbo  0.31 -0.31 -0.27  0.29 -0.10  0.30  0.45  0.36  0.54 -0.69  1.00

2.9 本章汇总

参数	类别	功能
ade4	包	可绘制最小生成树进行对应分析
vegan	包	提供非度量多维尺度分析的包装
gclus	包	针对聚类的散点图和平行坐标图
cluster	包	实现聚类
FD	包	计算功能多样性

2.10 参考文献

[1] 赖江山译. 《数量生态学——R语言的应用》.

[2]https://www.jianshu.com/p/744d71f235c4

源文件下载在线调试

返回：数据分析教程案例, R 语言

关联度测量