# Solving simultaneous linear equations with R

Am trying to solve a system of linear equations involving three variables using R. However, I get an error as shown below:

> A = matrix(data=c(-0.4, 0.4, 0.15, 0.3, -0.55, 0.6, 0.1, 0.15, -0.75), nrow=3, ncol=3, byrow=TRUE)
> A
[,1]  [,2]  [,3]
[1,] -0.4  0.40  0.15
[2,]  0.3 -0.55  0.60
[3,]  0.1  0.15 -0.75
>
> b = matrix(data=c(0, 0, 0), nrow=3, ncol=1, byrow=TRUE)
> b
[,1]
[1,]    0
[2,]    0
[3,]    0
>
> solve(A, b)
Error in solve.default(A, b) :
system is computationally singular: reciprocal condition number = 1.04615e-17


When I work them out using substitution method, I get the variable values as: 0.456, 0.4027, 0.1413 respectively (not exactly precise values). The sum of all the variables should be 1.

How could I solve such a system?

Background: These equations come from solving a Discrete Time Markov Chain. In reality, there could be around 20 such simultaneous equations on 20 variables.

• There are infinitely many solutions in this case. Check out A %*% c(1, .8837, .31) and A %*% c(3.225, 2.85, 1) for two examples. Commented Jun 26, 2012 at 19:45
• Ah, I missed out the constraint! Updating the post. Commented Jun 26, 2012 at 19:50

A1 = matrix(data=c(-0.4, 0.4, 0.15, 0.3, -0.55, 0.6, 0.1, 0.15, -0.75), nrow=3, ncol=3, byrow=TRUE)
A2<-eigen(A)\$vector
A%*%A2[,3]
[,1]
[1,]  2.461139e-17
[2,] -1.083796e-16
[3,]  3.932943e-17
1> A2[,3]
[1] 0.7298860 0.6450156 0.2263212


I suppose you mean the norm of the variables should be 1.

Otherwise:

1> A2[,3]/sum(A2[,3])
[1] 0.4558304 0.4028269 0.1413428

• +1 MASS::Null should also get this answer.
– whuber
Commented Jun 26, 2012 at 20:20
• This looks brilliant, although I'm not sure about the theory involved here. Actually, all the variables indicate certain probabilities (of being in a state) -- so their sum should be 1. Commented Jun 27, 2012 at 6:21
• @Barun: yes and, geometrically, this constraint implies that your LHS matrix is rank deficient and therefore that it's last eigenvalue is 0 and that the corresponding eigenvector gives (a basis for) the vector you are looking for (up to a rescaling/change of signs). Commented Jun 27, 2012 at 7:49