annotate DL/Majorization Minimization DL/mm1.m @ 214:b9b4dc87f1aa luisf_dev

Additional help comments in ~/DL/Majorization Minimization DL/wrapper_mm_DL.m.
author Aris Gretsistas <aris.gretsistas@elec.qmul.ac.uk>
date Wed, 21 Mar 2012 18:27:23 +0000
parents b14209313ba4
children
rev   line source
ivan@155 1 function [unhat,er] = mm1(Phi,x,u0,to,lambda,maxIT,eps,map)
ivan@155 2 %% Iterative Soft Thresholding (with optional debiasing)
ivan@155 3 %
ivan@155 4 % Phi = Normalized Dictionary
ivan@155 5 % x = Signal(x). This can be a vector or a matrix
ivan@155 6 % u0 = Initial guess for the coefficients
ivan@155 7 % to = 1/(step size) . It is larger than spectral norm of dictionary Phi
ivan@155 8 % lambda = Lagrangian multiplier. (regulates shrinkage)
ivan@155 9 % eps = Stopping criterion for iterative softthresholding and MM dictionary update
ivan@155 10 % map = Debiasing. 0 = No, 1 = Yes
ivan@155 11 % unhat = Updated coefficients
ivan@155 12 % er = Objective cost
ivan@155 13 %%
ivan@155 14 cont = 1;
ivan@155 15 in = 1;
ivan@155 16 % un = zeros(size(u0,1),size(u0,2));
ivan@155 17 un = u0;
ivan@155 18 c1 = (1/to^2)*Phi'*x;
ivan@155 19 c2 = (1/to^2)*(Phi'*Phi);
ivan@155 20 %%%%
ivan@155 21 while (cont && (in<=maxIT))
ivan@155 22 unold = un;
ivan@155 23 %%%%%% Soft Thresholding %%%%%%%
ivan@155 24 alphap = (un + c1 - c2*un);
ivan@155 25 un = (alphap-(lambda/(2*to^2))*sign(alphap)).*(abs(alphap)>=(lambda/(2*to^2)));
ivan@155 26 in = in+1;
ivan@155 27 cont = sum(sum((unold-un).^2))>eps;
ivan@155 28 end
ivan@155 29 %%%%%%%%%%
ivan@155 30 if map == 1,
ivan@155 31 %% Mapping on the selected space %%%%
ivan@155 32 [uN,uM] = size(un);
ivan@155 33 unhat = zeros(uN,uM);
ivan@155 34 for l = 1:uM,
ivan@155 35 unz = (abs(un(:,l))>0);
ivan@155 36 M = diag(unz);
ivan@155 37 PhiNew = Phi*M;
ivan@155 38 PhiS = PhiNew(:,unz);
ivan@155 39 unt = inv(PhiS'*PhiS+.0001*eye(sum(unz)))*PhiS'*x(:,l);
ivan@155 40 unhat(unz,l) = unt;
ivan@155 41 end
ivan@155 42 else
ivan@155 43 unhat = un;
ivan@155 44 end
ivan@155 45 %%% Cost function calculation
ivan@155 46 if map == 1,
ivan@155 47 er = sum(sum((Phi*unhat-x).^2))+lambda*(sum(sum(abs(unhat)>0))); %% l_0 Cost function
ivan@155 48 else
ivan@155 49 er = sum(sum((Phi*unhat-x).^2))+lambda*(sum(sum(abs(unhat))));
ivan@155 50 end