R中的积累曲线(非Vegan)

Question

我想要创建累积曲线，特别是度量累积曲线，使用引导和循环。我感兴趣的是抽样(并替换)示例数据集中的地块总数，从1开始直到总数(n=1…)(马克斯)每个样本将被取样1000次。

我不相信像Vegan这样的软件包会对此有所帮助，因为我不是在寻找物种积累曲线，而是需要根据丰度数据和植物物种的保守系数来计算度量标准(如果我错了，请纠正我！)

我的示例数据集是一个矩阵，包含图、植物物种名称、丰度值和c值(保守主义系数)：

https://docs.google.com/spreadsheets/d/1v-93sV4ANUXpObVbtixTo2ZQjiOKemvQ_cfubZPq9L4/edit?usp=sharing

对于每个_n_th样本的1000个迭代，我需要构建一个矩阵，它将保存具有物种名称、丰度和c值的1000个迭代结果，然后从该样本中消除任何重复的物种。对于每一次迭代，我必须计算植被度量。重要的是，我不计算整个1000个迭代的度量，而是为每个单独的迭代计算度量。

最后，我将将这些结果输入到最终结果的矩阵中，其中行为n个n+1。max n和1000列，并为这1000次迭代中的每一次计算度量。然后，我将在迭代中进行平均，然后从这些平均值中创建我想要的度量的累积曲线。

下面包含了我认为有用的代码，其中包含了一个不同的示例数据集，包括我感兴趣的指标。

https://docs.google.com/spreadsheets/d/1GcH2aq3qZzgTv2YkN-uMpnShblgsuKxAPYKH_mLbbh8/edit?usp=sharing

d<-Example2
d<-data.matrix(d)

MEANC<-function(x){
  mean(x, na.rm=TRUE)
}


FQI<-function(x){
  mean(x, na.rm=TRUE)*sqrt(sum(!is.na(x)))
}

RICH<-function(x){
  totalsprich<-sum(x)
  sum(x!=0, na.rm=TRUE)
}

shannon <- function(x){
  totalCov <- sum(x, na.rm=TRUE)
  (sum(x / totalCov * log(x / totalCov), na.rm=TRUE)) * -1

}

#for this particular example, the only two functions (metrics) that will work will be RICH and shannon

nrep<-1000
totalQuads<-nrow(d)

bootResultSD<-data.frame(matrix(nrow=nrep, ncol=totalQuads) )

bootResultMean<-data.frame(matrix(nrow=nrep, ncol=totalQuads)  )

for(j in 1:totalQuads){
  for(i in 1:nrep){

    bootIndex<-sample(1:totalQuads, j, replace=FALSE)

    bootSample<-d[bootIndex, na.rm=TRUE, drop=FALSE]

    VALUES<-apply(bootSample, 1, shannon)
    bootResultSD[i, j]<-sd(VALUES, na.rm=TRUE)
    bootResultMean[i, j]<-mean(VALUES, na.rm=TRUE)
  }
}

VALUES
bootResultSD
bootResultMean

meanDATA <- apply(bootResultMean, 2, mean, na.rm=TRUE)
meanDATASD <- apply(bootResultSD[-1], 2, mean, na.rm=TRUE)

我之前所创建的问题是，它是根据每个地块计算度量，而不是根据每个累积样本积累图表和重新计算度量。到目前为止，我是根据上面的代码得出的，但我不认为这是我需要的：

for(j in 1:totalQuads){
  for(i in 1:nrep){

    bootIndex<-sample(1:totalQuads, 10, replace=TRUE) 

    bootSample<-d[bootIndex, na.rm=TRUE, drop=FALSE]

    booted<-bootSample[!duplicated(bootSample[,2]),]

    bootResultSD[i, j]<-sd(booted, na.rm=TRUE)
    bootResultMean[i, j]<-mean(booted, na.rm=TRUE)

  }
}

我不知道该怎么做才能越过这一点。提前感谢！

Answer 1

更新：

我和一位同事合作，想出一个答案来回答我上面的问题。这是他创建的代码。

#d<-data file 
nrep <- 1000

totalQuads <- length(unique(d$Plot))

boot.result.fqi <- matrix(nrow=nrep, ncol=totalQuads)
boot.result.meanc <- matrix(nrow=nrep, ncol=totalQuads)
boot.result.shannon <- matrix(nrow=nrep, ncol=totalQuads)
boot.result.rich <- matrix(nrow=nrep, ncol=totalQuads)

for(j in 1:totalQuads){

for(i in 1:nrep){

    bootIndex <- sample(1:totalQuads, j, replace=TRUE)

    for(k in 1:length(bootIndex)){

        if(k == 1){
            bootSample <- subset(d, Plot %in% bootIndex[k])
        } else {
            bootSample <- rbind(bootSample, subset(d, Plot %in% bootIndex[k]))
        }

    }
    # bootSample <- subset(d, Plot %in% bootIndex)

    bootSampleUniqSp <- unique(bootSample[c("Species", "C_Value")])

    ## Calculate and store the results

    #Richness
    boot.result.rich[i,j] <- nrow(bootSampleUniqSp)

    #Mean C
    if(boot.result.rich[i,j] > 0){
        boot.result.meanc[i,j] <- mean(bootSampleUniqSp$C_Value)
    } else {
        boot.result.meanc[i,j] <- 0
    # This is the rule I've set up for when there are no species in the quads. Change the rule as you like
    }

    #FQI
    boot.result.fqi[i,j] <- boot.result.meanc[i,j] * sqrt(boot.result.rich[i,j])

    #Shannon
    covers <- aggregate(bootSample$CoverA_1.4mplot, 
by=list(bootSample$Species), sum)
    # covers <- bootSample$CoverA_1.4mplot
    total.cov <- sum(covers$x)
    boot.result.shannon[i,j] <- (sum(covers$x / total.cov * 
log(covers$x / total.cov), na.rm=TRUE)) * -1

}

}

par(mfcol=c(2,2))
boxplot(boot.result.rich, main="Richness")
boxplot(boot.result.meanc, main="Mean C")
boxplot(boot.result.fqi, main="FQI")
boxplot(boot.result.shannon, main="Shannon's index")


# the means across number of quadrats
apply(boot.result.shannon, 2, mean)

summary.dfr <- data.frame(quads=1:totalQuads, 
richness=apply(boot.result.rich, 2, mean), 
meanc=apply(boot.result.meanc, 2, mean), 
fqi=apply(boot.result.fqi, 2, mean), 
shannon=apply(boot.result.shannon, 2, mean)
)