function [y]=choose(n,k)
% binomial coefficient, n chose k
%
y=exp(gammaln(n+1)-gammaln(k+1)-gammaln(n-k+1)); 
%