6

私は生態学者で、主にビーガン R パッケージを使用しています。

私は2つの行列(サンプルx存在量)を持っています(以下のデータを参照):

マトリックス 1/ nrow= 6replicates*24sites、ncol=15 種の存在量 (魚) マトリックス 2/ nrow= 3replicates*24sites、ncol=10 種の存在量 (無脊椎動物)

サイトは両方のマトリックスで同じです。サイトのペア間で全体的なブレイ・カーティスの非類似度 (両方の行列を考慮) を取得したいと考えています。2 つのオプションが表示されます。

オプション 1、レプリケート (サイト スケールで) の魚類と大型無脊椎動物の存在量を平均し、2 つの平均存在量行列 (nrow=24sites、ncol=15+10 平均存在量) をバインドし、bray-curtis を計算します。

オプション 2、集合体ごとに、サイトのペア間のブレイカーティス非類似度を計算し、サイト重心間の距離を計算します。次に、2 つの距離行列を合計します。

よくわからない場合は、以下の R コードでこれら 2 つの操作を行いました。

選択肢 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() を使用して 2 つの距離行列を合計する

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

2 つのユークリッド距離行列の合計 (Jari Oksanen 補正のおかげ)

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)
4

3 に答える 3

4

より簡単な解決策は、各行列に重みを付けて、2 つの非類似度行列を柔軟に結合することです。重みの合計は 1 になる必要があります。2 つの非類似度行列の場合、融合非類似度行列は次のようになります。

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

ここwで、 は数値スカラー (長さ 1 のベクトル) の重みです。非類似度のセットの 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 ネットワーク (教師あり SOM) で使用され、複数層のデータを 1 つの分析に取り込みます。

analogはveganと緊密に連携するので、2 つのパッケージを並べて実行しても問題はありません。

于 2014-01-25T15:32:18.790 に答える
0

この回答のシンプルさが気に入っていますが、2 つの距離行列の追加にのみ適用されます。

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

そのため、標準の R パッケージを使用して、複数の距離行列 (2 つだけではない) の配列を組み合わせる独自のスニペットを作成しました。

# 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]]
于 2017-06-06T10:04:27.710 に答える