spp.area=function(c , z, ...){ curve(expr = c*x^z , from=1, to=1e5, #ylim=c(1,10000), xlab="Area", ylab="Numero de especies", font.lab=2, lwd=2, col=2, main=paste("c = ",c,"; z = ",z), ...) curve(expr = c*x^z , from=1, to=1e5, #ylim=c(1,10000), xlab="Area", ylab="Numero de especies", font.lab=2, lwd=2, col=2, main=paste("c = ",c,"; z = ",z), log="xy", ...) } chove.chuva=function(Nspp , chuva , abund , tempo){ spp=paste("sp.",1:Nspp) ilha=numeric() riq=numeric() for(i in 1:tempo){ ilha=union(unique(sample(spp,chuva,replace=T,prob=abund/sum(abund))),ilha) riq[i]=length(ilha) } plot(0:tempo,c(0,riq),type="l",lwd=2,bty='n',xlab="tempo",ylab="numero de especies", font.lab=2,col=2,ylim=c(0,Nspp),main=c("Riqueza da ilha",paste("(Nspp=",Nspp," ; chuva=",chuva,")"))) abline(h=Nspp,lty=3) return(riq) } arquip=function(areas, Nspp, chuva.total, abund, tempo){ n.ilhas=length(areas) spp=paste("sp.",1:Nspp) ilhas=paste("Ilha",1:n.ilhas) chuva=round(chuva.total*areas/sum(areas)) riq=numeric() par(mfrow=c(ceiling(sqrt(n.ilhas)),ceiling(sqrt(n.ilhas)))) for(i in 1:n.ilhas){ riq[i]=chove.chuva(Nspp,chuva[i],abund,tempo)[tempo] } names(riq)=ilhas mod=lm(log10(riq)~log10(areas)) x11() plot(areas,riq,log="xy",font.lab=2,pch=16,col=2,bty="l", main=paste("N ilhas=",n.ilhas,"; N spp=",Nspp,"; tempo=",tempo), xlab="Area da ilha",ylab="Numero de especies",ylim=c(1,Nspp)) abline(mod,lty=2) cat("c=",mod[[1]][1],"z=",mod[[1]][2],"\n") return(riq) } dancando.na.chuva=function(Nspp, chuva, abund, tempo, tx.ext){ spp=paste("sp.",1:Nspp,sep="") ilha=numeric() riq=numeric() ilha=unique(sample(spp,chuva,replace=T,prob=abund/sum(abund))) riq[1]=length(ilha) for(i in 2:tempo){ rr=sample(c("V","M"),length(ilha),replace=T,prob=c((1-tx.ext),tx.ext)) roleta=data.frame(ilha,rr) vivos=roleta[roleta$rr=="V",1] ilha=union(sample(spp,chuva,replace=T,prob=abund/sum(abund)),vivos) riq[i]=length(unique(ilha)) } plot(0:tempo,c(0,riq),type="l",lwd=2,bty='n',xlab="tempo",ylab="numero de especies", font.lab=2,col=2,ylim=c(0,Nspp), main=c("Riqueza da ilha",paste("Nspp =",Nspp," ; chuva =",chuva,"; tx.ext = ",tx.ext))) abline(h=Nspp,lty=3) return(riq) } grafeq=function(E , I , P){ S = I*P/(I+E) ; T = I*E/(I+E) plot(0,0,bty="n",xaxt="n",yaxt="n",xlab="Numero de especies", ylab="Taxas",font.lab=2,lwd=2,xlim=c(0,P), ylim=c(0,1)) segments(0,I,P,0, lwd=2) segments(0,0,P,E, lwd=2) abline(v=0) abline(h=0) mtext("P",side=1,at=P,font=2) mtext("I",side=2,at=I,font=2,las=1) mtext("E",side=4,at=E,font=2,las=1) mtext("S",side=1,at=S,font=2,col=2) mtext("T",side=2,at=T,font=2,las=1,col=2) points(S,T,col=2,pch=16,cex=1.3) segments(S,T,S,0,lty=3,col=2) segments(S,T,0,T,lty=3,col=2) } MW=function(areas , dist , P , a=1, b=-.01, c=1, d=-.01){ par(mfrow=c(1,2)) E=a+b*areas I=c+d*dist S=numeric() for(i in 1:length(areas)){S[i]=I[i]*P/(I[i]+E[i])} T=numeric() for(i in 1:length(areas)){T[i]=I[i]*E[i]/(I[i]+E[i])} curve(I[1]-I[1]*x/P,0,P,bty="n",xaxt="n",yaxt="n",xlab="Numero de especies", ylab="Taxas",font.lab=2,lwd=2,ylim=c(0,1)) curve((E[1]/P)*x,0,P,lwd=2,add=TRUE) abline(v=0) abline(h=0) mtext("P",side=1,at=P,font=2) mtext("I1",side=2,at=I[1],font=2,las=1) mtext("E1",side=4,at=E[1],font=2,las=1) mtext("S1",side=1,at=S[1],font=2,col=2) mtext("T1",side=2,at=T[1],font=2,las=1,col=2) points(S[1],T[1],col=2,pch=16,cex=1.3) segments(S[1],T[1],S[1],0,lty=3,col=2) segments(S[1],T[1],0,T[1],lty=3,col=2) for(i in 2:length(areas)){ curve(I[i]-I[i]*x/P,0,P,lwd=2,ylim=c(0,1),add=TRUE,lty=i) curve((E[i]/P)*x,0,P,lwd=2,add=TRUE,lty=i) mtext(paste("I",i,sep=""),side=2,at=I[i],font=2,las=1) mtext(paste("E",i,sep=""),side=4,at=E[i],font=2,las=1) mtext(paste("S",i,sep=""),side=1,at=S[i],font=2,col=i+1) mtext(paste("T",i,sep=""),side=2,at=T[i],font=2,las=1,col=i+1) points(S[i],T[i],col=i+1,pch=16,cex=1.3) segments(S[i],T[i],S[i],0,lty=3,col=i+1) segments(S[i],T[i],0,T[i],lty=3,col=i+1) } ex=data.frame(areas=areas,spp=S,dist=dist) y=lm(S~areas)[[1]][1] z=lm(S~areas)[[1]][2] plot(spp~areas,data=ex,log="xy",font.lab=2,pch=as.character(1:length(areas)),col=2,bty="l", main=c("Equilíbrio de M&W",paste("c = ",round(y,2),"; z = ",round(z,2))), xlab="Area da ilha",ylab="Numero de especies",ylim=c(1,P)) if(length(unique(ex$spp))>1&length(unique(ex$areas))>1){ abline(lm(log10(spp)~log10(areas),data=ex),lty=3) } return(ex) par(mfrow=c(1,2)) }