clear all; close all;
y = csvread('mm1.out');

% --- simulacao finita ---

% amostra fixa
n = 10;
beta = .02;
alpha = .05;

x = y(1:n,end);
while tinv(1-alpha/2,n-1)*sqrt(var(x)/n) > beta;
    n = n + 1;
end;
x = y(1:n,end);
disp(sprintf('amostra fixa:\nn = %d\nmedia = %f\nerro absoluto = %f', n, mean(x), tinv(1-alpha/2,n-1)*sqrt(var(x)/n)));
disp(sprintf('intervalo de confianca: [%f, %f]', mean(x) - tinv(1-alpha/2,n-1)*sqrt(var(x)/n), mean(x) + tinv(1-alpha/2,n-1)*sqrt(var(x)/n)));

% procedimento sequencial
n = 10;
gama = .05;
disp(sprintf('\nsequencial:\ngama = %f\ngama/(1-gama) = %f', gama, gama / (1 - gama)));

x = y(1:n,end);
u = mean(x); s2 = var(x); delta = tinv(1-alpha/2,n-1)*sqrt(var(x)/n);
while delta / abs(u) > gama;
    n = n + 1;
    x = y(1:n,end);
    u = mean(x); s2 = var(x); delta = tinv(1-alpha,n-1)*sqrt(var(x)/n);
end;
x = y(1:n,end);
disp(sprintf('n = %d\nmedia = %f\nerro relativo = %f', n, mean(x), tinv(1-alpha,n-1)*sqrt(var(x)/n) / abs(mean(x))));

% --- simulacao infinita ---

n = 5; m = 60;
x = y(1:n,1:m);
plot(mean(x));                                                axis([1 m .3 .6]); grid on;
k = 5; v = tsmovavg(mean(x),'s',k); figure; plot(v(floor(k/2):end)); axis([1 m .3 .6]); grid on;

% replicacao/remocao
l = 15;
n = 5; x = y(1:n,l:end); x = mean(x');
disp(sprintf('intervalo de confianca: [%f, %f], n = %d', mean(x) - tinv(1-alpha/2,n-1)*sqrt(var(x)/n), mean(x) + tinv(1-alpha/2,n-1)*sqrt(var(x)/n), n));

n = 30; x = y(1:n,l:end); x = mean(x');
disp(sprintf('intervalo de confianca: [%f, %f], n = %d', mean(x) - tinv(1-alpha/2,n-1)*sqrt(var(x)/n), mean(x) + tinv(1-alpha/2,n-1)*sqrt(var(x)/n), n));

% estimando proporcoes
z = x > .5; % mais de 30 segundos
disp(z);
p = sum(z)/n;
disp(sprintf('p = %f', p));
disp(sprintf('intervalo de confianca: [%f, %f]', p - norminv(1-alpha/2)*sqrt(p*(1-p)/n), p + norminv(1-alpha/2)*sqrt(p*(1-p)/n)));