chp4_4.m

来自「置信区间的模拟程序」· M 代码 · 共 60 行

%The Matlab code for Chp2.15
clear ; 
clc; 
close;

rand('seed',sum(100*clock));

%set params
p = 0.4;
nsamples = 1000;
num_exps = 100;
alpha = 0.05;
n_hits = zeros(nsamples,1);     %置信区间包含p的次数
interval_lengths = zeros(nsamples, 1);

for n=1 : nsamples
    %1-α置信区间
    e = log(2/alpha);
    e = e/2/n;
    e = sqrt(e);
    
    n_hits(n, 1) = 0;
    for k = 1: num_exps         %random generate samples
        tmp = rand(n, 1);
        for i= 1: n
            if(tmp(i, 1) <= p) 
                X_n(i) = 1;
            else
                X_n(i) = 0;
            end
        end
   
        p_cap = mean(X_n);
     
        %confidence interval 
        C_lower = p_cap - e;
        C_lower = max(C_lower, 0);
        C_upper = p_cap + e;
        C_upper = min(C_upper, 1);
        
        interval_lengths(n, 1) = C_upper - C_lower;
    
        if(p > C_lower && p < C_upper)   % in the confidence interval
            n_hits(n, 1) = n_hits(n, 1) + 1;
        end
    end
end
 



%plot the results of the last expriment
figure(1);
n_hits = n_hits/num_exps;
plot(n_hits);

figure(2);
plot(interval_lengths);