% Fluctuation of histogram bin counts over repeated experiments
% Textbook: Part I, "Statistical Plots", the variance of a bin count
%
% The count n_i in bin i is binomial, n_i ~ Bin(N, p_i), so repeating the
% experiment redraws every bar within about +-sqrt(N p_i (1 - p_i)) of its
% mean N p_i. Central bars look stable; tail bars swing by a large fraction
% of their own height (relative error ~ 1/sqrt(n_i)).

clearvars; close all
rng(0)

R = 30;         % repetitions of the experiment
N = 500;        % observations per experiment
k = 15;         % bins over [-3, 3]
edges = linspace(-3, 3, k + 1);
centers = edges(1:end-1) + diff(edges)/2;

Phi = @(x) erfc(-x/sqrt(2))/2;              % N(0,1) CDF
p = Phi(edges(2:end)) - Phi(edges(1:end-1));% bin probabilities
mean_th = N*p;
std_th = sqrt(N*p.*(1 - p));

counts = zeros(R, k);
figure
hold on
for r = 1:R
    counts(r, :) = histcounts(randn(1, N), edges);
    stairs([edges(1:end-1) edges(end)], [counts(r, :) counts(r, end)], 'Color', [0.7 0.7 0.7])
end
errorbar(centers, mean_th, std_th, 'ko', 'LineWidth', 1.2, 'MarkerFaceColor', 'k')
hold off
xlabel('x')
ylabel('n_i')
title(sprintf('%d repetitions of N = %d, theory N p_i \\pm sqrt(N p_i(1-p_i))', R, N))
grid on

fprintf('%8s %10s %10s %10s %10s\n', 'center', 'mean, th', 'mean, sim', 'std, th', 'std, sim')
fprintf('%8.2f %10.2f %10.2f %10.2f %10.2f\n', [centers; mean_th; mean(counts); std_th; std(counts)])