function [ci2] = bino_wald(n, x, aa)

%   Wald 2-sided interval of binomial

p_h = x./n;      za = norminv(1-aa/2,0,1);
se = za*sqrt(x.*(n-x)/n^3);
p1 = p_h - se;          p2 = p_h + se;

if sum(p1 < 0) >=1
    for k = 1 : length(x)
        if p1(k,1) < 0
            p1(k,1) = 0;
        else
        end
    end
else if sum(p2 > 1) >= 1
        for k = 1 : length(x)
            if p1(k,1) > 1
                p1(k,1) = 1;
            else
            end
        end
    else
    end
end
ci2 = [p1 p2];