% PROGRAM iidenv
% 		Simulates growth of a structured population in an iid
%		stochastic environment

%************* USER-SPECIFIED PARAMETERS ***************
% Enter name of m file containing matrices; e.g., the file  
% hudmats.m contains 4 matrices (A85, A86, A87, A88),
% one for each year of Frost's demographic study
hudmats;						
% Change names of matrices below to correspond with those 
% in the file named in the preceding command:
matrices=[A85(:) A86(:) A87(:) A88(:)];
% Enter probabilities of choosing each matrix below;
% NOTE: sum of penv must equal 1
penv=[.25 .25 .25 .25];	
% Enter time to predict future population size
tmax=50;						
% Enter number of realization of population growth to simulate
numreps=5000;				
% Enter initial population vector
n0=[4264; 3; 30; 16; 25; 5];				
%********************************************************

rand('state',sum(100*clock));		% seed random number generator
s=sqrt(size(matrices,1));			% s=number of classes
cumdist=cumsum(penv);				% CDF for environmental states

Nend=[];
for i=1:numreps                     % For each realization,
   n=n0;                            % starting with the initial pop. vector,
   for t=1:tmax                     % for each time step,
      x=sum(rand>=cumdist)+1;		% draw a matrix at random,
      A=reshape(matrices(:,x),s,s);	% extract the chosen matrix 
                                    %   from "matrices" and
      n=A*n;						% project the population 1 year ahead
   end;
   Nend=[Nend sum(n)];				% store population size at tmax
end;

hist(Nend,50)			% create histogram of population size at tmax
MeanN=mean(Nend)		% mean population size at tmax
MedianN=median(Nend)	% median population size at tmax