This repository has been archived on 2021-05-17. You can view files and clone it, but cannot push or open issues or pull requests.
scientificComputing/statistics/exercises/centrallimit.m

76 lines
2.0 KiB
Matlab

%% central limit theorem
n = 10000;
m = 10; % number of loops
%% (b) a single data set of random numbers:
x = rand(n, 1);
%% (c) plot probability density:
%histogram(x, 'Normalization', 'pdf');
[h,b] = hist(x, 20);
h = h/sum(h)/(b(2)-b(1)); % normalization
bar(b, h)
title('A uniform distribution')
xlabel('x')
ylabel('probability density')
pause(2.0)
%% (d) sum of two random numbers:
y = rand(n, 1);
x = x + y;
%histogram(x, 'Normalization', 'pdf');
[h,b] = hist(x, 20);
h = h/sum(h)/(b(2)-b(1)); % normalization
bar(b, h)
title('Sum of two uniform distributions')
xlabel('x')
ylabel('probability density')
pause(2.0)
%% (f) sum up more distributions:
x = zeros(n, 1);
means = zeros(m, 1);
stds = zeros(m, 1);
for i=1:m
y = rand(n, 1); % new uniform distributed numbers
x = x + y; % add them to the sum
mu = mean(x); % compute mean
sd = std(x); % compute standard deviation
means(i) = mu; % store mean
stds(i) = sd; % store standard deviation
%xx = min(x):0.01:max(x);
xx = -1:0.01:i+1; % x-axis values for plot of pdf
p = exp(-0.5*(xx-mu).^2/sd^2)/sqrt(2*pi*sd^2); % pdf
plot(xx, p, 'r', 'linewidth', 3 )
ns = sprintf('N=%d', i);
text(0.1, 0.9, ns, 'units', 'normalized')
hold on
%histogram(x, 20, 'Normalization', 'pdf');
[h,b] = hist(x, 20);
h = h/sum(h)/(b(2)-b(1)); % normalization
bar(b, h)
hold off
xlim([-0.5, i+0.5])
xlabel('x')
ylabel('summed pdf')
savefigpdf(gcf, sprintf('centrallimit-hist%02d.pdf', i), 6, 5);
if i < 6
pause(3.0)
end
end
%% (h) mean and standard deviation in dependence on number of summed up distributions:
nx = 1:m;
plot(nx, means, 'b', 'linewidth', 4);
hold on
plot(nx, stds, 'r', 'linewidth', 4);
xx = 0:0.01:m;
sdu = 1.0/sqrt(12); % standarad deviation of the uniform distribution
plot( xx, sqrt(xx)*sdu, 'k' )
legend('mean', 'std', 'theory')
xlabel('N')
hold off
savefigpdf(gcf, 'centrallimit-samples.pdf', 6, 5);