NC：多样性引发的确定性细菌组装限制了群落功能

rm(list=ls()) 
data.set.name = 'test'  
setwd("E:/NTI_files") 
library(picante) 
# read in a rarefied OTU table with all OTUs have abundance of 1 or greater. 
otu = as.data.frame(read.table(paste(data.set.name,"_rarefied_otu.txt",sep=""),header=T,row.names=1));  
# species as rows, samples as columns for otu table
# read in phylo and match 
phylowb = read.tree(paste(data.set.name,"_phylo_tree.nwk",sep="")) 
phylowb 
match.phylowb.otu = match.phylowb.data(phylowb, t(otu)) 
str(match.phylowb.otu)  
write.tree(match.phylowb.otu$phy,paste(data.set.name,"_matched_tree_to_rarified_ot u.tre",sep=""))
phylowbMatch = read.tree(paste(data.set.name,"_tree_matched_to_rarified_otu_table.tre",sep="")) 
phyloMatch 
match.phylowbMatch.otu = match.phylowbMatch.data(phylowbMatch, t(otu)) str(match.phylowbMatch.otu) 
beta.mntd.weighted = as.matrix(comdistnt(t(match.phylowbMatch.otu$data),cophenetic(match.phylowbMat ch.otu$phy),abundance.weighted=T)); 
dim(beta.mntd.weighted); 
beta.mntd.weighted[1:5,1:5]; 
write.csv(beta.mntd.weighted,'betaMNTD_weighted.csv',quote=F); 
identical(colnames(match.phylowbMatch.otu$data),colnames(beta.mntd.weighted)); 
# Just checking and should be TRUE 
identical(colnames(match.phylowbMatch.otu$data),rownames(beta.mntd.weighted)); 
# Just checking and should be TRUE 
# calculate randomized betaMNTD 
beta.reps = 999; # number of randomizations 
random.weighted.bMNTD.comp = array(c(-999),dim=c(ncol(match.phylowbMatch.otu$data),ncol(match.phylowbMatch. otu$data),beta.reps)); 
dim(random.weighted.bMNTD.comp); 
for (rep in 1:beta.reps) { 
  random.weighted.bMNTD.comp[,,rep] = as.matrix(comdistnt(t(match.phylowbMatch.otu$data),taxaShuffle(cophenetic(match. phylowbMatch.otu$phy)),abundance.weighted=T,exclude.conspecifics = F)); 
  print(c(date(),rep)); 
} 
weighted.bNTI = matrix(c(NA),nrow=ncol(match.phylowbMatch.otu$data),ncol=ncol(match.phylowb Match.otu$data)); 
dim(weighted.bNTI); 
for (columns in 1:(ncol(match.phylowbMatch.otu$data)-1)) {   
for (rows in (columns+1):ncol(match.phylowbMatch.otu$data)) { 
random.vals = random.weighted.bMNTD.comp[rows,columns,];     
weighted.bNTI[rows,columns] = (beta.mntd.weighted[rows,columns] - mean(random.vals)) / sd(random.vals);    
rm("random.vals"); 
  }; 
}; 
rownames(weighted.bNTI) = colnames(match.phylowbMatch.otu$data); 
colnames(weighted.bNTI) = colnames(match.phylowbMatch.otu$data); 
weighted.bNTI; 
write.csv(weighted.bNTI,"weighted_bNTI.csv",quote=F);