
MBW:StageBased Population Model Analysis
From MathBio
Article review and analysis by Luke Pederson
D. T. Crouse, L. B. Crowder, and H. Caswell. A stagebased population model for loggerhead sea turtles and implications for conservation. Ecology, 68(5):14121423, October 1987. *All Information, Unless Otherwise Specified, Comes From This Source
Summary
This wiki is based on a paper discussing population models for Loggerhead sea turtles. The study was conducted to determine the most effective method for turtle population conservation. This paper utilizes Leslie matrices to categorize turtles into 7 lifecycle ages. The paper utilizes differential equations techniques to find steady states and to analyze sensitivity of solutions to population changes at different lifecycle stages. The paper concludes that for turtle conservation, it is most efficient to try and maximize the populations of large juvenile turtles.
This same topic is covered in a paper on conservation of endangered plants. This paper points out a huge disconnect in communication between conversationalists and scientists. This paper takes a similar system and performs sensitivity analysis to determine whether a plant should or should not be considered endangered. It is interesting to see the broad range that sensitivity analysis can cover in the world of mathematical biology.
Introduction
 As people become more aware of the negative effects of development without regulations to balance the marginal social cost and marginal social benefits, there has been a shift in people’s mindset toward preservation or cost efficient use of rare and/or endangered animal species, natural resources, clean air, etc. Before the adoption of legislation to require cost benefit analysis of proposed preservation projects, efforts were typically been focused on simple cost analysis and in many cases have yielded poor results. With the adoption of costbenefit analysis, conservation efforts have been reconsidered to determine if the appropriate action was taken to optimally allocate resources to achieve the desired results. Crouse, Crowder, and Caswell's 1987 article on the stagebased population models for sea turtles affirms the necessity for examining the benefits of current and proposed conservation efforts and the implications that the findings pose on future efforts.
Methodology
 This article analyzed the impact of historic conservation methods and proposed efforts on the population of loggerhead sea turtles. To accomplish this, a sensitivity analysis on a stagebased population model was performed in order to determine where, in the turtles' life cycle, conservation efforts should be focused to prevent extinction. Stagebased population models are formed by a Leslie or Lefkovitch Matrix and make the assumptions that all members within a group or stage are affected equally by the fecundity, mortality, and growth rates. The overall population growth/decline for these models is determined by examining the magnitude of the largest eigenvalue and the population distribution, in the long run, among stages, is determined by normalization of the eigenvector corresponding to the largest eigenvalue.
Sensitivity Analysis
 Sensitivity of the model to perturbations was analyzed both by simulating changes to specific traits of a stage and examining the variation in the growth rate and population distribution and by determinig the sensitivity and elasticity matrices. The findings of the article suggested that the population of loggerhead turtles is most sensitive to perturbations in the large juvenile and subadult stages and significantly less sensitive to changes in the fecundity rates and egg/hatchling survivability. These findings supported the observation at the time that, even with very careful egg protection, the populations were still falling. Today, after regulation requiring Trawler Efficiency Devices/Turtle Exclusion Devices (TEDs) in the US for shrimp nets, the population of loggerhead turtles has decreased in its rate of decline from ~3% in 1987 [1] to ~1.8 % today [2]. To see another form of sensitivity analysis go to APPM4390:Modelling the dynamics of a complex lifecycle parasite.
Conclusions
 Overall, the Crouse, Crowder, and Caswell's article demonstrated the importance of examining the benefits of conservation efforts on a population by examining its sensitivity to changes in a stage of development. Although current conservation methods may be insufficient to prevent the decline in sea turtle populations, the sensitivity analysis illustrates the significance of attempting to direct conservation efforts at the stages of development that are the most elastic, or sensitive to change.
Context & History
 Crouse, Crowder, and Caswell's article describes a model of the loggerhead turtle to compare the effects of various management techniques on the species survival. Conservation efforts for the loggerhead turtles started largely in the 1960s and continue today. Most efforts until the early 1990s were focused on protecting nest and beach management for nest survival. In the mid to late 1980s biologist began to consider that juvenile and adult mortality may be equally, if not more, important to the survival of the population. This hypothesis was one of the driving factors for the 1987 article. Many biologists have made observations to try to determine the rates of remigration, survival, mortality, etc. This article primarily utilized data from Richardson's 20 year Little Cucumber Island project with Frazer's analysis of the data [1]. Many articles have been written since this article was published about conservation methods and population models based on similar techniques and often citing the Crouse, Crowder, and Caswell article. One of the more recent articles written was Enneson and Litzgus's article Using longterm data and a stage classified matrix to assess conservation strategies for an endangered turtle (Clemmys guttata) in 2008. This 2008 article used the framework and many of the same techniques as the 1987 paper to model the spotted turtle population using thirty years of data from Georgian Bay, Ontario [3].
Mathematical Model
StageBased Population Models
 Stagebased models are based on the concept that a large group can be broken into smaller relatively homogeneous groups that have similar parameters for each group and are somehow interconnected with the other groups. For a population model, it is common for the model to split up a group by age or characteristics. For example, a group of people could be broken into age groups such as infants (02), children (312), adolescents (1317), adults (1849), and elderly adults (50+). As time progresses, there must be an ability for each group to move to another group, remain in the same group, or produce offspring. If the probability of remaining in the same group is P, moving to the group is G, and producing offspring is F, the population model can be formed. Recalling the group of people example, if all groups older than the children group produce one child per couple per ten year span, the probability F is 0.05 of a member of that group having a child. If the age cap is 100, people progress at a proportion of the age span, and no one dies in the non elderly adult stage, then the probability P would be [ 0.5, 0.111, 0.25, 0.032] and the probability G would be [ 0.5, 0.889, 0.75, .968, 0.98]. This would yield the matrix shown in equation 1.
Equation 1: Example StageBased Model
 Many assumptions must be made for this model to hold. The population must be completely homogeneous within each group. That is, there must be an equal amount of males and females in each state and each person’s capability to survive, reproduce, stay in the group, and leave the group must be equal for all group members. The initial condition for modeling should be chosen based on historical distributions of the population. The model can then be validated by examining whether the predicted longrun distribution generated by the model closely resembles the initial state. With the model generated, the dominant right eigenvalue of the system can be extracted to determine the total population growth/decay rate. If the eigenvalue is less than one, the population will die out and conversely, if it is larger than one, the population will grow. In the example described by equation 1, the dominant eigenvalue is 1.013. This implies that as time progresses, the rate of population growth will approach ~1.3 %. Scaling the eigenvector to one, the population distribution in the long run would be 3.95 % infants, 15.93 % children, 6.72 % adolescent, 37.30 % adult, and 36.09 % elderly adult. The dominant left eigenvector provides key information about the relative reproductive impact of the various classes. In the example, the dominant left eigenvector (scaled to infants=1) demonstrates that the children contribute 2.0261 times more than the infants, the adolescents contribute 2.2647 times more than the infants, the adults contribute 2.1830 times more than the infants and the elderly contribute 1.5120 times more than the infants to the fecundity of the population.
 The left and right eigenvectors and eigenvalues were calculated by the following equations.
Equation 2: Right Eigenvalues/Eigenvectors
Equation 3: Left Eigenvalues/Eigenvectors
 Sensitivity analysis can be performed on the model by examining the variation in the dominant eigenvalue caused by a change in one of the parameters. For example, the reproduction rate of one of the stages in equation 1 could be varied and the new eigenvalue determined. This could be repeated for all parameters to determine which group and which parameter of that group changes the dominant eigenvalue the most or least depending on what is desired. Suppose one wished to examine the effect of an all out overseas war on population growth where, effectively, all the male adult population is recruited (Ex. A worse version of WWII). Since many may die, the proportion of adults moving into the next category and the proportion staying in the adult category would fall (25% for example). Also, since stages are independent in regards to the parameters, the women in the adult stage would not be able to reproduce (by the model’s assumptions) and thus their contribution to the infant group would fall to zero. The updated model of the population can be seen in equation 4.
Equation 4: Example StageBased Model With War
 Analyzing this matrix (Equation 4), the dominant eigenvalue was found to be 0.9834 meaning that under this war time scenario, the population growth would cease and the population would actually begin falling. Assuming it was a long war, the population distribution would also have changed. The new steady state distribution would be 3.5 % infants, 18.64 % children, 8.86 % adolescents, 8.61 % adults, and 60.37 % elderly adults. This is a crude example; however, it does show how population dynamics change when affected by positive/negative externalities. The dominant left eigenvalue states that the contribution distribution to fecundity will shift to [1.0000 1.9668 1.6731 1.3622 14.6104]. This suggest that the continuation of the population fell almost entirely on the elderly age group (those who were too old to fight in the war).
 The sensitivity of a system can be determined by examining a relationship between the sage contribution to fecundity and the stage population distribution. That is, the left and right eigenvectors. A matrix for the sensitivity can be set up by equation 5.
Equation 5: Sensitivity
 The elasticity of an element in the matrix can be determined by the product of the sensitivity with the models value scaled by the dominant eigenvalue as described by equation 5.
Equation 6: Elasticity
 For the example problem, the stability and elasticity matrices for the prewar case can be seen below. They show the model sensitivity to each entry of the matrix and the elasticity of the model to each individual parameter used in its construction. This allows for focusing efforts on the most significant parameters.
Equation 7: Sensitivity
Equation 8: Elasticity
 The plots of the elasticity for each type of parameter in the matrix model is shown below in figures 1 and 2. From these plots, the relative importance can be observed. From the plots, in the prewar case, the adult population survivability is the most elastic, whereas in the war case the elderly population survivability was the most elastic.
Figure 1: Elasticity: PreWar Case
Figure 2: Elasticity: War Case
SetUp
 One of the principle goals of this article was to determine what stage was most likely causing the downfall of the turtles. The concept was that if the most sensitive factor could be isolated and conservation resourced focused on helping that stage, the decline of the turtle population could potentially be stopped. The model utilized in the article broke the turtle population into seven groups: 1st year (eggs and hatchlings), small juveniles, large juveniles, subadults, novice breeders, 1styr remigrants, and mature breeders. A stagebased model was used instead of an age based model as the data typically grouped the animals by size and very little data was available about the age of the turtles measured. Using the values derived primarily from Richardson's 20 year Little Cucumber Island project for reproductive output, survival probabilities, and probabilities of remaining in the same stage, a 7 by 7 Lefkovitch Matrix was generated to model the turtle population shown in equation 9.
Equation 9: Lefkovitch Model of Loggerhead Population [1]
Analysis of Base Model
 Examination of the model revealed that the dominant right eigenvalue was 0.945. Since this is less than one, the population would decrease over time as suggested by the power method. Since the data showed that the population was diminishing towards extinction, it was not surprising that the dominant eigenvalue was less than one. To determine the population distribution over time, the eigenvector corresponding to the dominant right eigenvalue was determined and scaled so that the sum of all elements was one. This yielded a long term population distribution as described in table 1.
Table 1
1st Year

Small Juveniles

Large Juveniles

SubAdults

Novice Breeders

1styr Remigrants

Mature Breeders

20.651 % 
66.975 % 
11.460 % 
0.6623 % 
0.0363 % 
0.0005 % 
0.0029 %

 Examination of the dominant left eigenvalue describes the stage contribution to continued growth. Table 2 describes this as a factor of the contribution of the first stage.
Table 2
1st Year

Small Juveniles

Large Juveniles

SubAdults

Novice Breeders

1styr Remigrants

Mature Breeders

1.0000 
1.4007 
5.9955 
115.8445 
568.7809 
507.3730 
587.6693

 The response to a small initial perturbation from a steady state distribution was simulated using MATLAB to demonstrate the population response and distribution response. This can be seen in figures 3 and 4.
Figure 3: Population Response given Small Perturbation
Figure 4: Distribution Response given Small Perturbation
Sensitivity Analysis
 The benefit of a population matrix is that the growth, distribution, and fecundity contributions can be examined for specific parameters to determine where to optimally allocate conservation resources.
Parameter Variations: Fecundity & Survivability
Figure 5 : Variation of Fecundity [1]
 The parameters chosen for the model were based on a range of estimates. Therefore, to ensure that the model was valid, Crouse, Crowder, and Caswell simulated the model with a 50% decrease in fecundity for each stage and propagated the changes throughout the model. This caused the rate of population to decay to increase (decay faster) for all stages. The resulting change in the population decay was greatest for the juvenile and subadult categories. This can be seen graphically in the top figure, figure 5 taken from the Crouse, Crowder, and Caswell article.
 The model was also simulated for a perfect survivability case. This is not realistic, but it serves to demonstrate the effect survivability plays on the populations. Similar to the 50% decrease case, the greatest change was noted in the juvenile and subadult stages. It is interesting to note that even with 100% survivability of all eggs, the population, according to the model, would still decrease, suggesting that the conservation efforts at the time the article was written (almost solely egg/beach management) could not stop the population decay regardless of effort exerted.
 This property can be seen through sensitivity and elasticity analysis of the model. The sensitivity matrix for the population model is shown in equation 10 . The sensitivity matrix shows the relative importance of each parameter. To fully describe the relative response of changing a parameter, the elasticity matrix must be examined (equation 11). The elasticity matrix describes the response of the system relative to its dominant right eigenvalue with respect to a change in a parameter. The higher the elasticity, the greater the impact of a change in the particular matrix entry shown in equation.
Equation 10: Sensitivity Matrix
Equation 11: Elasticity Matrix
 The variations of all parameters, fecundity (F), survival with growth (G), and survival in the same class (P) for each stage was performed for the article and illustrated in figure 6 and reproduced in figure 7 to validate the plot in the paper. Elasticity or proportional sensitivity of the dominant eigenvalue is shown on the vertical axis. This figure demonstrates that fecundity has little effect on the population growth/decay, advancement to future states has a moderate effect, and that the most sensitive parameter was survival in the same stage. Once again, the juvenile, subadult, and mature adult populations demonstrated the most significant sensitivity to changing the growth/decay rate. The two figures shown are the original from the article as well as one generated to reproduce the results of the paper.
Figure 6: Elasticity Plot[1]
Figure 7: Reproduction of Elasticity Plot
Population Management Efforts
 The primary cause of death to turtles in the southeastern US is believed to be caused by fishing nets. Crouse, Crowder, and Caswell's article was written before Trawl Efficiency Devices or Turtle Exclusion Devices (TEDs) were mandated. The article suggest that implementation of these devices may significantly reduce the deaths of juvenile and adult turtles which, due to the high elasticity, contribute most significantly to the survivability of the species. Simulations were performed to demonstrate the significant impact these devices could play in improving survivability and generating growth of the population.
 Today, the population of loggerhead turtles continues to fall even after the implementation of the TEDs. Recent articles by the Fish and Wildlife Service [2] suggest that this decline may be at least partially caused by the following: the first few generations of TEDs, especially those designed for small boats, were ineffective or were disabled by fisherman who believed that the device would reduce their catch: the TEDs were designed primarily for shrimp boats and not mandated for all other net utilizing boats so many turtles were still caught in those types of nets: and finally, turtles are intentionally caught in places like Mexico to use their shells as souvenirs or decorations. Hopefully as technology improves and enforcement/implementation of TED like devices stepped up, the conservation efforts can halt the drop in sea turtle populations and allow them to replenish. [2]
Conclusions
 Stagebased population models are a useful tool for evaluating conservation and management efforts. From the findings, it is clear that exact values are not always required to determine the general trends in the population and to which parameters it is most sensitive to. The data used to generate the turtle model varied significantly between estimates; however, even factoring in significant variations 50+%, for some values, the model still showed that the most sensitive stages to the turtles' long term population growth was the juvenile, subadult, and mature adult stages. These findings were supported by the sensitivity analysis performed on the model. The article found that there was very little effect on the long term growth/decay rate from even 100% effective egg and beach management. This finding was supported by the data at the time showing that, although efforts had been underway to protect nest for over 20 years at the time the article was written, the population was still falling. The article suggested that more emphasis on management be placed on increasing the survivability of the stages that demonstrated the highest sensitivity, possibly through the use of TEDs. Although TEDs are now in use and the populations are still declining, the concept of using a sensitivity analysis of a stagebased population model is still valid. Even if it does not perfectly represent the population, it still gives planners and policy makers an idea of how to most effectively allocate resources and conservation efforts to maximize desired results.
External Links & References
 D. T. Crouse, L. B. Crowder, and H. Caswell. A stagebased population model for loggerhead sea turtles and implications for conservation. Ecology, 68(5):14121423, October 1987.
 National Marine Fisheries Service and U.S. Fish and Wildlife Service. 2008. Recovery Plan for the Northwest Atlantic Population of the Loggerhead Sea Turtle (Caretta caretta), Second Revision. National Marine Fisheries Service, Silver Spring, MD.
 "J. J. Enneson and J. D. Litzgus. Using longterm data and a stageclassified matrix to assess conservation strategies for an endangered turtle (Clemmys guttata). Biological Conservation, 141(2008):15601568, May 2008."
MATLAB Code
Function: Main
%% By: Luke Pederson
%% University of Colorado at Boulder
%% This code evaluates a matrix population model and determines the final
%% distribution, sensitivity, and elasticity as well as shows how the
%% system responds to an initial perturbation and the long term behavior
%% Note: The Tutle example was extracted from Crouse, Crowder, and
%% Caswell's 1987 article in Ecology "A stagebased population model for
%% loggerhead sea turtles and its implications on conservation"
%% Note 2: The function "func" is required and is shown at the bottom of
%% this code.
clc
clear all
close all
format compact
format long;
%% Model Choice
% People Example
% A=[0 0 0.05 0.05 0.05;
% 0.5 0.889 0 0 0;
% 0 .111 .75 0 0;
% 0 0 0.25 0.968 0;
% 0 0 0 0.032 0.98];
%People Example WAR
% A=[0 0 0.05 0 0.05;
% 0.5 0.889 0 0 0;
% 0 .111 .75 0 0;
% 0 0 0.25 0.726 0;
% 0 0 0 0.024 0.98];
%Paper Example (Turtles)
A=[0 0 0 0 127 4 80;
0.6747 0.7370 0 0 0 0 0;
0 0.0486 0.6610 0 0 0 0;
0 0 0.0147 0.6907 0 0 0;
0 0 0 0.0518 0 0 0;
0 0 0 0 0.8091 0 0;
0 0 0 0 0 0.8091 0.8089];
%% Eigen Decomposition, Distribution, Population, and Small Perturbation
[B pop Dist Evals Evect]=func(A) ;
%% Plots Section 1 : Distribution and Population
figure
hold on
plot(B(1,:),'c');
plot(B(2,:),'g');
plot(B(3,:),'y');
plot(B(4,:),'k');
plot(B(5,:),'r');
plot(B(6,:),'m');
plot(B(7,:),'b');
plot(pop,'bp');
legend('1st Year','Small Juveniles','Large Juveniles','SubAdults','Novice Breeders','1styr Remigrants','Mature Breeders','Total Population')
title('Population Progression with Small Initial Perturbation')
xlabel('Time')
ylabel('Population Size (Initially Scaled to 100)')
figure
hold on
plot(100*Dist(1,:),'c');
plot(100*Dist(2,:),'g');
plot(100*Dist(3,:),'y');
plot(100*Dist(4,:),'k');
plot(100*Dist(5,:),'r');
plot(100*Dist(6,:),'m');
plot(100*Dist(7,:),'b');
legend('1st Year','Small Juveniles','Large Juveniles','SubAdults','Novice Breeders','1styr Remigrants','Mature Breeders','Total Population')
title('Population Distributions with Small Initial Perturbation')
xlabel('Time')
ylabel('Population Size (Initially Scaled to 100)')
%% Sensitivity analysis
%find left eigenvalues/eigenvectors
[Levec Leval]=eig(A');
clear Leval
Leval=eig(A');
indi=find(abs(Leval)==max(abs(Leval)));
v=Levec(:,indi);
v=v/v(1);
%find right eigenvalues/eigenvectors
indi=find(abs(Evals)==max(abs(Evals)));
w=Evect(:,indi);
%initialize Matrices for Speed
sen=zeros(length(A),length(A));
El=zeros(length(A),length(A));
% Determine Sensitivity and Elasticity
for i=1:length(A)
for j=1:length(A)
sen(i,j)=v(i)*w(j)/(w'*v);
El(i,j)=A(i,j)/max(Evals)*sen(i,j);
end
end
%Determine elasticity per element type, Fecundity, survivability, etc.
F=El(1,:);
P(1)=El(1,1);
for i=1:length(El)1
P(i+1)=El(i+1,i+1);
G(i)=El(1+i,i);
end
G(length(El))=0;
%% Plot Section 2 : Elasticity by stage by parameter type
figure
hold on
plot(1:length(El),F,'r')
plot(1:length(El),P,'g')
plot(1:length(El),G,'b')
legend('Fecundity','Probability of Remaining in Stage','Probability of Advancing Stages')
title('Elasticity')
ylabel('Elasticity')
xlabel('Stage')
Function: func
function [B pop Dist Evals Evect]=func(A)
% This code is required to run the main program above
[Evect dummy]=eig(A);
clear dummy
Evals=eig(A);
B(:,1)=abs(Evect(:,length(A)));
pop(1)=sum(B(:,1));
t=rand(1,length(A));
B(:,1)=B(:,1)./pop(1).*100.*(2transpose(t))./2;
Dist(:,1)=B(:,1)/sum(B(:,1));
pop(1)=sum(B(:,1));
B(:,1)=B(:,1)/pop(1)*100;
pop(1)=sum(B(:,1));
for i=2:25
B(:,i)=A*B(:,i1);
pop(i)=sum(B(:,i));
Dist(:,i)=B(:,i)/pop(i);
end
