%Stochastic Competition
clear all 
a10=2; a11=0.03; a12=0.02;
a20=1.5; a21=0.01; a22=0.04;
t(1)=0; x(1)=50; y(1)=25; 
i=1; 
while (x(i)>0 && y(i)>0 && t(i)<5)
   b1=a10*x(i); b2=a20*y(i);
   d1=x(i)*(a11*x(i)+a12*y(i));
   d2=y(i)*(a21*x(i)+a22*y(i));
   tot=b1+b2+d1+d2;
   u1=rand; u2=rand;
   
   t(i+1)=-log(u1)/tot+t(i);
   x(i+1)=x(i); y(i+1)=y(i);
   if (u2<=b1/tot)
       x(i+1)=x(i)+1;
   elseif (u2>b1/tot && u2<=(b1+b2)/tot)
       y(i+1)=y(i)+1;
   elseif (u2>(b1+b2)/tot && u2<=(b1+b2+d1)/tot)
       x(i+1)=x(i)-1;
   else
       y(i+1)=y(i)-1;
   end
   i=i+1;
end

xx=linspace(0,100,11);
isox1=a10/a12-a11/a12*xx;
isox2=a20/a22-a21/a22*xx;
l1=stairs(x,y,'r'); set(l1, 'Linewidth', 2);
hold on
plot(xx, isox1, 'k', xx, isox2, 'k', 'Linewidth', 2);
xlabel('X'); ylabel('Y');

figure
plot(t,x, t, y)
xlabel('time'); ylabel('Population Size');
legend('X','Y')