求出第三个“总体”距离矩阵的2个距离矩阵之和(生态语境)

Question

我是生态学家，主要使用纯素R包装。

我有两个矩阵(样本x丰度)(见下文数据)：

矩阵1/ nrow= 6重复*24个位点，ncol=15物种丰度(fish)矩阵2/ nrow= 3重复*24个位点，ncol=10物种丰度(无脊椎动物)

两个矩阵中的位点是相同的。我想得到的整体bray不同(考虑两个矩阵)对的网站。我认为有两种选择：

备选方案1，对复制的鱼类和大型无脊椎动物的丰度进行平均处理，将两个平均丰度矩阵(nrow=24sites，ncol=15+10平均丰度)绑定起来，并计算bray。

选项2，对于每个组合，计算bray在两个站点之间的不同，计算站点质心之间的距离。然后对2距离矩阵进行求和。

如果我不清楚，我在下面的R代码中做了这两个操作。

请告诉我选项2是否正确，是否比选项1更合适。

提前谢谢你。

皮埃尔

下面是R代码示例

生成数据

library(plyr);library(vegan)

#assemblage 1: 15 fish species, 6 replicates per site
a1.env=data.frame(
  Habitat=paste("H",gl(2,12*6),sep=""),
  Site=paste("S",gl(24,6),sep=""),
  Replicate=rep(paste("R",1:6,sep=""),24))

summary(a1.env)

a1.bio=as.data.frame(replicate(15,rpois(144,sample(1:10,1))))

names(a1.bio)=paste("F",1:15,sep="")

a1.bio[1:72,]=2*a1.bio[1:72,]

#assemblage 2: 10 taxa of macro-invertebrates, 3 replicates per site

a2.env=a1.env[a1.env$Replicate%in%c("R1","R2","R3"),]

summary(a2.env)

a2.bio=as.data.frame(replicate(10,rpois(72,sample(10:100,1))))

names(a2.bio)=paste("I",1:10,sep="")

a2.bio[1:36,]=0.5*a2.bio[1:36,]


#environmental data at the sit scale

env=unique(a1.env[,c("Habitat","Site")])

env=env[order(env$Site),]

备选案文1，平均丰度和cbind

a1.bio.mean=ddply(cbind(a1.bio,a1.env),.(Habitat,Site),numcolwise(mean))

a1.bio.mean=a1.bio.mean[order(a1.bio.mean$Site),]

a2.bio.mean=ddply(cbind(a2.bio,a2.env),.(Habitat,Site),numcolwise(mean))

a2.bio.mean=a2.bio.mean[order(a2.bio.mean$Site),]

bio.mean=cbind(a1.bio.mean[,-c(1:2)],a2.bio.mean[,-c(1:2)])

dist.mean=vegdist(sqrt(bio.mean),"bray")

选项2，计算质心之间的每个组合距离，并将2距离矩阵相加。

a1.dist=vegdist(sqrt(a1.bio),"bray")

a1.coord.centroid=betadisper(a1.dist,a1.env$Site)$centroids

a1.dist.centroid=vegdist(a1.coord.centroid,"eucl")

a2.dist=vegdist(sqrt(a2.bio),"bray")

a2.coord.centroid=betadisper(a2.dist,a2.env$Site)$centroids

a2.dist.centroid=vegdist(a2.coord.centroid,"eucl")

用Gavin Simpson的fuse()求和两个距离矩阵

dist.centroid=fuse(a1.dist.centroid,a2.dist.centroid,weights=c(15/25,10/25))

总结两个欧氏距离矩阵(感谢Jari修正)

dist.centroid=sqrt(a1.dist.centroid^2 + a2.dist.centroid^2)

下面的“coord.centroid”用于进一步的基于距离的分析(是正确的吗?)

coord.centroid=cmdscale(dist.centroid,k=23,add=TRUE)

比较备选案文1和2

pco.mean=cmdscale(vegdist(sqrt(bio.mean),"bray"))

pco.centroid=cmdscale(dist.centroid)

comparison=procrustes(pco.centroid,pco.mean)

protest(pco.centroid,pco.mean)

Answer 1

一个更容易的解决方案是通过对每个矩阵进行加权，灵活地组合两个不同的矩阵。权重之和为1。对于两个不同的矩阵，融合的不同矩阵是

d.fused = (w * d.x) + ((1 - w) * d.y)

其中w是一个数值标量(长度1向量)权重。如果您没有理由对其中一组不同之处的权重大于另一组，只需使用w = 0.5。

在我的模拟包；fuse()中，我有一个函数可以帮您完成这个任务。?fuse的示例是

 train1 <- data.frame(matrix(abs(runif(100)), ncol = 10))
 train2 <- data.frame(matrix(sample(c(0,1), 100, replace = TRUE),
                      ncol = 10))
 rownames(train1) <- rownames(train2) <- LETTERS[1:10]
 colnames(train1) <- colnames(train2) <- as.character(1:10)

 d1 <- vegdist(train1, method = "bray")
 d2 <- vegdist(train2, method = "jaccard")

 dd <- fuse(d1, d2, weights = c(0.6, 0.4))
 dd
 str(dd)

这一思想被用于监督Kohonen网络(监督SOMs)，将多层数据整合到一个单一的分析中。

模拟与纯素密切合作，因此不会出现并行运行这两个包的任何问题。

Answer 2

平均距离的正确性取决于你如何处理这些距离。在某些应用程序中，您可能会期望它们实际上是距离。也就是说，它们满足一些度量性质，并且与原始数据有一个定义的关系。合并的不同可能不能满足这些要求。

这个问题与部分Mantel型分析的争议有关，即矩形数据在β多样性研究中的研究热点(我指的是红热)的差异分析和分析。我们素食主义者为两者提供了工具，但我认为在大多数情况下，对矩形数据的分析更加健壮和强大。对于矩形数据，我指的是正常抽样单位乘以物种矩阵。素食者地图中最喜欢的基于不同相似的方法在矩形形状上的不同。这些方法包括db-RDA (capscale)、置换MANOVA (adonis)和群内离散度(betadisper)分析.方法包括mantel、anosim、mrpp、meandis等。

不同或距离的平均值通常与原始矩形数据没有清晰的对应关系。也就是说，差异的平均值与数据的平均值不相对应。我认为，一般说来，最好是平均或处理数据，然后得到与转换后的数据不同的地方。

如果您想要结合不同之处，analogue::fuse()风格的方法是最实用的。但是，您应该理解fuse()也将不同的矩阵缩放为相等的最大值。如果在0..1标度中有不同的度量，这通常是一个小问题，除非其中一个数据集更同质，并且与其他数据集的最大差异较小。在fuse()中，它们都是均衡化的，这样就不是简单的平均，而是距离均衡化后的平均。此外，您必须记住，平均异同点通常会破坏几何学，如果您对矩形化数据使用分析方法(adonis、betadisper、capscale，在纯素中)，这将非常重要。

最后是关于几何组合的异同点。标度0.1中的不同指数是A/B型分数，只有当分母相等时，才能直接加两个分数(然后除以求平均值)。如果忽略这一点，直接对分数进行平均，则结果将不等于平均数据中的相同分数。这就是我所说的破坏几何学。一些开放式指数不是分数，可能是可加性的.曼哈顿距离是相加的。欧氏距离是平方差的平方根，它们的平方是相加的，而不是直接的距离。

我演示了这些东西，显示了将两个不同之处相加的效果(平均意味着将结果除以两个或适当的权重)。我以巴罗科罗拉多岛的纯素数据为例，将其划分为两个大小略有不同的子集。保持几何不变的数据子集距离的相加将产生与分析完整数据相同的结果：

library(vegan) ## data and vegdist
library(analogue) ## fuse
data(BCI)
dim(BCI) ## [1]  50 225
x1 <- BCI[, 1:100]
x2 <- BCI[, 101:225]
## Bray-Curtis and fuse: not additive
plot(vegdist(BCI), fuse(vegdist(x1), vegdist(x2), weights = c(100/225, 125/225)))
## summing distances is straigthforward (they are vectors), but preserving
## their attributes and keeping the dissimilarities needs fuse or some trick
## like below where we make dist structure dtmp to be replaced with the result
dtmp <- dist(BCI) ## dist skeleton with attributes
dtmp[] <- dist(x1, "manhattan") + dist(x2, "manhattan")
## manhattans are additive and can be averaged
plot(dist(BCI, "manhattan"), dtmp)
## Fuse rescales dissimilarities and they are no more additive
dfuse <- fuse(dist(x1, "man"), dist(x2, "man"), weights=c(100/225, 125/225))
plot(dist(BCI, "manhattan"), dfuse)
## Euclidean distances are not additive
dtmp[] <- dist(x1) + dist(x2)
plot(dist(BCI), dtmp)
## ... but squared Euclidean distances are additive
dtmp[] <- sqrt(dist(x1)^2 + dist(x2)^2)
plot(dist(BCI), dtmp)
## dfuse would rescale squared Euclidean distances like Manhattan (not shown)

我只考虑了上面的加法，但是如果你不能加，你就不能平均。这是一个味道的问题，如果这是重要的。勇敢的人会把不可能平均的事情平平，但有些人更胆小，想要遵守规则。我宁愿去第二组。

Answer 3

我喜欢这种简单的this answer，但它只适用于添加两个距离矩阵：

d.fused = (w * d.x) + ((1 - w) * d.y)

因此，我编写了自己的代码片段来组合一个多距离矩阵数组(而不仅仅是2)，并使用标准的R包：

# generate array of distance matrices
x <- matrix(rnorm(100), nrow = 5)
y <- matrix(rnorm(100), nrow = 5)
z <- matrix(rnorm(100), nrow = 5)
dst_array <- list(dist(x),dist(y),dist(z))

# create new distance matrix with first element of array
dst <- dst_array[[1]]

# loop over remaining array elements, add them to distance matrix
for (jj in 2:length(dst_array)){
  dst <- dst + dst_array[[jj]]
}

您还可以使用类似于dst_array大小的向量来定义缩放因子。

dst <- dst + my_scale[[jj]] * dst_array[[jj]]