diff options
Diffstat (limited to 'sourcecodes/bnt-master/KPMstats/dirichletrnd.m')
| -rw-r--r-- | sourcecodes/bnt-master/KPMstats/dirichletrnd.m | 32 |
1 files changed, 32 insertions, 0 deletions
diff --git a/sourcecodes/bnt-master/KPMstats/dirichletrnd.m b/sourcecodes/bnt-master/KPMstats/dirichletrnd.m new file mode 100644 index 00000000..f11074db --- /dev/null +++ b/sourcecodes/bnt-master/KPMstats/dirichletrnd.m @@ -0,0 +1,32 @@ +function x = dirichletrnd(alpha) +%DIRICHLETRND Random vector from a dirichlet distribution. +% x = dirichletrnd(alpha) returns a vector randomly selected +% from the Dirichlet distribution with parameter vector alpha. +% +% The algorithm used is the following: +% For each alpha(i), generate a value s(i) with distribution +% Gamma(alpha(i),1). Now x(i) = s(i) / sum_j s(j). +% +% The above algorithm was recounted to me by Radford Neal, but +% a reference would be appreciated... +% Do the gamma parameters always have to be 1? +% +% Author: David Ross +% $Id: dirichletrnd.m,v 1.1.1.1 2005/05/22 23:32:12 yozhik Exp $ + +%------------------------------------------------- +% Check the input +%------------------------------------------------- +error(nargchk(1,1,nargin)); + +if min(size(alpha)) ~= 1 | length(alpha) < 2 + error('alpha must be a vector of length at least 2'); +end + + +%------------------------------------------------- +% Main +%------------------------------------------------- +gamma_vals = gamrnd(alpha, ones(size(alpha)), size(alpha)); +denom = sum(gamma_vals); +x = gamma_vals / denom; \ No newline at end of file |