clear all
n = 10;     n1 = 10000;      ap = 0.1;
pp = 0.5;       x = binornd(n, pp, n1, 1);

w1 = bino_wilson(n, x, ap);
dcs1 = sum((pp <= w1)')';   rst1 = sum(dcs1 == 1)/n1;

pp = 0.9;       x = binornd(n, pp, n1, 1);
w2 = bino_wilson(n, x, ap);
dcs2 = sum((pp <= w2)')';   rst2 = sum(dcs2 == 1)/n1;


fprintf(1, 'The proportions are %5.4f at p=0.5 and %5.4f at p=0.9\n',...
    rst1, rst2);