Solving Optimisation Problem with Rank Constraint in MATLAB

I have a typical least squares problem, i.e I have to find the value of x that minimizes norm of C∗x(:)−d. C is 180x16 matrix(C is rank deficient,i.e rank(C)=7), x is 4x4 matrix & d is 180x1 vector. However, I have a constraint that rank(x)=1. If x was a 16x1 vector and didn't have rank constraint, this problem could be easily solved by using y = pinv(C)*d in MATLAB. But since x is a matrix and has rank constraint, I am not able to proceed further. I would be grateful if someone provides me hint or suggestion to tackle this problem.

Answer 1

[Torsten]

1. Solve y=pinv(C)*d

2. Determine the best rank-1 - approximation x to y as discussed in the previous thread:

http://www.mathworks.com/matlabcentral/answers/255420-representing-a-matrix-as-a-product-of-two-matrices-in-matlab

My guess is that x solves your original problem, but I'm not 100% certain.

Best wishes

Torsten.

[rihab]

I had exactly done the same before, but I have noticed that estimate of x obtained after svd (like in the thread you posted) is not a perfect match of y. So I was wondering if there is some other technique to somehow obtain a better match!

[Torsten]

To be sure, you could proceed as follows:

Write x as u*v' for two unknown 4x1-vectors u=(u1,u2,u3,u4) and v=(v1,v2,v3,v4).

Then use fminsearch to minimize |C*u*v'-d|_2 in the 8 unknowns u1,u2,u3,u4,v1,v2,v3,v4.

Best wishes

Torsten.

[rihab]

I appreciate your suggestion. I am using fminsearch now to solve it. But what value of "x0" should I give for this problem(x = fminsearch(fun,x0))?

I have written the function as follows:

function fun = optim(x,A,b)

residuals = (A*x) - b;

fun = sum(residuals.^2);

end

[Torsten]

function main
A=...;
b=...;
x0=ones(8,1);
x=fminsearch(@)(x)optim(x,A,b),x0);

function fun=optim(x,A,b)
u(1:4,1) = x(1:4);
v(1:4,1) = x(5:8);
rank1 = u*v';
y = reshape(rank1,[16,1]);
residuals = A*y-b;
fun = sum(residuals.^2);

Best wishes

Torsten.

[rihab]

thank you :)

[Torsten]

Try A and b real-valued first.

For complex A and b, the program has to be adapted.

Best wishes

Torsten.

[rihab]

Well real values work but I think complex values would provide a better estimate (because A and b in general are complex numbers ;) ). What should be added in the case of complex valued A and b?

[Torsten]

I did not test it, but maybe something like

function main
A=...;
b=...;
x0=ones(16,1);
x=fminsearch(@)(x)optim(x,A,b),x0);

function fun=optim(x,A,b)
realu(1:4,1) = x(1:4);
imagu(1:4,1) = x(5:8);
realv(1:4,1) = x(9:12);
imagv(1:4,1) = x(13:16);
u=complex(realu,imagu);
v=complex(realv,imagv);
rank1 = u*v';
y = reshape(rank1,[16,1]);
residuals = A*y-b;
fun = abs(residuals);

Best wishes

Torsten.

[rihab]

thank you so much :)

[Torsten]

In the last line of the code, you will have to replace

fun = abs(residuals)

by

fun = residuals'*residuals
or
fun = norm(residuals)

Best wishes

Torsten.