%matplotlib inline px = linspace(0,1,50) plot( px, px**5) xlabel(r'$\theta$',fontsize=18) ylabel(r'$\beta$',fontsize=18) from sympy.abc import p,k # get some variable symbols import sympy as S expr=S.Sum(S.binomial(5,k)*p**(k)*(1-p)**(5-k),(k,3,5)).doit() p0=S.plot(p**5,(p,0,1),xlabel='p',show=False,line_color='b') p1=S.plot(expr,(p,0,1),xlabel='p',show=False,line_color='r',legend=True,ylim=(0,1.5)) p1.append(p0[0]) p1.show() # implement 1st test where all five must be heads to declare underlying parameter > 0.5 from scipy import stats b = stats.bernoulli(0.8) # true parameter is 0.8 (i.e. hypothesis H_1) samples = b.rvs(1000).reshape(-1,5) # -1 means let numpy figure out the other dimension i.e. (200,5) print 'prob of detection = %0.3f'%mean(samples.sum(axis=1)==5) # approx 0.8**3 # here's the false alarm case b = stats.bernoulli(0.3) # true parameter is 0.3 (i.e. hypothesis H_0) samples = b.rvs(1000).reshape(-1,5) print 'prob of false alarm = %0.3f'%mean(samples.sum(axis=1)==5) # implement majority vote test where three of five must be heads to declare underlying parameter > 0.5 b = stats.bernoulli(0.8) # true parameter is 0.8 (i.e. hypothesis H_1) samples = b.rvs(1000).reshape(-1,5) print 'prob of detection = %0.3f'%mean(samples.sum(axis=1)>=3) # here's the false alarm case b = stats.bernoulli(0.3) # true parameter is 0.3 which means it's hypothesis H_0 samples = b.rvs(1000).reshape(-1,5) print 'prob of false alarm = %0.3f'%mean(samples.sum(axis=1)>=3)