femalemort=0.01:0.01:0.99; % extra mortality alpha=1.8; for k=1:1000 k for i=1:length(femalemort) alphamales=2*0.5*10000; alphamalesleft=alphamales; betamales=2*0.5*10000; betamalesleft=betamales; females=round(0.5*4*10000*(1-femalemort(i))); if females0 if rand(1)<(alpha*alphamalesleft)/(alpha*alphamalesleft+betamalesleft) alphamalesleft=alphamalesleft-1; else betamalesleft=betamalesleft-1; end females=females-1; end rvalpha(i,k)=1-alphamalesleft/alphamales; % one minus probability of remaining unmated rvbeta(i,k)=1-betamalesleft/betamales; else % different story if more females than males femalematingprob=1-(females-alphamales-betamales)/females; rvalpha(i,k)=1/femalematingprob; rvbeta(i,k)=1/femalematingprob; end end rva=mean(rvalpha,2); rvb=mean(rvbeta,2); figure(1); plot(femalemort,rva./rvb,'.'); drawnow; end %% figure(1); plot(1-femalemort,rva20./rvb20,'k-', 'LineWidth', 1.5); hold on; plot(1-femalemort,rva50./rvb50,'k:', 'LineWidth', 1.5); hold on; plot(1-femalemort,rva80./rvb80,'k--', 'LineWidth', 1.5); hold on; axis('square') axis([0 1 1 1.9]); yticks([1:0.2:1.9]); xlabel('Female survival {\it s}') ylabel('\itv_1\rm / \itv\rm_+') set(gca, 'Box','on','FontSize', 15) legend('\alpha = 0.2', '\alpha = 0.5', '\alpha = 0.8') text(0+0.04*(1-0), 1.9-0.04*(1.9-1), 'A', 'FontSize', 15, 'FontWeight','Bold'); print -depsc RVearlyLater