#Lab06.r  = basic script for density-dependent models

################################
# Discrete-time Logistic model #
################################

# Start with classic model: N(t+1)=N(t) + R*N(t)
#  and ... N(t+1)-N(t)=R*N(t)
# Modify this to allow R to vary with N(t)
# R=R0*(1-N(t)/K) ... where R0=intrinsic growth rate
#   when N(t) <<<< K, R ~= R0
#   when N(t) = K, R=0
R0= 0.5; N0=2; K=150; t=25
N=matrix(0,t,1)
N[1]=N0
R=matrix(0,(t-1),1)
deltaN=matrix(0,(t-1),1)
for (i in 1:(t-1)){
  R[i]=R0*(1-N[i]/K)
  N[i+1]=N[i]+R[i]*N[i]
  deltaN[i]=N[i+1]-N[i]}

par(mfrow=c(2,2))
plot(1:t,N,xlab="time",type='o')
plot(N[1:(t-1)],deltaN/N[1:(t-1)],xlab="N",ylab="per capita change in N", "o")
plot(N[1:(t-1)],deltaN,xlab="N",ylab="N[t+1]-N[t]", "o",
    ylim=c(min(deltaN),max(R0*K/4)))
 abline(h=0)
 abline(h=R0*K/4)
 abline(v=K/2)
plot(N[1:(t-1)],R,xlab="N",type='o')

##############
# Time lags  #
##############

R0= 0.092; N0=2; K=150; tau=4; yrs=150;
N=matrix(0,yrs,1)
N[1:(tau+1)]=N0
R=matrix(0,(yrs-1),1)
R[1:tau]=R0
deltaN=matrix(0,(yrs-1),1)
for (i in (tau+1):(yrs-1)){
  R[i]=R0*(1-N[i-tau]/K)
  N[i+1]=N[i]+R[i]*N[i]
  deltaN[i]=N[i+1]-N[i]}

windows()
par(mfrow=c(2,2))
plot(N[1:(yrs-1)],R,xlab="N")
plot(1:yrs,N,xlab="time",type="o")
plot(N[1:(yrs-1)],deltaN,xlab="N",ylab="N[t+1]-N[t]")
plot(N[1:(yrs-1)],deltaN/N[1:(yrs-1)],xlab="N",ylab="Per Capita Rate of Change")
R0*tau

##################
# theta-logistic #
##################
R0= 2; K=100; 
N=1:K
theta=1
R=R0*(1-(N/K)^theta)

windows()
plot(N,R,xlab="N",type='l',col='black')

theta=5
R=R0*(1-(N/K)^theta)
points(N,R,type='l',col='red')

theta=0.5
N=1:K
R=R0*(1-(N/K)^theta)
points(N,R,type='l',col='blue')
legend(20,.5,lty=1,
  col=c('red','black','blue'),c('theta=5','theta=1','theta=0.5'))