%    Markov Chain model for landscape succession
%            from Ecobeaker - Math 151 
% This is to illustrate eigenvalues and 
% eigenvectors arising from a simple community
% succession model from Ecobeaker in which there are 
% 3 species (Bare soil, Grass and Shrub) with
% several different possible transition matrices used
%   The program calculates for each matrix the dominant 
% eigenvalue and shows it is  = 1 and then finds the 
% eigenvector associated with this. 
% P is the projection matrix and n0 is the
% initial landscape structure, ni is the population
% structure at time I which you could think of as the 
% total area in each vegetation type at each time, 
% nitot is the total area of the region
% at time i, ninorm is the normalized landscape structure
% at time i (sum of the elements of ninorm is 1)
P=[.8 0 .01; .2 .5 .49; 0 .5 .5 ],pause
n0=[100;0;0],pause
n0tot=sum(n0)
n1=P*n0,pause
n1tot=sum(n1),pause
n1norm=n1/n1tot,pause
% Now compute the landscape structure at time 100
n100=P^100*n0
n100tot=sum(n100)
n100norm=n100/n100tot,pause
% Now compute the eigenvalue and the eigenvectors and compare these
% to the landscape structure at time 100
[X,e]=eig(P),pause
% Now get the eigenvector for the dominant eigenvalue, 
% normalize it to get ev1norm, and compare this to 
% the normalized landscape structure at time 21
ev1=X(:,2),pause
ev1norm=ev1/sum(ev1),pause
n100norm,pause
%
% Now repeat for a different matrix
%
P=[.8 0 .2; .2 .5 .3; 0 .5 .5 ],pause
n0=[100;0;0],pause
n0tot=sum(n0)
n1=P*n0,pause
n1tot=sum(n1),pause
n1norm=n1/n1tot,pause
% Now compute the landscape structure at time 100
n100=P^100*n0
n100tot=sum(n100)
n100norm=n100/n100tot,pause
% Now compute the eigenvalue and the eigenvectors and compare these
% to the landscape structure at time 100
[X,e]=eig(P),pause
% Now get the eigenvector for the dominant eigenvalue, 
% normalize it to get ev1norm, and compare this to 
% the normalized landscape structure at time 21
ev1=X(:,2),pause
ev1norm=ev1/sum(ev1),pause
n100norm,pause
%
% Now repeat for a different matrix
%
P=[.8 0 .1; .2 .4 .3; 0 .5 .5 ],pause
n0=[100;0;0],pause
n0tot=sum(n0)
n1=P*n0,pause
n1tot=sum(n1),pause
n1norm=n1/n1tot,pause
% Now compute the landscape structure at time 100
n100=P^100*n0
n100tot=sum(n100)
n100norm=n100/n100tot,pause
% Now compute the eigenvalue and the eigenvectors and compare these
% to the landscape structure at time 100
[X,e]=eig(P),pause
% Now get the eigenvector for the dominant eigenvalue, 
% normalize it to get ev1norm, and compare this to 
% the normalized landscape structure at time 21
ev1=X(:,2),pause
ev1norm=ev1/sum(ev1),pause
n100norm,pause