这是群落生态学中的一个经典（也是难题）问题。这picante包具有基于 Miklos & Podani (2004) 的蒙特卡罗算法，该算法在保留边距的同时随机化二进制矩阵。However，该算法假设您start使用保留约束的二进制矩阵；然后，它会根据您的需要提供任意数量的具有这些约束的随机矩阵。
我想不出一种简单的算法来生成二进制矩阵（甚至是非随机的），并满足用作起始值的约束，所以我使用了蛮力 - 我使用r2dtable()生成lots矩阵并发现其中一些是二元矩阵。
开始了 ...
r2dtable 生成起始矩阵
rt <- c(4, 2, 3, 5, 3)
ct <- c(5, 1, 5, 2, 4)
set.seed(101)
system.time(tt <- r2dtable(100000, rt, ct)) ## 0.06 seconds
w <- which(sapply(tt, max) == 1)
length(w) ## 36/100000 (0.036%) are binary
m0 <- tt[[w[1]]] ## pick the first one
如果你不关心效率，你可以到此为止。但是，如果您需要数千个满足条件的矩阵，那么第二阶段更好......
尝试交换洗牌
Miklos 和 Podani 的算法进行“试验交换”（交换行/列对以保留行/列总计）；默认情况下，它会进行 1000 次交换来随机化矩阵。
library(picante)
results <- list(m0)
nm <- 100000
pb <- txtProgressBar(max = nm, style = 3)
system.time(
for (i in 1:nm) {
setTxtProgressBar(pb, i)
results[[i+1]] <- randomizeMatrix(results[[i]], null.model = "trialswap")
}
)
这在我的机器上会在 7 秒内生成 10^5 矩阵（这将花费大约 20 倍的时间）r2dtable()，根据我的计算）。
然而 ...
所有的10^5矩阵都是这样生成的，and全部r2dtable满足约束的矩阵（只有 36 个）是相同的！ （我们通过将整个矩阵折叠成一个二进制字符串来进行比较......）可能只有一个矩阵满足这些约束，或者可能在一个非常大的空间中有一个非常小的数字，因此很难从其中洗牌一个对另一个...
table(sapply(results, paste, collapse = ""))
## 1111100010111111001010111
## 100001
table(sapply(tt[w], paste, collapse = ""))
## 1111100010111111001010111
## 36
最后的想法
- 如果您已经有一个满足约束的（可能是非随机的）矩阵，或者您有一个合理的算法来生成单个实例，那么您就可以开始
- 试验交换算法应该能够很好地扩展（到更大的尺寸）； Miklos 和 Podani (2004) 提出的两个例子分别是 56x28 和 118x80。这
r2dtable()然而，hack 的规模将会非常大。如果你的实际问题比 5x5 大得多，那就不切实际了。 （模拟退火或遗传算法可能会找到一个初始矩阵——从满足行约束的随机矩阵开始，并使用等于列约束平方偏差的目标函数进行洗牌——但我懒得花时间搞清楚了。）
Miklos I. & Podani J. 2004。存在-不存在矩阵的随机化：评论和新算法。生态学 85：86-92。
