# X is number of taxa to sample from each community to to account for differences in sequencing depth
# Y is number of times the communities are resampled
# matrix is a community matrix -see picante manual
# DD is a phylogenetic distance matrix. You can calculate it from a phylogeny using the cophenetic() command in APE.
# if you want to calculate aMPD, you need to change the abundance.weighted=FALSE, to =TRUE
# if you get this error you need to update your PICANTE package: Error in mpd(samp, dis, abundance.weighted = abundance.weighted) : 
	#unused argument(s) (abundance.weighted = FALSE)
	
# be warned, this script can take awhile depending on the size of your communities!
#
# required R packages: picante, APE, vegan


MPDrarefy=function (matrix,DD,X,Y){

df.list<-vector("list",Y)

for(i in 1:Y){
print(i)

matrixA<-rrarefy(matrix,X)
a<-colSums(matrixA)
b<-a !=0

matrixB<-matrixA[,b==TRUE]

c<-as.character(colnames(matrixB))
DDA<-DD[c,c]
#print(DDA)
#print(matrixB)
df.list[[i]]<-ses.mpd(matrixB,DDA,null.model="taxa.labels",runs=999,abundance.weighted=TRUE)

#print(df.list)

} # close for loop

Res<-Reduce("+",df.list)/length(df.list)

return(Res)

} # close program