Matrix left division with constraints?
5 views (last 30 days)
Show older comments
I have an underdetermined system of equations A*x=b, and I am looking for any non-negative solution. I.e. A is an n-by-m matrix with n < m, and I need to solve the system for x subject to the following two constraints:
- sum(x) == 1
- all(x >= 0)
I can solve this using linprog, but it is fairly slow (I have to solve 1e6 equations of this form, and A is relatively large, on my machine this was going to take about 6 days).
On the other hand I just discovered that matrix left division will solve this problem orders of magnitude faster, but I don't know how to enforce the non-negativity constraint (the summing to 1 constraint is easy as an additional row of 1's in A and b). On my machine it took about 2 min to solve all 1e6 systems of equations (but obviously the solutions were not non-negative).
I suppose I can use lsqlin or lsqnonneg, but these seem to be unnecessary for my needs. I know that there are an infinite number of solutionns for every system of equations I want to solve, and I don't care which solution I get (i.e. doesn't have to have the smallest norm, or anything), so long as the non-negativity and summing to 1 constraints are satisfied. Is there a way to take advantage of the speed of matrix left division and also enforce non-negativity?
0 Comments
Answers (2)
Shashank Prasanna
on 29 Aug 2013
Your best bet is LSQLIN if you have constraints for a linear system.
Since your have a large system, do you know if it is sparse? If it is then you can make the matrix sparse and LSQLIN will handle it.
Why do you say that they are "unnecessary for your needs" ? Did LSQLIN not work or not provide a solution?
Backslash or MLDIVIDE is optimized for unconstrained linear system solution. I am not aware of a way to modify \ to honor these constraints.
1 Comment
John D'Errico
on 7 May 2023
And of course, you cannot simply modify backslash(mldivide) to honor constraints. Except that the tool which does solve that problem is lsqlin OR linprog. For the huge problem size posed by the OP, this is not a surprise. The constraints mean that any speed gains found in backslash are no longer possible.
Big problems take big time.
Mohamed Abdelhameed
on 24 Mar 2018
Edited: Walter Roberson
on 7 May 2023
I am trying to solve an overdetermenind system of equations (45 equations with 18 unknowns) in the form of A*X = B; where:
A = [8.19 0.3; 12.39 0.86; 16.15 1.68; 17.7 2.16; 24.6 8.33];
A is 5 rows, 2 columns
B = [0.85 3.8 0.225 0.175 0.05 0.05 0 1.919 0.165; 1.5 6.825 0.45 0.4 0.175 0.125 0.025 2.907 0.315; 3 7.625 0.9 1.05 0.475 0.225 0.1 3.823 0.544; 3.05 7.9 1.1 1.8 1.225 0.25 0.45 4.153 0.881; 5.1 12.95 1.975 3.65 3.075 0.25 1.075 4.045 1.217];
B is 5 rows, 9 columns
X = [a1 b1 c1 d1 e1 f1 g1 h1 i1;a2 b2 c2 d2 e2 f2 g2 h2 i2];
X is 2 rows, 9 columns
I need to solve X with the following constraints: sum (a1 to i1) = 1 & sum (a2 to i2) = 1 & all elements in X shall be higher than or equal to zero (non negative).
4 Comments
John D'Errico
on 7 May 2023
@Benjamin Kelm - Why have you not found an answer? This is a trivial probem to solve using lsqlin, at least for a small problem. And sereral people have suggested lsqlin. So saying you have not found an answer just means you have not looked.
It is even easier to formulate using the problem based tools in the optimization toolbox, where you don't even need to understand how to call lsqlin.
Torsten
on 7 May 2023
Edited: Torsten
on 7 May 2023
X = sym('X',[2 9]);
A = [8.19 0.3; 12.39 0.86; 16.15 1.68; 17.7 2.16; 24.6 8.33];
B = [0.85 3.8 0.225 0.175 0.05 0.05 0 1.919 0.165; 1.5 6.825 0.45 0.4 0.175 0.125 0.025 2.907 0.315; 3 7.625 0.9 1.05 0.475 0.225 0.1 3.823 0.544; 3.05 7.9 1.1 1.8 1.225 0.25 0.45 4.153 0.881; 5.1 12.95 1.975 3.65 3.075 0.25 1.075 4.045 1.217];
eqn = A*X==B;
[As,Bs] = equationsToMatrix(eqn);
X_without_constraints = double(As)\double(Bs);
norm(double(As)*X_without_constraints-double(Bs))
Aeq = [ones(1,9) zeros(1,9);zeros(1,9) ones(1,9)];
beq = [1;1];
lb = zeros(18,1);
ub = Inf(18,1);
X_with_constraints = lsqlin(double(As),double(Bs),[],[],Aeq,beq,lb,ub);
norm(double(As)*X_with_constraints-double(Bs))
X = reshape(X_with_constraints,[9,2]).'
sum(X,2)
See Also
Categories
Find more on Linear Least Squares in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!