lbincoeff

(tmp_for_tar/misc/lbincoeff.m)

Function Synopsis

l = lbincoeff (n,k)

Help text

 l = lbincoeff (n,k)          - Natural logarithm of binomial coefficient

           / /   \ \
           | | n | |        /  n (n-1) (n-2) ... (n-k+1) \
 l =   log | |   | |  = log | -------------------------  |
           | | k | |        \      k!                    /
           \ \   / /

 Last modified: April 2001

Cross-Reference Information

test_lbincoeff tmp_for_tar/misc/test_lbincoeff.m

Listing of function file tmp_for_tar/misc/lbincoeff.m

## l = lbincoeff (n,k)          - Natural logarithm of binomial coefficient
## 
##
##           / /   \ \
##           | | n | |        /  n (n-1) (n-2) ... (n-k+1) \
## l =   log | |   | |  = log | -------------------------  |
##           | | k | |        \      k!                    /
##           \ \   / /
##

## Author:        Etienne Grossmann  <etienne@isr.ist.utl.pt>
## Last modified: April 2001

function l = lbincoeff (n,k)

sn = size (n); psn = prod (sn);
sk = size (k); psk = prod (sk);

if psn != psk && all (max ([sn;sk]')>1)
  printf ("lbincoeff : Args have size %ix%i and %ix%i\n",\
	  sn,sk);
  keyboard
end


if psn != psk
  if     all (sn == 1), n = n * ones (sk);
  elseif all (sk == 1), k = k * ones (sn); 
  end
elseif any (sn != sk)
  k = reshape (k, sn);
end

l = zeros (size (n));


l(find (k>n  | k<0)) = -inf;
l(find (k==n | k==0)) = 0;

ii = find (k<n & k>0);
ii = ii(:)';

if isempty (ii), return; end

N = max (n);
c = ones (N,1);
l0 = log(1:N);
				# Avoid doing a N x length(n) matrix ... at
				# the cost of a loop
for i = ii
  c(1:n(i)) = 1;
  c(n(i)+1:N) = 0;
  c(1:k(i))--;
  c(1:n(i)-k(i))--;
  l(i) = l0*c;
end