在"permute“R包中定义排列设计

Question

我在R中的"permute“包中准确定义我的排列设计/层次结构时遇到了麻烦。

给出一组假设的样地，在这些样地中，我记录了物种的出现情况，我希望在保持每个样地中的物种数量的同时，还保持整个种子库中个体物种的总体丰富度。

最终，我试图构建一个零分布，该分布在地块级别(每个地块有n个物种)和整个物种池级别(每个物种的总观测值)受到约束。

# build dataset representing the presence/absence of 10 species (columns)
# in 100 plots (rows)

set.seed(123)

dat = matrix(
    sample(c(0,1), size = 100*10, replace = T, prob = c(0.75, 0.25)),
    nrow = 100, 
    ncol = 10) # let this matrix represent the observed data

rowSums(dat) # represents the number of species present in each plot
colSums(dat) # represents the overall number of observations of each species

relative_abund = colSums(dat) / sum(dat) 
    # proportion of occurrences of each species in the entire species pool

# use "permute" package to shuffle species in plots
# while maintaining the total number of species in each plot
# and the relative abundance of all species in the species pool

library(permute)

# single permutation of "plot # 1" maintaining number of species per plot
dat[1, shuffle(dat[1,])] 

# single permutation maintaining total observations of "species # 1"
dat[shuffle(dat[,1]), 1 ] 


# use permutation design/control to shuffle data, such that 

rowSums(permuted_dat) == rowSums(dat)
colSums(permuted_dat) == colSums(dat) # at least approximately

Answer 1

谢谢-加文·辛普森。素食包在这方面效果很好。

对于那些希望在生成的空模型中同时保留行和列频率的人，这里有一个示例。"vegan“中的?commsim帮助文件详细介绍了用于指定null模型的不同算法。

# load library
library(vegan)

# build hypothetical dataset of presence/absence data,
# where rows are plots and columns are species

dat = matrix(
    sample(c(0,1), size = 100*10, replace = T, prob = c(0.75, 0.25)),
    nrow = 100, 
    ncol = 10) 

# build null models by reshuffling the dataset
nm_qs = nullmodel(dat, "quasiswap")
nm_qs$data # show presence/absence if original data were abundances

# 100 iterations/reshuffles
# using the "quasiswap" non-sequential algorithm for binary data
sm_qs = simulate(nm_qs, nsim = 100)

# confirm reshuffling preserved row frequencies (i.e., # of species per plot)
rowSums(nm_qs$data) == rowSums(sm_qs[1:100, 1:10, 1]) # first reshuffle
rowSums(nm_qs$data) == rowSums(sm_qs[1:100, 1:10, 2]) # second reshuffle

# confirm reshuffling preserved column frequencies (i.e., overall species abundances)
colSums(nm_qs$data) == colSums(sm_qs[1:100, 1:10, 1]) # first reshuffle
colSums(nm_qs$data) == colSums(sm_qs[1:100, 1:10, 1]) # second reshuffle

# view single reshuffled community
sm_qs[1:100, 1:10, 1]

# build distribution from null models as needed...