lsqr result strongly depends on weights

I need to solve argmin||W*FT^-1(Ax)-W*p|| using lsqr. p is image, x is k-space matrix, A is a matrix and W is a weighting matrix. In order to pass them to matlab lsqr, I vectorized p, x and W. This is my code:

 b=W.*p;
 x=lsqr(@(x1,modo)FUNC(x1,W_vector,A,modo),b,tol,maxit);
 X_k_lsqr=reshape(x,dim(1),dim(2),dim(3));
 X_lsqr=real(ifftn(X_k_lsqr)).*MASK(:,:,:,1);

%% Auxiliary function

    function [result modo]=FUNC(x1,W_vector,A,modo)
        %Computes y=A*x for modo='notransp'
        %Computes y=A'*x for modo='transp'

        switch modo
            case 'notransp'
                res=A*x1;
                R1=reshape(res,norient,dim(1)*dim(2)*dim(3));
                for co=1:norient
                    R2(:,:,:,co)=reshape(R1(co,:),dim(1),dim(2),dim(3));
                    FR(:,:,:,co)=ifftn(R2(:,:,:,co));
                    aux=FR(:,:,:,co);
                    R3(co,:)=aux(:).';
                end
                result=W.*R3(:); 

            case 'transp'
                RR1=reshape(x1./(W+eps),norient,dim(1)*dim(2)*dim(3));
                for co=1:norient
                    RR2(:,:,:,co)=reshape(RR1(co,:),dim(1),dim(2),dim(3));
                    FRR(:,:,:,co)=fftn(RR2(:,:,:,co));
                    aux=FRR(:,:,:,co);
                    RR3(co,:)=aux(:).'; 
                end
                result=A'*RR3(:); 
        end
    end
  end

As W appears in both terms of the minimization problem, I would have expected the image I obtain as a result to be almost independent on W values. The image looks qualitatively the same if I change W, but its values strongly depend on W. I don't know if something is wrong with my code. Should I actually obtain almost the same values for different W? Thank you very much for your help.

Answer 1

[Cam Salzberger]

Hello Matteo,

Your code is a bit involved and specific for me to determine if it would be expected that the image values would change as a result of the weight values changing. Also, do you mean scaling the weights, or modifying them?

However, I can say that I believe weighted least squares is usually done by trying to solve argmin(W*||FT^-1(Ax)-p||), rather than having the weights inside the error magnitude calculation.

Instead of using "lsqr", you may wish to check out other functions that are capable of doing weighted least squares , such as "lscov" in base MATLAB, or "robustfit" in the Statistics Toolbox.

I hope this helps with at least testing if you are receiving the expected results.

-Cam