camir-aes2014: toolboxes/FullBNT-1.0.7/KPMstats/dirichletpdf.m annotate

first hg version after svn

rev	line source
wolffd@0	1 function p = dirichletpdf(x, alpha)
wolffd@0	2 %DIRICHLETPDF Dirichlet probability density function.
wolffd@0	3 % p = dirichletpdf(x, alpha) returns the probability of vector
wolffd@0	4 % x under the Dirichlet distribution with parameter vector
wolffd@0	5 % alpha.
wolffd@0	6 %
wolffd@0	7 % Author: David Ross
wolffd@0	8
wolffd@0	9 %-------------------------------------------------
wolffd@0	10 % Check the input
wolffd@0	11 %-------------------------------------------------
wolffd@0	12 error(nargchk(2,2,nargin));
wolffd@0	13
wolffd@0	14 % enusre alpha is a vector
wolffd@0	15 if min(size(alpha)) ~= 1 \| ndims(alpha) > 2 \| length(alpha) == 1
wolffd@0	16 error('alpha must be a vector');
wolffd@0	17 end
wolffd@0	18
wolffd@0	19 % ensure x is is a vector of the same size as alpha
wolffd@0	20 if any(size(x) ~= size(alpha))
wolffd@0	21 error('x and alpha must be the same size');
wolffd@0	22 end
wolffd@0	23
wolffd@0	24
wolffd@0	25 %-------------------------------------------------
wolffd@0	26 % Main
wolffd@0	27 %-------------------------------------------------
wolffd@0	28 if any(x < 0)
wolffd@0	29 p = 0;
wolffd@0	30 elseif sum(x) ~= 1
wolffd@0	31 disp(['dirichletpdf warning: sum(x)~=1, but this may be ' ...
wolffd@0	32 'due to numerical issues']);
wolffd@0	33 p = 0;
wolffd@0	34 else
wolffd@0	35 z = gammaln(sum(alpha)) - sum(gammaln(alpha));
wolffd@0	36 z = exp(z);
wolffd@0	37
wolffd@0	38 p = z * prod(x.^(alpha-1));
wolffd@0	39 end

Mercurial > hg > camir-aes2014