%FEATSELLR Plus-L-takeaway-R feature selection for classification
%
% [W,RES] = FEATSELLR(A,CRIT,K,L,R,T)
%
% INPUT
% A Dataset
% CRIT String name of the criterion or untrained mapping
% (optional; default: 'NN', i.e. 1-Nearest Neighbor error)
% K Number of features to select
% (optional; default: return optimally ordered set of all features)
% L Number of features to select at a time (plus-L, default: 1), L ~= R
% R Number of features to deselect at a time (takeaway-R, default: 0)
% T Tuning set (optional)
% N Number of cross-validations (optional)
%
% OUTPUT
% W Output feature selection mapping
% RES Matrix with step-by-step results of the selection
%
% DESCRIPTION
% Floating selection of K features using the dataset A, by iteratively
% selecting L optimal features and deselecting R. Starts from the full
% set of features when L < R, otherwise from the empty set. CRIT sets the
% criterion used by the feature evaluation routine FEATEVAL. If the dataset
% T is given, it is used as a tuning set for FEATEVAL. Alternatively
% a number of cross-validations N may be supplied. For K = 0, the optimal
% feature set (maximum value of FEATEVAL) is returned. The result W can be
% used for selecting features by B*W. In this case, features are ranked
% optimally.
% The selected features are stored in W.DATA and can be found by +W.
% In R, the search is reported step by step as:
%
% RES(:,1) : number of features
% RES(:,2) : criterion value
% RES(:,3:3+max(L,R)) : added / deleted features
%
% SEE ALSO (<a href="http://37steps.com/prtools">PRTools Guide</a>)
% MAPPINGS, DATASETS, FEATEVAL, FEATSEL
% FEATSELO, FEATSELB, FEATSELF, FEATSELI, FEATSELP, FEATSELM, PRPROGRESS
% Copyright: D. de Ridder, D.deRidder@ewi.tudelft.nl
% Faculty of Electrical Engineering, Mathematics and Computer Science
% Delft University of Technology, Mekelweg 4, 2628 CD Delft, The Netherlands
% $Id: featsellr.m,v 1.6 2009/11/27 08:53:00 duin Exp $
function [w,res] = featsellr(a,crit,ksel,l,r,t,fid)
if (nargin < 2) | isempty(crit)
prwarning(2,'No criterion specified, assuming 1-NN.');
crit = 'NN';
end
if (nargin < 3) | isempty(ksel)
ksel = 0; % Consider all the features and sort them.
end
if (nargin < 4) | isempty(l)
l = 1;
end;
if (nargin < 5) | isempty(r)
r = 0;
end;
if (nargin < 6)
prwarning(3,'No tuning set supplied.');
t = [];
end
if (nargin < 7)
fid = [];
end
if (l==r)
error('In ''Plus-L-takeaway-R feature selection'' L should be unequal to R')
end
% No inputs arguments provided, return an untrained mapping.
if (nargin == 0) | (isempty(a))
w = prmapping('featsellr',{crit,ksel,t});
w = setname(w,'+L-R FeatSel');
return
end
isvaldfile(a,1,2); % At least 1 object per class, 2 classes.
a = testdatasize(a);
a = setprior(a,getprior(a));
if ~is_scalar(t), iscomdset(a,t); end
[m,k,c] = getsize(a);
featlist = getfeatlab(a);
% Initialise counters.
state.l = l; state.r = r;
if (state.l > state.r) % Start from empty set.
state.I = []; % Pool of selected feature indices.
state.critval_opt = 0; % Maximum criterion value found so far.
state.res = zeros(0,2+max(l,r)); % Report of selection process.
else
state.I = [1:k]; % Pool of the feature indices to be used.
state.critval_opt = feateval(a,crit,t); % Max. criterion found so far.
state.res = [k,state.critval_opt,zeros(1,max(l,r))]; % Report of selection process.
end;
state.I_opt = state.I; state.critval = [];
% Evaluate the criterion function for the entire dataset A.
prprogress(fid,['\nfeatsellr: Plus-' num2str(l) '-takeaway-' num2str(r) ...
' feature selection\n'])
if ksel == 0
if (l > r)
s = sprintf('Forward search of optimal feature set out of %i: ',k);
else
s = sprintf('Backward search of optimal feature set out of %i: ',k);
end
nsize = k;
else
s = sprintf('Selection of %i features: ',ksel);
nsize = abs(length(state.I)-ksel);
end
prwaitbar(nsize,s);
while (ksel == 0 & (((state.l > state.r) & (length(state.I) <= k-state.l)) | ...
((state.l < state.r) & (length(state.I) > state.r)))) | ...
(ksel > 0 & state.l > state.r & length(state.I) < ksel) | ...
(ksel > 0 & state.l < state.r & length(state.I) > ksel)
% Calculate criterion values C for the features in I and find the feature
% giving the maximum.
if (l > r)
state = plusl(a,crit,t,state,fid);
if (r > 0), state = takeawayr(a,crit,t,state,fid); end;
prwaitbar(nsize,length(state.I),[s int2str(length(state.I))]);
else
state = takeawayr(a,crit,t,state,fid);
if (l > 0), state = plusl(a,crit,t,state,fid); end;
prwaitbar(nsize,k-length(state.I),[s int2str(length(state.I))]);
end;
end;
prwaitbar(0);
% Make sure we end up with exactly KSEL features.
if (ksel == 0) & (state.l > state.r) & (length(state.I) < k)
state.l = k - length(state.I);
state = plusl(a,crit,t,state,fid);
elseif (ksel == 0) & (state.l < state.r) & (length(state.I) > 1)
state.r = length(state.I) - 1;
state = takeawayr(a,crit,t,state,fid);
elseif (ksel > 0) & (state.l > state.r) & (length(state.I) > ksel)
state.r = length(state.I) - ksel;
state = takeawayr(a,crit,t,state,fid);
elseif (ksel > 0) & (state.l < state.r) & (length(state.I) < ksel)
state.l = ksel - length(state.I);
state = plusl(a,crit,t,state,fid);
end;
if (ksel ~= 0), state.I_opt = state.I; end;
w = featsel(k,state.I_opt);
w = setmapping_type(w,'trained');
w = setsize(w,[k length(state.I_opt)]);
if ~isempty(featlist)
w = setlabels(w,featlist(state.I_opt,:));
end
w = setname(w,'+L-R FeatSel');
res = state.res;
prprogress(fid,'featsellr finished\n')
return;
function state = plusl(a,crit,t,state,fid)
% J contains the indices of the nonselected features.
[m,k,c] = getsize(a); J = setdiff(1:k,state.I);
critval_opt = 0; Jsub_opt = [];
% Find all possible choices of L indices out of J.
[Jsub,ind] = npickk(J,state.l);
while (~isempty(Jsub))
% Calculate the criterion when subset JSUB is added.
if (isempty(t))
critval = feateval(a(:,[state.I Jsub]),crit);
elseif is_scalar(t)
critval = feateval(a(:,[state.I Jsub]),crit,t);
else
critval = feateval(a(:,[state.I Jsub]),crit,t(:,[state.I Jsub]));
end;
% Store the best subset and its criterion value.
if (critval > critval_opt)
critval_opt = critval;
Jsub_opt = Jsub;
end;
[Jsub,ind] = npickk(J,state.l,ind);
end;
state.I = [state.I Jsub_opt];
if (critval_opt > state.critval_opt) | ...
(critval_opt == state.critval_opt & ...
length(state.I) < length(state.I_opt)),
state.critval_opt = critval_opt;
state.I_opt = state.I;
end;
line = [length(state.I),critval_opt,...
[Jsub_opt zeros(1,size(state.res,2)-2-length(Jsub_opt))]];
prprogress(fid,' %d %f \n',line(1:2));
state.res = [state.res; line];
return
function state = takeawayr(a,crit,t,state,fid)
% J contains the indices of the selected features.
[m,k,c] = getsize(a); J = state.I;
critval_opt = 0; Jsub_opt = [];
% Find all possible choices of L indices out of J.
[Jsub,ind] = npickk(J,state.r);
while (~isempty(Jsub))
% Calculate the criterion when subset JSUB is removed.
if (isempty(t))
critval = feateval(a(:,[setdiff(state.I,Jsub)]),crit);
elseif is_scalar(t)
critval = feateval(a(:,[setdiff(state.I,Jsub)]),crit,t);
else
critval = feateval(a(:,[setdiff(state.I,Jsub)]),crit,...
t(:,[setdiff(state.I,Jsub)]));
end;
% Store the best subset and its criterion value.
if (critval > critval_opt)
critval_opt = critval;
Jsub_opt = Jsub;
end;
[Jsub,ind] = npickk(J,state.r,ind);
end;
state.I = setdiff(state.I,Jsub_opt);
if (critval_opt > state.critval_opt) | ...
(critval_opt == state.critval_opt & ...
length(state.I) < length(state.I_opt)),
state.critval_opt = critval_opt;
state.I_opt = state.I;
end;
line = [length(state.I),critval_opt,...
[-Jsub_opt zeros(1,size(state.res,2)-2-length(Jsub_opt))]];
prprogress(fid,' %d %f \n',line(1:2));
state.res = [state.res; line];
return
% NPICKK - A kind-of-recursive implementation of NCHOOSEK.
%
% Picks K elements out of vector V and returns them in R, with their indices
% in I. Initialise with [R,I] = NPICKK(V,K); subsequent calls should be
% [R,I] = NPICKK(V,K,I). If the set is exhausted, returns R = I = [].
function [r,i] = npickk(v,k,i)
n = length(v);
if (nargin < 3) | (isempty(i))
i = 1:k; j = 1;
else
% Increase the last index.
j = k; i(j) = i(j) + 1;
% While there's an overflow, increase the previous index
% and reset the subsequent ones.
while (i(j) > n-(k-j)) & (j >= 2)
i(j-1) = i(j-1) + 1;
for l = j:k
i(l) = i(l-1) + 1;
end;
j = j - 1;
end;
end;
% If the last index is too large, we're done.
if (i(end) > n)
r = []; i = [];
else
r = v(i);
end;
return