mfbox_rel_paramsweep_run.m

toolbox for spm 5 for data, model free analysis

function [mfbss,params]=mfbox_rel_paramsweep_run(mfbss,params,runflag)
% run bss algorithm with varying parameters
%
% Usage:
%  [mfbss,params]=mfbox_rel_paramsweep_run(mfbss,params,runflag)
%
%  mfbss   - mfbss structure (see mfbox_init)
%            name - bss algorithms
%  params  - struct with
%            threshold_type - threshold components by
%            threshold      - threshold value
%            parameter1     - parameter1 to vary
%            range_from1    - from
%            range_to1      - to
%            parameter2     - parameter2 to vary
%            range_from2    - from
%            range_to2      - to
%            parameter3     - parameter3 to vary
%            range_from3    - from
%            range_to3      - to
%            sortby         - temporalspatial
%                             temporal
%                             spatial
%                             variance
%  runflag - -1 get default parameter
%             0 interactive ask parameters
%             1 interactive ask parameters and run
%             2 run
%
% Copyright by Peter Gruber and Fabian J. Theis
% Signal Processing & Information Theory group
% Institute of Biophysics, University of Regensburg, Germany
% Homepage: http://research.fabian.theis.name
%           http://www-aglang.uni-regensburg.de
%
% This file is free software, subject to the 
% GNU GENERAL PUBLIC LICENSE, see gpl.txt

error(nargchk(1,3,nargin));
error(nargchk(1,2,nargout));

if (nargin<2), params = []; end
if (nargin<3), runflag = 1; end

params = mfbox_checkparam(params,'rel','paramsweep', ...
    struct('steps',10,'range_from1',0,'range_to1',0, ...
    'parameter1','','range_from2',0,'range_to2',0, ...
    'parameter2','','range_from3',0,'range_to3',0, ...
    'parameter3','','sortby','temporalspatial','threshold',0.01, ...
    'threshold_type','none'),'experimental');

if (abs(runflag-0.5)<1)
    params = mfbox_rel_paramsweepg(mfbss,params,runflag);
end

if (runflag>0 && isstruct(params))
    if (runflag<3)
        prgs = mfbox_progress([],'title','BSS analysis','string', ...
            sprintf('Running BSS\nanalysis...'),'progress',[1,4+params.steps]);
    end

    [tg,mask,smask,xmask] = mfbox_getmask(mfbss.grid,mfbss.mask);

    if (~isempty(mfbss.reference)), reference = mfbss.reference(mask);
    else reference = [];
    end

    cmd_func = str2func(sprintf('mfbox_%s_run',mfbss.name));
    if (runflag<3)
        mfbox_progress(prgs,'string',sprintf('Running BSS\nanalysis...'),'progress',[2,4+params.steps]);
    end

    n = params.steps;
    As = cell(1,n); Ws = cell(1,n); Ss = cell(1,n);

    paramevo = struct();
    for i=1:n
        rparams = mfbss.params.(mfbss.name);

        v = (i-1)/(n-1);
        boolind = 0;
        if (~isempty(params.parameter1) && isfield(rparams,params.parameter1))
            p = 1;
            if (strcmp(rparams.([params.parameter1,'_type']),'int'))
                p = round((1-v)*params.range_from1+v*params.range_to1);
                rparams.(params.parameter1) = p;
            elseif (strcmp(rparams.([params.parameter1,'_type']),'float'))
                p = (1-v)*params.range_from1+v*params.range_to1;
                rparams.(params.parameter1) = p;
            elseif (strcmp(rparams.([params.parameter1,'_type']),'boolean'))
                p = mod(floor((i-1)/(2^boolind)),2)+1;
                t = {'off','on'};
                rparams.(params.parameter1) = t{p};
                boolind = boolind+1;
            end
            if (isfield(paramevo,params.parameter1)), x = paramevo.(params.parameter1);
            else x = [];
            end
            x(i) = p;
            paramevo.(params.parameter1) = x;
        end
        if (~isempty(params.parameter2) && isfield(rparams,params.parameter2))
            p = 1;
            if (strcmp(rparams.([params.parameter2,'_type']),'int'))
                p = round((1-v)*params.range_from2+v*params.range_to2);
                rparams.(params.parameter2) = p;
            elseif (strcmp(rparams.([params.parameter2,'_type']),'float'))
                p = (1-v)*params.range_from2+v*params.range_to2;
                rparams.(params.parameter2) = p;
            elseif (strcmp(rparams.([params.parameter2,'_type']),'boolean'))
                p = mod(floor((i-1)/(2^boolind)),2)+1;
                t = {'off','on'};
                rparams.(params.parameter2) = t{p};
                boolind = boolind+1;
            end
            if (isfield(paramevo,params.parameter2)), x = paramevo.(params.parameter2);
            else x = [];
            end
            x(i) = p;
            paramevo.(params.parameter2) = x;
        end
        if (~isempty(params.parameter3) && isfield(rparams,params.parameter3))
            p = 1;
            if (strcmp(rparams.([params.parameter3,'_type']),'int'))
                p = round((1-v)*params.range_from3+v*params.range_to3);
                rparams.(params.parameter3) = p;
            elseif (strcmp(rparams.([params.parameter3,'_type']),'float'))
                p = (1-v)*params.range_from3+v*params.range_to3;
                rparams.(params.parameter3) = p;
            elseif (strcmp(rparams.([params.parameter3,'_type']),'boolean'))
                p = mod(floor((i-1)/(2^boolind)),2)+1;
                t = {'off','on'};
                rparams.(params.parameter3) = t{p};
            end
            if (isfield(paramevo,params.parameter3)), x = paramevo.(params.parameter3);
            else x = [];
            end
            x(i) = p;
            paramevo.(params.parameter3) = x;
        end

        if (isempty(smask))
            [f,A,W,S] = cmd_func(mfbss.X,mask,rparams,max(2,runflag));
        else
            tX = mfbox_getmaskdata(mfbss.X,tg,mask);
            [f,A,W,S] = cmd_func(tX,mask,rparams,max(2,runflag));
        end
        vAp = mfbox_calcthreshold(A,S,params.threshold_type,params.threshold, ...
            mfbss.design,reference);
        As{i} = A(:,vAp); Ws{i} = W(vAp,:); Ss{i} = single(S(:,vAp));
        clear A W S
        if (runflag<3)
            mfbox_progress(prgs,'string',sprintf('Running BSS\nanalysis...'),'progress',[2+i,4+params.steps]);
        end
    end

    if (runflag<3)
        mfbox_progress(prgs,'string',sprintf('Running BSS\nanalysis...'),'progress',[2+params.steps,4+params.steps]);
    end

    H = cell(1,length(As));
    for i=1:length(As)-1
        [e,hc] = mfbox_cumindependence(As{i},As{i+1});
        [e,ho] = mfbox_overlap(Ss{i},Ss{i+1});
        v = complex(hc,ho);
        y = mfbox_getsort(complex(var(As{i+1}),var(Ss{i+1})),v,params.sortby);
        As{i+1} = As{i+1}(:,y);
        Ws{i+1} = Ws{i+1}(y,:);
        Ss{i+1} = Ss{i+1}(:,y);
        H{i} = v(:,y);
    end
    H{length(As)} = paramevo;
    mfbss.extraplot = {@mfbox_plotsweepwin,H};

    if (runflag<3)
        mfbox_progress(prgs,'string',sprintf('Running BSS\nanalysis...'),'progress',[3+params.steps,4+params.steps]);
    end

    part = cell(1,length(As)); n = 1;
    for i=1:length(As), part{i} = n:(n+size(As{i},2)-1); n = n+size(As{i},2); end
    
    Ss = cat(2,Ss{:});
    m = repmat(1:length([part{:}]),size(mfbss.X,2),1);
    v = 1./std(Ss);
    mfbss.A = cat(2,As{:})./(v(m));
    mfbss.W = cat(1,Ws{:}).*(v(m))';
    mfbss.S = struct('mask',xmask,'map',m,'dat',Ss*diag(v), ...
        'part',{part});

    if (runflag<3), mfbox_progress(prgs,'close',[]); end
end

return