Python Linear Equations - Gaussian Elimination

Question

Goal

Given a set of points, I'm trying to find the coefficients of the linear equation that satisfies all the points provided.

For example, if I wanted to find the linear equation (ax + by + c = z):

3x + 2y + 2 = z

I would need a minimum of three, three-dimensional points:

(2, 2, 12)
(3, 4, 19)
(4, 5, 24)

Given enough points with coordinates (x, y, z), I should be able to find (a, b, c) using Gaussian Elimination.

However, I think I'm having issues back solving the matrix in special cases. You can view my first stab at a python implementation here: https://gist.github.com/anonymous/8188272

Let's take a look at a few examples...

Data Set 1

Use the following "hand made" points (x, y, z):

(2, 2, 12)
(3, 4, 19)
(4, 5, 24)

Perform LU decomposition on the following matrix:

[[  2.   2.   1.  12.]
 [  3.   4.   1.  19.]
 [  4.   5.   1.  24.]]

Back solve the U matrix:

[[  4.    5.    1.   24. ]
 [  0.   -0.5   0.5   0. ]
 [  0.    0.    0.5   1. ]]

The returned result (a, b, c):

[3.0, 2.0, 2.0]

Correct! Everything seems fine...

Data Set 2

Use the following "hand made" points (x, y, z):

(3, 4, 19)
(4, 5, 24)
(5, 6, 29)

Perform LU decomposition on the following matrix:

[[  3.   4.   1.  19.]
 [  4.   5.   1.  24.]
 [  5.   6.   1.  29.]]

Back solve the U matrix:

[[  5.00000000e+00   6.00000000e+00   1.00000000e+00   2.90000000e+01]
 [  0.00000000e+00   4.00000000e-01   4.00000000e-01   1.60000000e+00]
 [  0.00000000e+00   0.00000000e+00   4.44089210e-16   0.00000000e+00]]

The returned result (a, b, c):

[1.0, 4.0, 0.0]

While technically it is a solution, not what I was looking for!

Data Set 3

Use the following "hand made" points (x, y, z):

(5, 6, 29)
(6, 7, 34)
(7, 8, 39)

Perform LU decomposition on the following matrix:

[[  5.   6.   1.  29.]
 [  6.   7.   1.  34.]
 [  7.   8.   1.  39.]]

Back solve the U matrix:

[[  7.00000000e+00   8.00000000e+00   1.00000000e+00   3.90000000e+01]
 [  0.00000000e+00   2.85714286e-01   2.85714286e-01   1.14285714e+00]
 [  0.00000000e+00   0.00000000e+00   0.00000000e+00   3.55271368e-15]]

Implementation crashes...

Thoughts

In Data Set 2 and 3 the last row and second to last row are "special". The second to last row has the same value for "b" and "c" (which is true in my special example!). Unfortunately, I lack the mathematical knowledge to make heads or tails from it.

Is there something special case I need to handle when the last row is all zeros and the row above it has equal values across?

Thanks in advance!

Answer 1

Yes, it is a special case you need to handle differently. In cases 2 and 3 you have a rank deficient matrix. Generally it can mean there are infinitely many solutions, or no solutions.

You can identify whether these cases are going to occur by checking the determinant of the matrix you make by stacking those 3-vectors.

>>> import numpy as np
>>> from scipy.linalg import det
>>> data1 = np.array([(2, 2, 12), (3, 4, 19), (4, 5, 24)])
>>> data2 = np.array([(3, 4, 19), (4, 5, 24), (5, 6, 29)])
>>> data3 = np.array([(5, 6, 29), (6, 7, 34), (7, 8, 39)])
>>> det(data1)
-1.9999999999999982
>>> det(data2)
5.551115123125788e-17
>>> det(data3)
8.881784197001213e-16

Example 1 was a full rank matrix, which geometrically tells you that the 3 points are linearly independent.

Examples 2 and 3 make matrices with zero determinant, which tells you that the points are linearly dependent.