% L5Bdirichlet   Demo of Dirichlet distribution

% Counts of observations.
counts = [5 3 2];

% Hint: Try bigger values, e.g. [50 30 20] for a 100-person sample
% with those counts of parties 1, 2, 3.
% The posterior mode will be the same (0.5, 0.3) but the
% credible region will be much smaller.

N    = 100;     % Image resolution
grid = linspace(0,1,N);
[p,q] = meshgrid(grid,grid);
r = 1-p-q;
P = [p(:) q(:) r(:)];

% Prior
subplot(1,3,1)
A = [1 1 1];    % parameters for uniform Dirichlet distribution
fprior = dirpdf(P, A);
fprior = reshape(fprior, N, N);
surf(p, q, fprior);
colorbar
title('Prior density')
xlabel('p')
ylabel('q')

% Posterior
subplot(1,3,2)
B = 1+counts;
fpost = dirpdf(P, B);
fpost = reshape(fpost, N, N);
surf(p, q, fpost);
colorbar
title('Posterior density');
xlabel('p')
ylabel('q')


% Posterior 50%, 95% and 99% credible regions
subplot(1,3,3)
I = credibleregion(fpost, 0.50);
J = credibleregion(fpost, 0.95);
K = credibleregion(fpost, 0.99);
imagesc(grid, grid, I+J+K)
set(gca,'ydir','normal')
title('50%/95%/99% credible regions');
xlabel('p')
ylabel('q')