Correct and incorrect answer from linsolve

The following gives a correct answer:

A=[-2/5,1/5;2/5,-1/5;1,1]
B=[0;0;1]
linsolve(A,B)

 ans =0.3333 0.6667

The following, however, gives an incorrect answer:

C=[-0.2,0.3;0.25,-0.3;1,1]
D=[0;0;1]
C\D

 ans = 0.5699 0.4297

The correct answer, found by substitution, is 0.6000 0.4000

How can I get a correct answer by linsolve?

Answer 1

[Walter Roberson]

x = sym('x',[2 1]);
C*x

ans =
(3*x2)/10 - x1/5
x1/4 - (3*x2)/10
        x1 + x2

You have a system of 3 equations in 2 unknowns. It is overdetermined, and might not have any exact solution. The \ operation will do a least-squared fit to find an answer that is least bad in some sense.

Solving (3*x2)/10 - x1/5 = 0 for x1 gives x1 = (3*x2)/2. Substituting that back into C*x gives

         0
  (3*x2)/40
  (5*x2)/2

solving (3*x2)/40 = 0 for x2 gives x2 = 0. Substituting that back gives (5*0)/2 = 1 which is 0 = 1 which has no solution.

>> C*[0.6;0.4]
ans =
                         0
                      0.03
                         1

So 0.6 0.4 is not a solution after-all.

>> sum((D-C*(C\D)).^2)
ans =
      0.000407074042245239
>> sum((D-C*[0.6;0.4]).^2)
ans =
                    0.0009

so the solution found by C\D gives less of an error than [0.6 0.4] does.

[Oddur Bjarnason]

Thank you for answering my question. I come from the Mathematica World and am new to MATLAB. I am working hard at understanding your answer. After having looked at your MATLAB Answers I expect it to be right. By substitution I got the answer 0.6 and 0.4. Now you tell me that this is not an answer at all and that I should put my faith in the answer that I got from linsolve, that is 0.5699 0.4297. When I multiply the initial state vector by a high power of the transition matrix the vector seems to approach a stationary vector of [0.6,0.4].

[Steven Lord]

The last sentence of your comment makes no sense. The matrices you posted are non-square so you can't compute the matrix power -- it's not defined.

C=[-0.2,0.3;0.25,-0.3;1,1]
C^2 % throws an error

Using element-wise power doesn't work either, that doesn't converge to [0.6; 0.4].

Q = zeros(2, 20);
for k = 1:20
  Q(:, k) = (C.^k)\[0; 0; 1];
end
[1:20; Q]

Please explain in more detail the exact problem you're trying to solve and how you're trying to solve it in MATLAB and we may be better able to explain what's going on.

[Walter Roberson]

You have a system of 3 equations in 2 unknowns. Such systems only have a solution if at least one of the equations can be derived from the others, which is not the case here. There are no exact solutions to the system given: there are only solutions that are wrong. The \ operator returns the solution that is least wrong in some sense.

What it would take for a system with those C coefficients to be consistent would be if instead of D = [0; 0; 1], it were D = [3/110; 0 1], for which the solution would be [6/11, 5/11]

Alternately, if you want to hold D to [0; 0; 1], and you want [0.6 0.4] to be the solution, then one possible change would be to change C from

[ -1/5,  3/10]
[  1/4, -3/10]
[    1,     1]

to

[ -1/5, 3/10]
[  1/4, -3/8]
[    1,    1]

[Oddur Bjarnason]

Dear Walter and Steven

Thank you for your answer/comments.

A Markov chain has the initial state vector S0= [0.25,0.75] and the transition matrix T=[0.8,0.2;0.3,0.7].

The stationary state vector is approached by [0.25,0.75]*[0.8,0.2;0.3,0.7]^N if N is sufficiently high. For N is = 100 the stationary state vector is [0.6,0.4].

Let the stationary state vector be SS and let SS=[SS1,SS2]. The vector can be found by solving the following equation:

SS*T=SS

or

[SS1,SS2]*[0.8,0.2;0.3,0.7]=[SS1,SS2]

This yields the following system of linear equations:

 0.8*SS1+0.3*SS2=SS1
 0.25*SS1+0.7*SS2=SS2
 SS1+SS2=1

This equation can be solved by substitution and the stationary state vector SS = [SS1,SS2] = [0.6,0.4].

It should also be possible to solve the system of equations by linsolve. Before trying to do that SS1 must be subtracted from both sides of the first equation and SS2 from both sides of the second equation, yielding the following equations:

 -0.2*SS1+0.3*SS2=0
 0.25*SS1-0.3*SS2=0
 SS1+SS2=1

The coefficient matrix is

A=[-0.2,0.3;0.25,-0.3]

and the vector is

B=[0;0;1]

The solution by linsolve is

linsolve(A,B) = [0.5699,0.4279] Which is approximately "correct".

Can you help me to understand why [0.6,0.4] should not be a solution?

Steven asked what I was trying to do.

I am interested in developing mathematical models of migration. One of the methods I should like to use is Markov chains. I thought that it might be good to use MATLAB seing that it is based on matrices. I therefore got MATLAB,started learning it and began to try to use it for studying Markov chains. As you see I am still crawling.

Another reason, quite unrelated, for getting MATLAB, was that I bought a book called Cosmology with Matlab and thought that I should like to have access to MATLAB while reading it.

[Steven Lord]

Now the cause of the problem is clear. You wrote:

 This yields the following system of linear equations:

 0.8*SS1+0.3*SS2=SS1 0.25*SS1+0.7*SS2=SS2 SS1+SS2=1

You have a typo in your second equation. You wrote 0.25 instead of 0.2. So linsolve is solving the problem you told it to solve, not the problem you wanted it to solve. When I write the correct coefficient matrix, both linsolve and backslash (\) return the answer you expected.

A = [0.8-1 0.3; 0.2 0.7-1; 1 1];
b = [0; 0; 1];
x1 = linsolve(A, b)
x2 = A\b

Since you're interested in Markov chains you may find the PageRank and Markov Chains section of the chapter on Linear Equations in Cleve's Numerical Computing in MATLAB interesting and informative.

[Oddur Bjarnason]

Thank you Steven, When I was filling in the coefficient Matrix I read 0.2S1 to be 0.25. So you are right - it was a typo.