LESSON 22: Population dynamics with the logistic equation

In this lesson you will:

What is the logistic equation? The logistic equation is a well-known model of population size as a function of generation number. In the real world, populations cannot grow exponentially for an indefinite period. The logistic model takes this into account by modeling the production per capita as: r*(1 - x(t)). Here the population x(t) at generation t is expressed a fraction of the maximum number possible. Hence, x(t) is a value between 0 and 1. The parameter r is the maximum per capita production possible. The model says that as the population gets closer to the maximum capacity (i.e., x(t) gets close to 1), the production rate per capita falls to 0.

The update equation for the logistic model written in MATLAB syntax is:

x(t+1) = r*(1 - x(t))*x(t)

See <http://en.wikipedia.org/wiki/Logistic_map> for a discussion of the logistic map.

Setup for this lesson:

Contents

>>>EXAMPLE 1: Calculate 100 generations of a logistic population

Create a new cell in your Lesson22Script.m Editor Window and execute:

   r = 1.5;                           % Maximum per capita production
   start = 0.5;                       % Starting value of population
   generations = 100;                 % Calculate 100 generations forward
   x = zeros(generations + 1, 1);     % Set up a vector to hold population
   x(1) = start;                      % Set the starting population
   for t = 1:generations
       x(t+1) = r*(1 - x(t))*x(t);      % Update equation
   end;

You should see r, start, generations, x, and t variables in your Workspace.

>>>EXAMPLE 2: Display the population evolution for r = 1.5

Create a new cell in your Lesson22Script.m Editor Window and execute:

   t = (0:generations)';
   figure('Name', 'Logistic population r = 1.5');
   stairs(t, x);

You should see a Figure Window with an unedited stair plot:

>>>EXAMPLE 3: Develop a calculateLogicistic function based on EXAMPLE 1

    |function x = calculateLogistic(r, start, generations)|

>>>EXAMPLE 4: Call the calculateLogistic function to calculate population

Create a new cell in your Lesson22Script.m Editor Window and execute:

    figure('Name', 'Logistic population r = 3.59');
    stairs(t, calculateLogistic(3.59, 0.5, 100));

You should see a Figure Window with an unedited stair plot:

>>>EXAMPLE 5: Calculate the population for r = 2.5, 3.5, 3.59, and 3.83

Create a new cell in your Lesson22Script.m Editor Window and execute:

   r = [2.5, 3.5, 3.59, 3.83];      % Calculate evolution for different r's
   start = 0.5;                      % Starting value of population
   generations = 100;                % Calculate 100 generations forward
   x = zeros(generations+1, length(r));  % Calculate populations for each r
   for k = 1:length(r)
       x(:, k) = calculateLogistic(r(k), start, generations);
   end;

You should see r, start, generations, x, and k variables in your Workspace.

>>>EXAMPLE 6: Display the population for r = 2.5, 3.5, 3.59, and 3.83

Create a new cell in your Lesson22Script.m Editor Window and execute:

   figure('Name', 'Logistic population for different values of r');
   subplot(2, 2, 1)
   stairs(t, x(:, 1));
   title('r=2.5');
   subplot(2, 2, 2)
   stairs(t, x(:, 2));
   title('r=3.5');
   subplot(2, 2, 3)
   stairs(t, x(:, 3));
   title('r=3.59');
   subplot(2, 2, 4)
   stairs(t, x(:, 4));
   title('r=3.83');

You should see a Figure Window with four unedited stair subplots:

>>>EXAMPLE 7: Calculate the evolution for random starts with r = 3.59

Create a new cell in your Lesson22Script.m Editor Window and execute:

   r = 3.59;
   generations = 100;
   numStarts = 10;
   x = zeros(generations+1, numStarts);
   for k = 1:numStarts
       x(:, k) = calculateLogistic(r, rand(1, 1), generations);
   end;

>>>EXAMPLE 8: Display the evolution for random starts with r = 3.59

Create a new cell in your Lesson22Script.m Editor Window and execute:

   t = (0:generations)';
   figure('Name', 'Logistic population for random starts with r = 3.59');
   plot(t, x, '.');

>>>EXAMPLE 9: Calculate the evolution after throwing out 1000 generations

Create a new cell in your Lesson22Script.m Editor Window and execute:

   r = 3.5;                           % Maximum per capita production
   start = 0.5;                       % Starting value of population
   generations = 100;                 % Calculate 100 generations forward
   transient = 1000;                  % Number of points before saving
   newStart = start;                  % Eliminate transient by discarding
   for t = 1:transient
       newStart = r*(1 - newStart)*newStart;      % Update equation
   end;
   x = zeros(generations + 1, 1);     % Set up a vector to hold population
   x(1) = newStart;                      % Set the starting population
   for t = 1:generations
       x(t+1) = r*(1 - x(t))*x(t);      % Update equation
   end;

>>>EXAMPLE 10: Modify the calculateLogistic function to take a transient parameter

function x = calculateLogistic(r, start, generations, transient)

|if (nargin == 3)|
    |transient = 0;|
|end;|

>>>EXAMPLE 11: Call the calculateLogistic function with a transient

Create a new cell in your Lesson22Script.m Editor Window and execute:

    x1 = calculateLogistic(3.5, 0.5, 100, 1000);

You should see an x1 variable in your Workspace. Compare x and x1 in the Variable Editor.

>>>EXAMPLE 12: Call the calculateLogistic function without a transient

Create a new cell in your Lesson22Script.m Editor Window and execute:

   x2 = calculateLogistic(3.5, 0.5, 100);

You should see an x2 variable in your Workspace. Compare x and x2 in the Variable Editor.

>>>EXAMPLE 13: Plot the logistic map: x(t) versus x(t+1)

Create a new cell in your Lesson22Script.m Editor Window and execute:

    x = calculateLogistic(3.59, 0.05, 100, 1000); % Calculate population
    figure('Name', 'x(t) versus x(t+1) for r = 3.59')
    plot(x(1: end-1, :), x(2: end, :), 'k.')  % Plot x(t) vs x(t+1) using black dots

You should see a Figure Window with an unedited scatter plot:

>>>EXAMPLE 14: Calculate x(t) versus x(t+1) for 10 random starts

Create a new cell in your Lesson22Script.m Editor Window and execute:

   r = 3.59;
   generations = 100;
   transient = 1000;
   numStarts = 10;
   x = zeros(generations + 1, numStarts);
   for k = 1:numStarts
       x(:, k) = calculateLogistic(r, rand(1,1), generations, transient);
   end;

You should see a Figure Window with an unedited scatter plot:

>>>EXAMPLE 15: Display x(k) versus x(k+1) for 10 random starts

Create a new cell in your Lesson22Script.m Editor Window and execute:

    y = linspace(0, 1.0, 100);  % Create a vector of even steps on [0, 1]
    y1 = r.*(1 - y).*y;         % Parabola of logistic map for reference
    figure('Name', 'x(t) versus x(t+1) for r = 3.59, 10 starts')
    hold on
    line([0, 1], [0, 1], 'Color', 'r');   % Plot a red 45% line
    plot(y, y1, 'LineWidth', 3, 'Color', [0.8, 0.8, 0.8]); % In gray
    plot(x(1: end-1, :), x(2: end, :), 'k.')  % Plot x(t) vs x(t+1) using black dots
    legend({'45^o line', 'equally-spaced iterates', ...
        'random starts after transient'});
    hold off

You should see a Figure Window with a scatter plot overlaid on a parabolic background curve:

>>>EXAMPLE 16: Edit the figure of EXAMPLE 15

>>>EXAMPLE 17: Display the logistic map in a new Figure Window

Create a new cell in your Lesson22Script.m Editor Window and execute:

    open LogisticMap.fig;

You should see a Figure Window with an edited scatter plot overlaid on a parabolic background curve:

>>>EXAMPLE 18: Calculate logistic population for 200 equally-space r's

Create a new cell in your Lesson22Script.m Editor Window and execute:

    r = linspace(0, 4, 200); % Create vector of evenly spaced r's in [0, 4]
    generations = 100;       % Run each calculation for 100 generations
    transient = 1000;        % Throw away 1000 steps before each calculation
    x = zeros(generations + 1, length(r)); % Hold values for plotting
    rv = zeros(size(x));     % rv is same size as x and has r values for x
    for k = 1:length(r)      % For each r
       rv(:, k) = r(k);      % Set the kth column to r(k) for plotting
       x(:, k) = calculateLogistic(r(k), rand(1,1), generations, transient);
    end;

You should see r, |generations, transient, x, rv, and k variables in your Workspace.

>>>EXAMPLE 19: Display a bifurcation diagram r versus x(k) for the logistic map

Create a new cell in your Lesson22Script.m Editor Window and execute:

    figure('Name', 'Bifurcation diagram for logistic equation')
    plot(rv, x, '.k')

You should see a Figure Window with an unedited bifurcation diagram:

>>>EXAMPLE 20: Edit the figure of EXAMPLE 19

>>>EXAMPLE 21: Display the logistic map bifurcation diagram

Create a new cell in your Lesson22Script.m Editor Window and execute:

    open LogisticBifurcationDiagram.fig;

You should see a Figure Window with an unedited bifurcation diagram:

SUMMARY OF SYNTAX

This lesson was written by Kay A. Robbins of the University of Texas at San Antonio and last modified on 27-Jan-2009. Please contact krobbins@cs.utsa.edu with comments or suggestions.