import numpy as n
import numpy.linalg as la
import numpy.random as ran
import pylab as p

# Simple least-squares fit of data to a linear function.
x = n.arange(0, 10., 0.1)
y = x + ran.normal(0, 1., len(x))    # A little noise has been added.  Instead of x, try x**2.
# Fit y to the form mx + c.
# Cast the problem as y = ap.
# Here, a = [[x 1]]  is just a matrix with two columns, the set x and a set of 1's.
# The quantity p is a vector p = [[m], [c]], where [m] is a set of the slope m and [c] is a set of the constant c.
a = n.vstack([x, n.ones(len(x))]).T
[m, c], r, rank, s = la.lstsq(a, y)
print r     # sum of residues
# z is the functional fit to the data
z = m*x + c

p.plot(x,y)
p.plot(x,z)
p.show()