LESSON 22: Population dynamics with the logistic equation
In this lesson you will:
- Calculate a population from the logistic model and display the the population as a function of generation.
- Investigate how the population depends on r, the maximum production per capita.
- Invesigate how the population depends on the starting population.
- Display the logistic map x(t) versus x(t+1).
- Create a bifurcation diagram for the logistic map.
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:
- Set the Current Directory to Z:\working\MATLAB\Lesson22. (You will need to make a new directory for Lesson22 .)
- Use File->New->M-File on the MATLAB Window menubar to create a new script called Lesson22Script.m.
Contents
- >>>EXAMPLE 1: Calculate 100 generations of a logistic population
- >>>EXAMPLE 2: Display the population evolution for r = 1.5
- >>>EXAMPLE 3: Develop a calculateLogicistic function based on EXAMPLE 1
- >>>EXAMPLE 4: Call the calculateLogistic function to calculate population
- >>>EXAMPLE 5: Calculate the population for r = 2.5, 3.5, 3.59, and 3.83
- >>>EXAMPLE 6: Display the population for r = 2.5, 3.5, 3.59, and 3.83
- >>>EXAMPLE 7: Calculate the evolution for random starts with r = 3.59
- >>>EXAMPLE 8: Display the evolution for random starts with r = 3.59
- >>>EXAMPLE 9: Calculate the evolution after throwing out 1000 generations
- >>>EXAMPLE 10: Modify the calculateLogistic function to take a transient parameter
- >>>EXAMPLE 11: Call the calculateLogistic function with a transient
- >>>EXAMPLE 12: Call the calculateLogistic function without a transient
- >>>EXAMPLE 13: Plot the logistic map: x(t) versus x(t+1)
- >>>EXAMPLE 14: Calculate x(t) versus x(t+1) for 10 random starts
- >>>EXAMPLE 15: Display x(k) versus x(k+1) for 10 random starts
- >>>EXAMPLE 16: Edit the figure of EXAMPLE 15
- >>>EXAMPLE 17: Display the logistic map in a new Figure Window
- >>>EXAMPLE 18: Calculate logistic population for 200 equally-space r's
- >>>EXAMPLE 19: Display a bifurcation diagram r versus x(k) for the logistic map
- >>>EXAMPLE 20: Edit the figure of EXAMPLE 19
- >>>EXAMPLE 21: Display the logistic map bifurcation diagram
- SUMMARY OF SYNTAX
>>>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
- Use File->New->M-File on the Editor Window menu bar to create a new M-file.
- The first line of this file should be:
|function x = calculateLogistic(r, start, generations)|
- Copy the code from EXAMPLE 1 inside the calculateLogistic function and make the necesesary adjustments.
- Save the file as calculateLogistic.m.
>>>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
- Modify the calculateLogistic function so that its header is:
function x = calculateLogistic(r, start, generations, transient)
- Modify the calculateLogistic function so that it includes the code from EXAMPLE 9 to throw away the transient.
- Add the following lines after the header line in calculateLogistic:
|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
- Edit the figure so that it is properly labeled and visually pleasing.
- Save the figure as LogisticMap.fig.
>>>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
- Edit the figure so that it is appropriately labeled and visually pleasing.
- Save the figure as LogisticBifurcationDiagram.fig.
>>>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
- linspace(a, b, n) creates a vector that consists of n values that are equally spaced in the interval [a, b].
- nargin can be called inside a function to determine the number of arguments that the user has actually passed to the function during the call.
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.