Clear Filters
Clear Filters

How do I solve over determined system of non linear equation?

6 views (last 30 days)
Suppose I have three matrices A,B,X. A and B being now and size 200x4. X being unknown with size 3x4. The system is AX=B.
syms x1 x2 x3 x4 y1 y2 y3 y4 z1 z2 z3 z4
A = rand(200,4);
X = [ x1 x2 x3 x4; y1 y2 y3 y4; z1 z2 z3 z4];
B = rand(200,4);
I want solve for X while keeping x1 x2 x3 y1 y2 y3 z1 z2 z3 values constrained to [-1,1]. What is the best approach?
  2 Comments
John D'Errico
John D'Errico on 27 May 2018
Edited: John D'Errico on 27 May 2018
I'm not sure why you think it is nonlinear. Purely linear. Admittedly, you have a size problem, so matrix multiplication does not conform.
Do I understand that you want the last row of the matrix X to be [0 0 0 1]?
Ahmad Adee
Ahmad Adee on 27 May 2018
Edited: Ahmad Adee on 27 May 2018
I agree to linear part (the way they are given) but in the original problem 3x3 of A and B are function of trigonometric identities. For X is 3x4 unknown and if you want to settle is to 4x4, 4th row can be [0 0 0 1], a constant one. and similar can be done to random data, all 4th and multiple of 4 rows can are similar and known [0 0 0 1].

Sign in to comment.

Accepted Answer

John D'Errico
John D'Errico on 27 May 2018
Edited: John D'Errico on 27 May 2018
Why do people think they need to use syms, just because they don't know the value of a result? (I think that somehow comes from the design of the symbolic toolbox and how it interacts with MATLAB.) Basic linear algebra suffices here. And you CAN solve the problem. You just need to accept that since it is over-determined, there is no exact solution in general.
A = rand(200,4);
B = rand(200,4);
Now, you want to solve the problem A*X=B, where X has a fixed, known row. Or apparently it may even be multiple rows in your actual problem. So I'll show how to solve it for an array of size 4x4, but with some set of known rows.
fixedRowIndex = 4;
X_fixed = [0 0 0 1];
estRowIndex = setdiff(1:size(A,2),fixedRowIndex);
Now, we can imagine partitioning the matrix X into two sub-matrices, they will look sort of like this:
% X = [X_est;X_fixed]
A is a 200x4 array. But when we multiply it by X, the problem splits into two pieces.
% A*X = A(:,estRowIndex)*X_est + A(:,fixedRowIndex)*X_fixed
So, does this make sense? If so, then you need to recognize that the second term is fully KNOWN. So we can move it to the right hand side. The problem reduces to solving for the unknown 3x4 array Xest:
% A(:,estRowIndex)*X_est = B - A(:,fixedRowIndex)*X_fixed
Now that we have eliminated those known rows from X, solve it simply, using backslash.
X_est = A(:,estRowIndex)\(B - A(:,fixedRowIndex)*X_fixed)
X_est =
0.45678 0.3027 0.2734 0.027561
0.25445 0.39719 0.3131 0.0085818
0.23474 0.28378 0.32242 -0.067524
X(estRowIndex,:) = X_est;
X(fixedRowIndex,:) = X_fixed
X =
0.45678 0.3027 0.2734 0.027561
0.25445 0.39719 0.3131 0.0085818
0.23474 0.28378 0.32242 -0.067524
0 0 0 1
The result is the least squares solution, consistent with the constraint that X has some fixed, known in advance rows. The nice thing is as I wrote it, if you have multiple fixed, known rows, it is trivially extensible.
Could I have done this using lsqlin? Yes, of course. But why bother?
  2 Comments
Ameer Hamza
Ameer Hamza on 27 May 2018
Edited: Ameer Hamza on 27 May 2018
This solution should be fine for a non-constrained case. But the question required that element of X are constrained to [-1, 1]. For example
A = rand(200,4);
B = 10*rand(200,4);
fixedRowIndex = 4;
X_fixed = [0 0 0 1];
estRowIndex = setdiff(1:size(A,2),fixedRowIndex);
X_est = A(:,estRowIndex)\(B - A(:,fixedRowIndex)*X_fixed);
X(estRowIndex,:) = X_est;
X(fixedRowIndex,:) = X_fixed
X =
3.7302 2.2907 3.1585 3.1444
1.9967 3.4358 2.7390 3.2424
3.4536 2.6884 3.0786 1.9789
0 0 0 1.0000
The X elements no longer belong to [-1, 1].
John D'Errico
John D'Errico on 27 May 2018
Edited: John D'Errico on 27 May 2018
Sigh. Still trivial. Of course, lsqlin will now be needed. Must I point out that you first wanted to "solve" the problem using X = A\B?
So how do we EFFICIENTLY solve this problem using LSQLIN and bound constraints?
LSQLIN applies to linear problems where you have bound constraints. Of course, LSQLIN will only solve one right hand side at a time. And really this problem has 4 right hand sides. backslash solved that in one line, but really, it was 4 basically parallel problems.
So the answer is one of two solutions. Simplest is to loop. This was the original solve.
X_est = A(:,estRowIndex)\(B - A(:,fixedRowIndex)*X_fixed)
But, IF we wanted to solve using LSQLIN and bound constraints, with a loop...
options = optimset('lsqlin');
options.Display = 'none';
Bhat = B - A(:,fixedRowIndex)*X_fixed;
A_est = A(:,estRowIndex);
nrhs = size(Bhat,2);
X = zeros(size(A,2),size(B,2));
X(fixedRowIndex,:) = X_fixed;
lb = repmat(-1,size(estRowIndex));
ub = repmat(1,size(estRowIndex));
for i = 1:nrhs
X(estRowIndex,i) = lsqlin(A_est,Bhat(:,i),[],[],[],[],lb,ub,[],options);
end
No use of fmincon needed.
Finally, could I have re-written this to avoid the loop completely? Well, yes, all in one call to lsqlin. The solution would be somewhat less efficient unless it is done carefully though.

Sign in to comment.

More Answers (2)

Ameer Hamza
Ameer Hamza on 27 May 2018
Edited: Ameer Hamza on 27 May 2018
Since the system is Overdetermined, the unique solution is highly unlikely. The best you can expect is the solution to the following least squares minimization problem
argmin_X mean2(abs(A*X-B))
You can get it like this
A = rand(200,4);
B = rand(200,4);
X = A\B
Edit: since you mentioned other constraints to be fulfilled, therefore you can solve the modified problem as follow
A = rand(200,4);
B = rand(200,4);
f = @(x) sum(sum((A*[x; 0 0 0 1]-B).^2)); % constraining the last row of x to be [0 0 0 1]
X = fmincon(f, zeros(3, 4), [], [], [], [], -1*ones(3, 4), ones(3, 4));
X = [X; 0 0 0 1];
This will force the last row of X to be [0 0 0 1] and also put constraints on other elements of X to be between -1 to 1.
Note the updated objective function use element-wise norm, which is same as used by MATLAB's lsqlin, thereby producing the same solution.
  3 Comments
Ameer Hamza
Ameer Hamza on 27 May 2018
The objective function of this fully linear least square problem is non-linear. I agree with your comment that lsqlin() is a possible solution, but fmincon is solving the same problem in a more compact and intuitive way to solve a non-linear objective function.
John D'Errico
John D'Errico on 28 Feb 2021
But the point is, the OBJECTIVE function is not nonlinear. It is quadratic. Using a nonlinear solver will force the tool to differentiate the problem, and then iterate to find a solution on a problem that requires no iteration. In the end, you will get a result that will be less accurate, since this is only an iterative method that will try then to converge, but only to within a tolerance.

Sign in to comment.


Walter Roberson
Walter Roberson on 26 May 2018
That is not possible. If A is something by 4, and X is 3 by something, then A*X would be trying to do algebraic matrix multiplication between arrays of inconsistent sizes.
In order to be able to calculate A*X with A being something by 4, then X would have to be 4 by something. If the result is to be 200x4 with A being 200 x 4, then X would have to be 4 x 4, definitely not 3 x 4.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!