source("comm.R")
source("distr.R")
source("SpecArea.R")

#type = 1 geometric (param is the proportion of neighbouring species)
#     = 2 lognormal (param is the standard deviation)
#     = 3 Tokeshi power fraction
#     = 4 broken stick (fully random)



#sample from theoretical distributions and plot

type = 2
param = 1
pool_size = 100

size = 5000

a=comm(size,type, param,pool_size,2)
plot(sort(a,decreasing=T),log="y",ylab="Number of individuals",xlab="Species",type="b",pch=20,lty=1,lwd=2,col="red")


# plot species-area relationship due to sampling

sizemin = 50
sizeincrease = 2
nosizes = 10
noreal = 20
type = 1
param = 0.5
pool_size = 100

SpecArea(sizemin,sizeincrease,nosizes,noreal,type, param,pool_size,2) 



#plot several S-A relationships together to compare

sizemin = 50
sizeincrease = 2
nosizes = 10
noreal = 10
type = 3
param = 0.5
pool_size = 20

a=SpecArea(sizemin,sizeincrease,nosizes,noreal,type, param,pool_size,2) 

type = 3
param = 0.5
pool_size = 200
b=SpecArea(sizemin,sizeincrease,nosizes,noreal,type, param,pool_size,2) 

outp = cbind(a[,2],b[,2])
matplot(a[,1],outp,log="xy",ylab="Number of species",xlab="Number of individuals",type="b",pch=20,lty=1,lwd=2,col=c("red","blue","green","black"))