[Master Index] [Index for Toolbox]

# muopt

## (Toolbox/muopt.m in BrainStorm 2.0 (Alpha))

### Function Synopsis

    [reg_min] = muopt(U,s,b,method)

### Help Text

    MUOPT - Quasi-optimality criterion for choosing the regularization parameter.
     function [reg_min] = muopt(U,s,b,method)
    QUASIOPT Quasi-optimality criterion for choosing the regularization parameter.

     [reg_min,Q,reg_param] = quasiopt(U,s,b,method)
     [reg_min,Q,reg_param] = quasiopt(U,sm,b,method)  ,  sm = [sigma,mu]

     Plots the quasi-optimality function Q for the following methods:
        method = 'Tikh' : Tikhonov regularization   (solid line )
        method = 'tsvd' : truncated SVD or GSVD     (o markers  )
        method = 'dsvd' : damped SVD or GSVD        (dotted line)
     If no method is specified, 'Tikh' is default.

     If any output arguments are specified, then the minimum of Q is
     identified and the corresponding reg. parameter reg_min is returned.

     Reference: T. Kitagawa, "A deterministic approach to optimal
     regularization - the finite dimensional case", Japan J. Appl.
     Math. 4 (1987), 371-391.

     Per Christian Hansen, UNI-C, 03/17/93.

     Adapted from Hansen for BrainStorm

### Cross-Reference Information

This function calls

    dsvd                               C:\BrainStorm_2001\PublicToolbox\regutools\dsvd.m
    tsvd                               C:\BrainStorm_2001\PublicToolbox\regutools\tsvd.m

### Listing of function C:\BrainStorm_2001\Toolbox\muopt.m

    function [reg_min] = muopt(U,s,b,method)
    %MUOPT - Quasi-optimality criterion for choosing the regularization parameter.
    % function [reg_min] = muopt(U,s,b,method)
    %QUASIOPT Quasi-optimality criterion for choosing the regularization parameter.
    %
    % [reg_min,Q,reg_param] = quasiopt(U,s,b,method)
    % [reg_min,Q,reg_param] = quasiopt(U,sm,b,method)  ,  sm = [sigma,mu]
    %
    % Plots the quasi-optimality function Q for the following methods:
    %    method = 'Tikh' : Tikhonov regularization   (solid line )
    %    method = 'tsvd' : truncated SVD or GSVD     (o markers  )
    %    method = 'dsvd' : damped SVD or GSVD        (dotted line)
    % If no method is specified, 'Tikh' is default.
    %
    % If any output arguments are specified, then the minimum of Q is
    % identified and the corresponding reg. parameter reg_min is returned.
    %
    % Reference: T. Kitagawa, "A deterministic approach to optimal
    % regularization - the finite dimensional case", Japan J. Appl.
    % Math. 4 (1987), 371-391.
    %
    % Per Christian Hansen, UNI-C, 03/17/93.
    %
    % Adapted from Hansen for BrainStorm

    %<autobegin> ---------------------- 26-May-2004 11:31:14 -----------------------
    % --------- Automatically Generated Comments Block Using AUTO_COMMENTS ---------
    %
    % CATEGORY: Inverse Modeling
    %
    % At Check-in: $Author: Mosher $  $Revision: 14 $  $Date: 5/26/04 10:00a $
    %
    % This software is part of BrainStorm Toolbox Version 2.0 (Alpha) 24-May-2004
    %
    % Principal Investigators and Developers:
    % ** Richard M. Leahy, PhD, Signal & Image Processing Institute,
    %    University of Southern California, Los Angeles, CA
    % ** John C. Mosher, PhD, Biophysics Group,
    %    Los Alamos National Laboratory, Los Alamos, NM
    % ** Sylvain Baillet, PhD, Cognitive Neuroscience & Brain Imaging Laboratory,
    %    CNRS, Hopital de la Salpetriere, Paris, France
    %
    % See BrainStorm website at http://neuroimage.usc.edu for further information.
    %
    % Copyright (c) 2004 BrainStorm by the University of Southern California
    % This software distributed  under the terms of the GNU General Public License
    % as published by the Free Software Foundation. Further details on the GPL
    % license can be found at http://www.gnu.org/copyleft/gpl.html .
    %
    % FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE
    % UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY
    % WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF
    % MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY
    % LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE.
    %<autoend> ------------------------ 26-May-2004 11:31:14 -----------------------

    % Set defaults.
    npoints = 80;  % Number of points for 'Tikh' and 'dsvd'.
    if (nargin==3), method = 'Tikh'; end   % Default method.

    % Initialization.
    [m,n] = size(U); [p,ps] = size(s);
    if (ps==2), s = s(p:-1:1,1)./s(p:-1:1,2); U = U(:,p:-1:1); end
    xi = (U'*b)./s;
    if (nargout > 0), find_min = 1; else find_min = 0; end

    % Compute the quasioptimality function Q.
    if (method(1:4)=='Tikh' | method(1:4)=='tikh')

      Q = zeros(npoints,1); reg_param = Q;
      reg_param(npoints) = s(p);
      ratio = (s(1)/s(p))^(1/(npoints-1));
      for i=npoints-1:-1:1, reg_param(i) = ratio*reg_param(i+1); end
      for i=1:npoints
        f = (s.^2)./(s.^2 + reg_param(i)^2);
        Q(i) = norm((1 - f).*f.*xi);
      end

    elseif (method(1:4)=='dsvd' | method(1:4)=='dgsv')

      Q = zeros(npoints,1); reg_param = Q;
      reg_param(npoints,1) = s(p);
      ratio = (s(1)/s(p))^(1/(npoints-1));
      for i=npoints-1:-1:1, reg_param(i) = ratio*reg_param(i+1); end
      for i=1:npoints
        f = s./(s + reg_param(i));
        Q(i) = norm((1 - f).*f.*xi);
      end

    elseif (method(1:4)=='tsvd' | method(1:4)=='tgsv')

      Q = abs(xi); reg_param = [1:p]';

    else, error('Illegal method'), end

    % Plot the function.
    if (method=='tsvd' | method=='tgsv')
      semilogy(reg_param,Q,'o'), xlabel('k')
      title('Quasi-optimality function')
      if (find_min)
        [minQ,minQi] = min(Q); reg_min = reg_param(minQi);
        HoldState = ishold; hold on;
        semilogy([reg_min,reg_min],[minQ,minQ/1000],'--')
        if (~HoldState), hold off; end
        title(['Quasi-optimality function, minimum at ',num2str(reg_min)])
      end
    else
    %  if (method(1:4)=='Tikh' | method(1:4)=='tikh' | ...
    %      method(1:4)=='dsvd' | method(1:4)=='dgsv'      )
    %    loglog(reg_param,Q), xlabel('lambda')
    %  else
    %    loglog(reg_param,Q,':'), xlabel('lambda')
    %  end
    %  title('Quasi-optimality function')
      if (find_min)
        [minQ,minQi] = min(Q); reg_min = reg_param(minQi);
    %   HoldState = ishold; hold on;
    %   loglog([reg_min,reg_min],[minQ,minQ/1000],'--')
    %    if (~HoldState), hold off; end
    %    title(['Quasi-optimality function, minimum at ',num2str(reg_min)])
      end
    end

Produced by color_mat2html, a customized BrainStorm 2.0 (Alpha) version of mat2html on Tue Oct 12 12:05:14 2004
Cross-Directory links are: ON