After using the code below to create a series of .png images, from the command line, use, e.g., ImageMagick convert to create an animated gif from the .pngs

convert -delay 20 -loop 0 img*.png bsim.gif
In [1]:
from scipy.misc import comb
r=exp(log(2)/10)
rest=map(int,32*cumprod(r*ones(60))+.5)  #log-spaced steps from 32 to 2048
In [3]:
trials=10000

f=2  # for p=1/2
#f=6 #for p=1/6

p=1./f
sigma=p*(1-p)

for k,n in enumerate(range(1,33)+rest):
    figure(figsize=(6,4))
    
    hist([sum(np.random.rand(n) < p) for i in xrange(trials)],bins=range(0,n+2),
               linewidth=0,align='left')
#simulation done entirely in the above line, rest just adds text annotations
    
    if n<=32:
       m=arange(n+1)
       ms=trials*max(p**m*(1-p)**(n-m)*comb(n,m))
    else: 
       ms = trials/sqrt(2*pi*sigma*n)
#ms is used to set the upper ylim, so that it decreases smoothly with n
        
    text(0.7, 0.5,'n={}'.format(n), horizontalalignment='left',
             verticalalignment='center',fontsize=24,
              transform=gca().transAxes)  #write the n=...
    text(0.7, 0.3,'$\sigma/\mu={}$'.format(round(sqrt((1-p)/(p*n)),2)), horizontalalignment='left',
             verticalalignment='center',fontsize=18,
              transform=gca().transAxes)  #write the sigma/mu=...
    annotate('n/{}'.format(f),(p*n,0),(p*n,-ms/6.4),arrowprops=dict(arrowstyle="->"),ha='center',fontsize=20,color='b')

    xlim(-.5,n+.5*(1-p)/p)
    ylim(0,1.1*ms)

    savefig('pngs/img{:03d}.png'.format(k))  #save 92 images in pngs subdir
    close()
In [ ]: