物种共现矩阵随机化 R答案

【问题标题】：species co-occurrence matrix randomization R物种共现矩阵随机化 R
【发布时间】：2014-11-07 00:08:50
【问题描述】：

我的数据是一个物种共现矩阵，我想生成随机矩阵来测试共现模式。

我发现用于此类分析的唯一函数是 R 包 picante 中的 randomizeMatrix 函数。它运行良好，但是此函数中可用的空模型类型数量有限。

当前实现的空模型（null.model 的参数）：frequency（保持物种出现频率）、richness（保持样本物种丰富度）、independentswap 和试玩

有没有人知道其他函数或对此函数的修改，这将允许我测试其他空模型，例如等概率或比例列总和。

这是我使用函数的方式

> test <- matrix(c(1,1,0,1,0,1,0,0,1,0,0,1,0,1,0,0),nrow=4,ncol=4)
> test
     [,1] [,2] [,3] [,4]
[1,]    1    0    1    0
[2,]    1    1    0    1
[3,]    0    0    0    0
[4,]    1    0    1    0
> randomizeMatrix(test,null.model = "richness",iterations = 1000)
     [,1] [,2] [,3] [,4]
[1,]    1    1    0    0
[2,]    1    1    0    1
[3,]    0    0    0    0
[4,]    0    1    0    1
> randomizeMatrix(test,null.model = "independentswap",iterations = 1000)
     [,1] [,2] [,3] [,4]
[1,]    1    0    1    0
[2,]    1    1    0    1
[3,]    0    0    0    0
[4,]    1    0    1    0
>

我在一个循环中运行该函数以获得多次迭代

提前谢谢你

【问题讨论】：

您看过vegan 包吗...？如果您发布了对您所追求的更具体的定义，那么有人可能会在这里为您编写（例如，生成等概率案例听起来很容易）。您也可以在r-sig-ecology@r-project.org 邮件列表中询问...
感谢您的建议。我想要的是一种对随机化进行更多控制的方法，即为单独的行和列设置不同的选项、等概率、比例或固定总和。我确实设法使用其他软件来做到这一点，但如果能在 R 中提供该选项会很好。
嗯，它当然可以在 R 中完成，但你需要更精确/具体地了解你想要什么。大多数 R 用户不是社区生态学家，但如果您能准确定义所需矩阵的特征，可能会帮助您获得答案。

标签： r

【解决方案1】：

我编写了一个函数来生成一个空的随机矩阵，其中每一列的概率都不同。这是针对物种出现矩阵的。

nullMatrix <- function(nrows, ncols, prob1) {
matrixVector <- c()
for(i in 1:ncols){
    columnVector <- sample(c(0,1), nrows, replace = T, prob = c((1-prob1[i]), prob1[i]))
    matrixVector <- c(matrixVector, columnVector)
}
matrix(matrixVector, nrow = nrows, ncol = ncols)
}

这会生成一个包含指定 n 行和列的矩阵以及每行为 1 的概率向量。

我使用以下代码对其进行了测试。

#generate probability vectore for column matrix, each 10 rows the probability increases by .1
p1 <- c(rep(.1, 10), rep(.2, 10), rep(.3, 10), rep(.4, 10), rep(.5, 10), rep(.6, 10), rep(.7, 10), rep(.8, 10), rep(.9, 10), rep(1, 10) )
#generate the matrix with the nullMatrix function
m1 <- nullMatrix(100, 100, p1)
#The average of every ten rows should roughly = .1, .2, ..., 1
sum(m1[,c(1:10)])/1000 #Should ~.1
sum(m1[,c(11:20)])/1000 #Should ~.2
sum(m1[,c(21:30)])/1000 #Should ~.3
sum(m1[,c(31:40)])/1000 #Should ~.4
sum(m1[,c(41:50)])/1000 #Should ~.5
sum(m1[,c(51:60)])/1000 #Should ~.6
sum(m1[,c(61:70)])/1000 #Should ~.7
sum(m1[,c(71:80)])/1000 #Should ~.8
sum(m1[,c(81:90)])/1000 #Should ~.9
sum(m1[,c(91:100)])/1000 #Should ~1

我的谈话有点晚了，但希望这对某人有所帮助。

【讨论】：