function W = nchoosecrit(S, FUN)
% NCHOOSECRIT - subsets (combinations of elements of a set) that fulfill a criterion
%
% W = NCHOOSECRIT(S, FUN) returns those combinations of one or more elements
% of the set S (called a subset) that fulfill a specific criterion.
% This criterion is specified by the function FUN. FUN is a function
% handle to a function that takes one input argument and returns a
% logical scalar value.
% W will be cell array of row vectors. Each cell of W holds one of the
% combinations C of S for which FH(C) is true.
%
% In total, there are 2^N-1 possible subsets of S (N being the number of
% elements of S). This number grows rapidly with increasing N. W is a
% selection of those subsets.
%
% S can be a cell array, and each cell of W will then contain a cell array.
%
% Examples:
% % find the subsets that sum op to 6
% nchoosecrit([1 2 3 4 5 6], @(x) sum(x)==6)
% % -> { [1 2 3], [2 4], [1 5], [6]}
%
% % find subgroups of 3 or more people that contain either James or Bob,
% % but not both!
% S = {'Bob' 'Tom' 'Joe' 'Bill' 'James'} ; % the whole group
% % criterion 1:
% fh1 = @(x) numel(x) >= 3 ;
% % criterion 2
% fhname = @(x,y) any(strncmp(y,x,numel(y))) ;
% fh2 = @(x) xor(fhname(x,'James'), fhname(x,'Bob')) ;
% % the 2 criterions combined:
% fhcomb = @(x) fh1(x) && fh2(x) ;
% W = nchoosecrit(S, fhcomb)
%
% Notes:
% - If S contain non-unique elements (e.g. S = [1 1 2]), NCHOOSECRIT will
% return non-unique cells. In other words, NCHOOSECRIT treats all elements
% of S as being unique. One could use NCHOOSECRIT(UNIQUE(S)) to avoid that.
% - The output is the same as
% Wtemp = nchoose(S) ; W = Wtemp(cellfun(Wtemp, fh)) ;
% but does not create the (possible large) temporary array Wtemp.
%
% See also NCHOOSEK, PERMS
% NCHOOSE, COMBN, ALLCOMB on the File Exchange
% for Matlab R13 and up
% version 1.0 (feb 2013)
% (c) Jos van der Geest
% email: jos@jasen.nl
% History
% 1.0, feb 2013 - inspired by a post on Matlab Answers
error(nargchk(nargin,2,2)) ;
N = numel(S) ;
if N == 0,
% nothing to do
W = {} ;
return ;
end
try
% these two tries catch many errors
X1 = FUN(S(1)) ; % try it on a single element
X2 = FUN(S) ; % try it on the whole set
catch ME
disp('This seems to be an invalid criterion function') ;
rethrow(ME)
end
if numel(X1) ~= 1 || numel(X2) ~= 1 || ~islogical(X1) || ~islogical(X2)
error('The criterion function should return a logical scalar') ;
end
M = (2^N)-1 ; % The total number of possible combinations
% The selection of elements is based on the binary representation of all
% numbers X between 1 and M. This binary representation is is retrieved by
% the formula: bitget(X * (2.^(N-1:-1:0)),N) > 0
% See NCHOOSE (available on the File Exchange) for details
p2=2.^(N-1:-1:0) ; % This part of the formula can be taken out of the loop
% We pre-allocate the output, but only to a certain extend
W = cell(2^16,1) ;
wk = 0 ; % We'll add the subsets to the list one by one
% loop over all subsets
for k=1:M,
% calculate the (reversed) binary representation of k
% select the elements of the set based on this representation
subset = S(bitget(k*p2,N) > 0) ;
if FUN(subset),
% does it fullfill the criterion?
wk = wk + 1 ; % go to the next element in W
W{wk} = subset ; % add to the list
end
end
W = W(1:wk) ; % remove the unused pre-allocated elements of W