function out = performICA(in)
% ICA computation and component selection function
% (c) Jakub Kuzilek, Lenka Lhotska
% http://bio.felk.cvut.cz/~kuziljak/ E-mail: jakub.kuzilek@gmail.com
% Version: 1.0 Last update: 07/01/2013.
% (Version: 1.0, 07/01/2013)
%
%======================================================
%
% PURPOSE: This function deploys JADE ICA algorithm and selects
% components containing ECG.
%
% MANDATORY INPUT ARGUMENTS
% in ..... input data MxN, M - length of data, N - number of leads
% OPTIONAL INPUT ARGUMENTS
% none
% OUTPUT ARGUMENTS
% out .... selected components MxO, where O<=N
%=========================================================
%
% Related Bibliography:
%
% [1] Kuzilek J, Lhotska L. Electrocardiogram beat detection
% enhancement using independent component analysis. Medical
% Engineering and Physics. 2012.
%
%=========================================================
M = size(in,2); % input ecg data size
% if there is small number of mixtures perform trick
% to double them (widely used in ICA) - move signals
% by one sample and append them to original data
if M < 4
in = [in, [in(2:end,:);zeros(1,M)]];
end
% perform JADE ICA algorithm
in = real((jadeR(in')*in')');
% compute kurtosis of components
x0 = in - ones(size(in,1),1)*mean(in);
s2 = mean(x0.^2); % this is the biased variance estimator
m4 = mean(x0.^4);
kurt = m4 ./ s2.^2;
% kurt = kurtosis(in);
% sort and normalize them
[skurt, order] = sort(kurt, 'descend');
skurt = cumsum(skurt);
skurt = skurt/skurt(end);
% choose componets whose sum of kurt is larger than 75 percents of sum(kurt)
order = order(1:find(skurt >= 0.76));
% output data
out = in(:,order);