0

What is the general procedure in solving systems of equations using R (as opposed to manual Gauss-Jordan/Gaussian elimination)?

Must I first determine if the system is determined/under/overdetermined?

If a system is determined, I just use 

solve(t(a)%*%a)%*%t(a)%*%b

to get $x$ in $Ax = b$

If it overdetermined or underdetermined, I am not quite sure what to do. I think the above sometimes gives an answer depending on the rank, but the solution is not always unique. How can I get all the solutions? I think if there's no solution, will R just give an error?

Context: I am planning to recommend to my Stochastic Calculus professor that we use R in our upcoming exam (as opposed to tedious calculators/by-hand computation) so I have a feeling only simple functions will do (e.g. solve) for over/underdetermined systems rather than lengthy programs/functions.

Edit: I tried using solve(a,b), but I think that still doesn't give me all the solutions.

Here is an underdetermined example (R cannot give an answer since a is not square):

a=matrix(c(1,1,1,3,2,1),byrow=T,nrow=2)
a
b=matrix(c(1,2),byrow=T,nrow=2)
b
solve(a,b)
David Arenburg
  • 91,361
  • 17
  • 137
  • 196
BCLC
  • 139
  • 1
  • 3
  • 13
  • 1
    Please add an [tag:r] tag. – jotik Aug 12 '14 at 14:59
  • In the case you show you should be using `solve(a,b)` instead of what you are doing. – Bhas Aug 12 '14 at 15:10
  • @Bhas Sorry. I forgot to mention that. I think I tried that but it doesn't give me all the solutions. I think it just gives me the particular of infinitely many solutions that minimizes least squares error or something. – BCLC Aug 12 '14 at 16:22
  • @Bhas Updated. Using solve gives me no answer. Using the function in the link in the OP gives me A solution but not all solutions. – BCLC Aug 12 '14 at 16:40
  • `solve` gives you an error message. Matrix `a` is not square which it should be for `solve`. You'd better give more details of what you want. Have a look here [System of linear equations](http://en.wikipedia.org/wiki/System_of_linear_equations) especially the section `Matrix solution` for non square systems. – Bhas Aug 12 '14 at 17:40
  • 1
    @BCLC updated. For R you can look at package `limSolve` which provides a function `Solve` (note the uppercase S) for solving non square linear equation systems. Experiment with that. – Bhas Aug 12 '14 at 17:43
  • @Bhas As I feared, what limSolve gave me was one of the infinitely many solutions.* What I'm hoping to get is the infinitely many solutions represented parametrically which are [t,1-2t,t]. What limSolve/R gave me was the solution when t=1/3. Bhas, do you know anything in R that gives such parametric solutions? Thanks much for your replies, btw. ^_^ *> a=matrix(c(1,1,1,3,2,1),byrow=T,nrow=2) > b=matrix(c(1,2),byrow=T,nrow=2) > library(limSolve) > Solve(a,b) [1] 0.3333333 0.3333333 0.3333333 – BCLC Aug 12 '14 at 21:14

2 Answers2

4

The link I gave in section Matrix solution in the Wikipedia article on linear systems shows how to get what you want.
Define matrix A and vector b like this 

A <- matrix(c(1,1,1,3,2,1),byrow=T,nrow=2)
A
b <- matrix(c(1,2),byrow=T,nrow=2)
b

The following code will give you the general solution to your underdetermined system 

library(MASS)

Ag <- ginv(A)
Ag
xb <- Ag %*% b
xb
Aw <- diag(nrow=nrow(Ag)) - Ag %*% A
Aw

You can check that this is correct with 

w <- runif(3)
z <- xb + Aw %*% w
A %*% z - b

where the vector w is any arbitrary vector.
You can simplify the solution further manually to what you gave; I leave that as an exercise for you. As far as I know you can't get that solution automatically but maybe package Ryacas can do it.

You can get what you want by using package MASS or package pracma. E.g. with MASS: 

library(MASS)
N <- Null(t(A))

Then the solution is 

xb + N * q

where q is an arbitrary scalar.

With pracma:

N <- null(A)  # or nullspace(A)

with the same expression as above for the solution.

Bhas
  • 1,844
  • 1
  • 11
  • 9
4

Try qr.solve(A,b). That should work for both under- and over-determined systems.

Sheila
  • 2,438
  • 7
  • 28
  • 37