
% selectio.m
% This program, selectio.m, simulates selection
% in a one locus two allele model, with random
% mating, infinite population size, no mutation
% You are prompted to read in three selection
% coefficients: s1 for AA, s2 for Aa, and 
% s3 for aa - these represent the conditional
% probability of surviving to reproduce for 
% each of these three genotypes.
% You also are asked to give the time interval
% over which the run is to be made.
!c:
q=1;
while q>0;
s1=input('selection coefficient for AA: s1    ')
s2=input('selection coefficient for Aa: s2    ')
s3=input('selection coefficient for aa: s3    ')
p0=input('initial frequency of A allele:     ')
n=input('end of time interval b:     ')
p=zeros(n+1,1);
q=zeros(n+1,1);
t=zeros(n+1,1);
p(1)=p0;
q(1)=1-p0;
for i=1:n
t(i)=i-1;
a=s1*p(i)^2+2*s2*p(i)*q(i)+s3*q(i)^2;
p(i+1)=(s1*p(i)^2+s2*p(i)*q(i))/a
q(i+1)=(s2*p(i)*q(i)+s3*q(i)^2)/a;
end
t(n+1)=n;
plot(t,p,t,p,'o'),pause
q=input('Do you want to stop - if so enter 0')
end

