MATLAB Answers

0

How to get all the solutions to an underdetermined system of linear equations using lsqnonneg()?

Asked by MUKTIKANTA on 13 Oct 2018 at 1:19
Latest activity Commented on by MUKTIKANTA on 14 Oct 2018 at 0:35

I am trying find out all the positive solutions to an underdetermined system of linear equations and my code is as follows:

 format rat
a=[3 22 119 620 3 1 3 3 2 20 20 22 22 119 119 3 3 3 3 3 3 1 1 1 20 20 20 20 2 1 2 20 20 3 3 3 3 3 3 1 1 1 2 2 1;
    3 0 0 0 3 1 3 3 2 0 0 0 0 0 0 3 3 3 3 3 3 1 1 1 0 0 0 0 2 1 2 0 0 3 3 3 3 3 3 1 1 1 2 2 1;
    0 22 0 0 3 1 0 0 0 20 20 22 22 0 0 3 3 3 3 0 0 1 1 1 20 20 20 20 2 1 0 20 20 3 3 3 3 3 3 1 1 1 2 2 1;
    0 0 119 0 0 0 3 0 0 20 20 0 0 119 119 3 3 0 0 3 3 1 0 0 20 20 20 20 0 0 2 20 20 3 3 3 3 3 3 1 1 1 2 2 1;
    0 0 0 620 0 0 0 3 2 0 0 22 22 119 119 0 0 3 3 3 3 0 1 1 20 20 20 20 2 1 2 20 20 3 3 3 3 3 3 1 1 1 2 2 1];
b=[620;3;22;119;620];
 w=lsqnonneg(a,b);

Each time I run this code, it gives me only one solution. How can I get all the solution? Is there any method available in MATLAB to achieve this?

  1 Comment

That is rather the wrong function to use.

An underdetermined system of linear equations describes hyperplane (or possibly a plane or a line.)

You have 45 variables and 5 equations. Given any arbitrary values for the last 40 variables, exact values of the first 5 can be calculated.

The set of all positive solutions could only be described by boundary equations. I would have to work out how many of those you would need. It is late at night here and one avenue of thought is suggesting that you might need 45!/5! so I would want to think about that when I am more awake. Maybe it is 45!/40! ?

Sign in to comment.

1 Answer

Answer by Bruno Luong on 13 Oct 2018 at 6:55
Edited by Bruno Luong on 13 Oct 2018 at 7:17
 Accepted Answer

The number of solutions is infinity. Explanation:

>> size(a)
ans =
       5    45

In your case you have 5 eqt and 45 unknowns, and

>> rank(a)
ans =
       5 

confirms your matrix is full row-rank. meaning all the 5 rows of A are independent vectors in R^45.

One particular solution can be obtained with backshlash (or lsqnonneg if you prefer)

>> x0=a\b;
>> fprintf('x0 = %s\n', mat2str(x0));
x0 = [0;0;5.93708057585062e-17;0.77741935483871;0;0;0;0;0;0;0;0.863636363636363;0;0.974789915966387;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0.999999999999998;0;0;0;0;0;0;0;0;0;0;0]

The kernel of a, meaning a basis of equation A*X = 0 can be obtained with NULL command.

>> X=null(a);
>> size(X)
ans =
      45    40  

It tells that the kernel has dimension 40 and a maths theorem tells it's equals to size(a,2)-rank(a).

To confirm X is the kernel, you can check with frobenus norm of max values (they are 0 at the numerical precision)

>> y=a*X;
>> norm(y,'fro')
ans =
     1.5079e-13
>> max(max(abs(y)))
ans =
     7.5533e-14

Now all solutions of a*x=b can be written as

x = x0 + X*r

where r is any arbitrary vector of size 40 x 1. The set of solutions is an affine space of dimension 40 and there are of course an infinity vectors in this huge space. To check it try:

r = rand(40,1); % any random vector
x = x0 + X*r; % general form of the solution a*x=b

check if a*x = b

>> max(abs(a*x-b))
ans =
     1.1369e-13

If you add the constraints (x >= 0), meaning

x(i) >= 0 for all i = 1,...45

Then the general solution is the intersection of 45 half-hyperspace and the affine space of dimension 40.

It will be a polytope of dimension 40 living in R^45. It is still huge. This polytope is actually a an interior of the convex combination of n vertices, n can be as high as:

>> nchoosek(45,40)
ans =
       1221759

each of them is one particular solution of

    a*x = b
    x >= 0

There is no way you can enumerate them all and make on top of that an arbitrary sub-convex combination if that is your intention.

  3 Comments

Thanks, Bruno. I would have found it tough to pull all those steps together.

Note: what I wrote is not entirely correct, sometime the polytope is unbounded if no other constraint on x is introduced and it can not express as the convex combination of finite set of vertexes.

However I hop it discourage people trying to beat the market with such brute force method.

Thanks, Bruno for such a well-explained answer.

Sign in to comment.