How can I solve a system of many equations in 2 variables ?

Hello everyone,

I have a system of 12 equations (6pairs of equations) with 2 variables A and B.

I want to find a value of A and B that can fit all the equations (optimizing the result)

First, I've resolved each system of 2 equations, then I found 6 values of A and B. Next, I've calculated the sum over 6 (the average value).

The results is not really good (it does'nt fit all the equations.

So then I defined two variables only A and B, and in the function solve in MATLAB I inserted 12equations.

It seems that MATLAB cannot resolve the problem in this case, I got this message (S is empty, no solution):

Warning: 4 equations in 2 variables.

> In C:\Program Files\MATLAB\R2014a\toolbox\symbolic\symbolic\symengine.p>symengine at 56

In mupadengine.mupadengine>mupadengine.evalin at 97

In mupadengine.mupadengine>mupadengine.feval at 150

In solve at 170

In systeme_equation at 64

Warning: Explicit solution could not be found.

> In solve at 179

In systeme_equation at 64

S =

[ empty sym ]

The aim is not to find a result (A and B) that satisfy all the equations, but at least a results that approaches the solutions of each systedm of 2 equations.

Do you have an idea about how to resolve this kind of problems with MATLAB ?

syms A B

S = solve(0.3*A + 0.3*B == 1.47e13, 0.3*B - 0.3*A == 8.3e12, 0.2*A + 0.2*B == 9.9e12, 0.2*B - 0.2*A == 5.3e12, 0.1*A + 0.1*B == 5.2e12, 0.1*B - 0.1*A == 2.5e12, -0.1*A - 0.1*B == -5.7e12, 0.1*A - 0.1*B == -2.1e12, -0.2*A - 0.2*B == -1.25e13, 0.2*A - 0.2*B == -3.9e12, -0.3*A - 0.3*B == -2.05e13, -0.3*B + 0.3*A == -5.5e12 );

Answer 1

[John D'Errico]

I think you don't understand the concept of linear least squares, or linear regression.

You have essentially a problem with 12 equations, in 2 unknowns. Usually, there will be no solution that satisfies all of the equations at once. There is no need to use symbolic methods to solve the problem either. As well, you need to learn to use arrays in MATLAB.

rhs = [1.47e13; 8.3e12; 9.9e12; 5.3e12; 5.2e12; 2.5e12; -5.7e12; -2.1e12; -1.25e13; -3.9e12; -2.05e13; -5.5e12];
M = [0.3 0.3;
     0.3 -0.3;
     0.2 0.2;
     0.2 -0.2;
     0.1 0.1;
     0.1 -0.1;
     -0.1 -0.1;
     0.1 -0.1;
     -0.2 -0.2;
     0.2 -0.2;
     -0.3 -0.3;
     -0.3 0.3];

I've created a vector rhs, and an array M.

To convince you this is the same problem, you might have written:

syms A B
vpa(M*[A;B] == rhs)
ans =
    0.3*A + 0.3*B == 14700000000000.0
     0.3*A - 0.3*B == 8300000000000.0
     0.2*A + 0.2*B == 9900000000000.0
     0.2*A - 0.2*B == 5300000000000.0
     0.1*A + 0.1*B == 5200000000000.0
     0.1*A - 0.1*B == 2500000000000.0
  - 0.1*A - 0.1*B == -5700000000000.0
    0.1*A - 0.1*B == -2100000000000.0
 - 0.2*A - 0.2*B == -12500000000000.0
    0.2*A - 0.2*B == -3900000000000.0
 - 0.3*A - 0.3*B == -20500000000000.0
    0.3*B - 0.3*A == -5500000000000.0

and you would see we have recovered the original set ot equations, only now in a simpler to use form.

However, now you want to solve for the unknown variables that minimizes the residual errors for the problem

M*X = rhs

where X is a 2x1 vector. Thus, we wish to find X, such that norm(M*X - rhs) is as small as possible. Typically the 2-norm is used. The solution in MATLAB is just

AB = M\rhs
AB =
          36767857142857.1
          20839285714285.7

So we would have the optimal values of A and B as

A = AB(1);
B = AB(2);

How well did we do? Again, remember this is an overall optimum set for A and B.

format short g
[rhs,M*AB, rhs - M*AB]/1e11
ans =
          147       172.82      -25.821
           83       47.786       35.214
           99       115.21      -16.214
           53       31.857       21.143
           52       57.607      -5.6071
           25       15.929       9.0714
          -57      -57.607      0.60714
          -21       15.929      -36.929
         -125      -115.21      -9.7857
          -39       31.857      -70.857
         -205      -172.82      -32.179
          -55      -47.786      -7.2143

Could be worse. The first column is actual, then predicted, and finally the error. I've divided by 1e11 to make the numbers easier to read. Or, I could just have plotted the results as:

plot(1:12,[rhs,M*AB],'-o')
legend('actual','predicted')