% =========================================================================
%  Biological_Models_Logistic.m
% 
%  Continuous and Discrete Logistic Models  
%
%  This is a MATLAB m-file for the Logistic model. 
%  -- first the continuous model for various growth rates is examined.
%  -- second the discrete model for various growth rates is examined.
% 
% Note the difference between the two models
%
% =========================================================================
 
% Clears history and command window 
close all; clear all;
 
% Continouus Logistic Model
% Investigate the long-time behavior of the model using
% different values of N0 and r
 
% model paramters
K = 1; % population is normalized, so capacity is 1
y0 = @(t) K;
% -- various initial populations
N0_1 = 0.2;
N0_2 = 0.5;
N0_3 = 1.5;
% -- various growth rates
r1 = 1.0;
r2 = 1.3;
r3 = 0.7;
 
% Continuous Malthus Model
% --examine various initial values
y1 = @(t) (N0_1*K*exp(r1*t))/(K+N0_1*(exp(r1*t)-1));
y2 = @(t) (N0_2*K*exp(r1*t))/(K+N0_2*(exp(r1*t)-1));
y3 = @(t) (N0_3*K*exp(r1*t))/(K+N0_3*(exp(r1*t)-1));
 
figure (1);
hold on;
fplot (y1, [0, 5], 'b');
fplot (y2, [0, 5], 'r'); 
fplot (y3, [0, 5], 'g'); 
fplot (y0, [0, 5], 'black');
legend('N(0)=0.2','N(0)=0.5','N(0)=1.5','K');
xlabel('Time');
ylabel('Population');
hold off;
 
% --examine various growth rates
y1 = @(t) (N0_2*K*exp(r1*t))/(K+N0_2*(exp(r1*t)-1));
y2 = @(t) (N0_2*K*exp(r2*t))/(K+N0_2*(exp(r2*t)-1));
y3 = @(t) (N0_2*K*exp(r3*t))/(K+N0_2*(exp(r3*t)-1));
 
figure (2);
hold on;
fplot (y1, [0, 5], 'b');
fplot (y2, [0, 5], 'r'); 
fplot (y3, [0, 5], 'g'); 
fplot (y0, [0, 5], 'black');
legend('r=1.0','r=1.5','r=2.0');
xlabel('Time');
ylabel('Population');
hold off;
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
% Discrete Logistic Model
% Investigate the long-time behavior of the model using
% different values of P0 and r
 
% create two vectors of zeros [t, p(t)]
N = 21; % number of all time periods to be checked
y = linspace(0,20,N);
%p = zeros(1,N);
 
% initial values for r, p
K = 1.0; % population is normalized, so capacity is 1
r1 = 1.0;
r2 = 1.5;
r3 = 2.0;
p(1:3,1) = 0.01;
 
% calculate population for all time periods
for i = 2:N
  p(1,i) = p(1,i-1)*(1 + r1*(1-p(1,i-1)/K));  
  p(2,i) = p(2,i-1)*(1 + r2*(1-p(2,i-1)/K));  
  p(3,i) = p(3,i-1)*(1 + r3*(1-p(3,i-1)/K));  
end
 
figure (3);
hold on;
plot (y,p(1,:),'b',y,p(2,:),'r',y,p(3,:),'g'); 
legend('r=1.0','r=1.5','r=2.0');
xlabel('Time Steps');
ylabel('Population');
hold off;
 
% examine chaotic behavior
N = 31; % number of all time periods to be checked
y = linspace(0,30,N);
K = 1; % population is normalized, so capacity is 1
r1 = 2.0;
r2 = 2.2;
r3 = 2.4;
r4 = 2.8;
p(1:4,1) = 0.01;
 
% calculate population for all time periods
for i = 2:N
  p(1,i) = p(1,i-1)*(1 + r1*(1-p(1,i-1)/K));  
  p(2,i) = p(2,i-1)*(1 + r2*(1-p(2,i-1)/K));  
  p(3,i) = p(3,i-1)*(1 + r3*(1-p(3,i-1)/K));  
  p(4,i) = p(4,i-1)*(1 + r4*(1-p(4,i-1)/K));  
end
 
figure (4);
hold on;
plot (y,p(1,:),'b',y,p(2,:),'r',y,p(3,:),'g',y,p(4,:),'black'); 
legend('r=2.0','r=2.2','r=2.4','r=2.8');
xlabel('Time Steps');
ylabel('Population');
hold off;

% Discrete Logistic Model - Bifurcation diagram
%
K = 1.0; % population is normalized, so capacity is 1
rmin = 1.0; % min value of growth rate
rmax = 2.8; % max value of growth rate
p0 = 0.02; % initial value for population
n = 1000; % max number of time steps
icheck = 750; % to check behavior at long times
jmax = 400; % number of different growth rate values
del = (rmax-rmin)/jmax; % growth rate step increase
r = zeros(jmax+1,1); % array of different growth rates to be examined
z = zeros(jmax+1,250); % population values at long times
for j=1:jmax+1
  p = zeros(n+1,1);
  p(1) = p0;
  r(j)=(j-1)*del+rmin;
  rcur = r(j); % current value of growth rate
  for i=1:n
    p(i+1) = p(i)*(1+rcur*(1-p(i)/K));
    % check population at long times
    if (i>icheck) 
      z(j,i-icheck) = p(i+1);
    end
  end
end
% plot bifurcation diagram 
figure (5);
plot(r,z,'b.','MarkerSize',4)
xlabel('r','FontSize',10), ylabel('Population','FontSize',10)
title('Bifurcation diagram for the logistic difference equation')
 