Null space solutions in the presence of noise

Usually in linear algebra I end up posing problems of the form

A * x = b

These are usually easy to solve in a noise-robust manner with standard least squares approaches. I often end up with very big arrays where I would use lsqminnorm to solve A^-1 * b.

However I am currently faced with a problem of the form

A * x = 0

I have tried a variety of approaches to solve this: null(), lsqlin(), fmincon().... these all work fine when there is no noise, but return either garbage (null) or a fairly poor solution in the presence of a modest amount of noise. It seems that with noise, the true solution no longer minimises A*x.

Does anyone have any advice for working in the null space, in the presence of noise?

Thanks in advance for your time.

[KALYAN ACHARJYA]

The following link mayhelp you

https://www.math.colostate.edu/~gerhard/M345/CHP/ch7_4.pdf

Answer 1

[Matthew Dvorsky]

This is known as the homogeneous linear least squares problem.

You are trying to minimize the value of . The issue is that we can arbitrarily scale any arbitrary vector to make the value as small as we like (including the trival solution ). Thus, to find a unique solution, we must restrict the value of . So we can instead ask to minimize by choosing such that .

This minimization can be done by computing the singular value decomposition of A, whcih decomposes A as . If we take to be a column of V whose corresponding singular value is λ, then . Thus, we can minimize by choosing to be the column of V with the smallest corresponding singular value. This is demonstrated in the following MATLAB example.

% Generate A and x such that Ax = 0
A = rand(4, 3) * rand(3, 4);
x = null(A)         % Returns a single null vector

x = 4×1

-0.7998 0.4316 -0.1600 0.3854

% Add noise to A
A_noise = A + 0.00001 * rand(size(A));
null(A_noise)       % Does not work due to noise

ans =

  4×0 empty double matrix

% Calculate least square solution using SVD.
[~, ~, V] = svd(A_noise);
x_lsq = V(:, end)   % Gives a good approximation for x (or sometimes -x).

x_lsq = 4×1

-0.7998 0.4316 -0.1600 0.3853

Note that this method even works for non-square matrices A.

Edit: As Bruno pointed out in his answer, the we can use the MATLAB 'svds' function, to compute only the column of V that corresponds to the smallest singular value. This is efficient even for very large sparse matrices, as noted in the documentation of 'svds'.

[~, ~, x_lsq] = svds(A_noise, 1, "smallest")    % Same result as the above code.

x_lsq = 4×1

-0.7998 0.4316 -0.1600 0.3853

[Phillip]

This is very handy and well described, thank you.

I haven't been using PCA for my problem as the arrays are very large (out of memory error). However, the array I am inverting is sparse and hence I have been using the Matlab function "lsqminnorm", which uses optimization to attemtpt to find the global minimum.

Would you happen to know of a suitable Matlab alternative to SVD in, which can handle large arrays in order to estimate V as in this example?

[John D'Errico]

You can use one of the options in SVD for larger problems.

A = randn(50000,299)*randn(299,300) + randn(50000,300)*1e-8;
tic,[~,~,V] = svd(A);toc
Elapsed time is 24.572171 seconds.

tic,[~,~,V] = svd(A,0);toc
Elapsed time is 0.432794 seconds.

I could think of some other ideas too, but SVD is best.

[Matthew Dvorsky]

As Bruno pointed out in his answer, we can use the svds function, which will be more efficient, especially when the input matrix is large. I have edited my answer to give an example.

[Phillip]

Thanks all, this is very useful for me -- and thanks for your many contributions to the file exchange over the years, especially Bruno and John! I have all manner of snippets of your code embedded into my workflow.

Answer 2

[Bruno Luong]

You can call

[x_lsq ,a] = eigs(@(x) (S'*S)\x, size(S,2), 1, 'smallestabs')

where S is your (sparse) matrix

sqrt(a) is the smallest singular value of S, and x_lsq is the corresponding singular vector.

Or

[~,~,x_lsq] = svds(S,1,0)

Span x_lsq approximates the null space of S, assuming its dimension is 1.