Chapter 10 - Curve fitting regression

Example No. 10_01 Pg No. 326

In [13]:
from numpy import array
from sympy.abc import x
#Fitting a Straight Line

X = range(1,6)
Y = array([[ 3, 4, 5 ,6 ,8 ]])
n = len(X)#
X=array(X)
b = ( n*sum(X*Y) - sum(X)*sum(Y) )/( n*sum(X*X) - (sum(X))**2 )
a = sum(Y)/n - b*sum(X)/n
print 'b = ',b
print 'a = ',a
y = a + b*x
print 'y=',y

b =  [-1 -1  0  0  1]
a =  [ 3  3  1  1 -2]
y= [-x + 3 -x + 3 1 1 x - 2]

Example No. 10_02 Pg No. 331

In [1]:
from numpy import array,log,exp
from sympy.abc import x

#Fitting a Power-Function model to given data

X = array(range(1,6))
Y = [ 0.5, 2 ,4.5 ,8 ,12.5 ]
Xnew = log(X)
Ynew = log(Y)
n = len(Xnew)
b = ( n*sum(Xnew*Ynew) - sum(Xnew)*sum(Ynew) )/( n*sum(Xnew*Xnew) - ( sum(Xnew) )**2 )
lna = sum(Ynew)/n - b*sum(Xnew)/n
a = exp(lna)
print 'b = ',b
print 'lna = ',lna
print 'a = ',a
print '\n The power function equation obtained is \n y = %.4Gx**%.4G'%(a,b)#

b =  2.0
lna =  -0.69314718056
a =  0.5

 The power function equation obtained is 
 y = 0.5x**2

Example No. 10_03 Pg No. 332

In [1]:
from numpy import array,log,exp
time = array(range(1,5))
T = [ 70 ,83 ,100 ,124 ]
t = 6
Fx = exp(time/4.0)
n = len(Fx)
Y = T #
b = ( n*sum(Fx*Y) - sum(Fx)*sum(Y) )/( n*sum(Fx*Fx) - (sum(Fx))**2 )
a = sum(Y)/n - b*sum(Fx)/n
print 'b = ',b
print 'a = ',a
print 'The relationship between T and t is \nT = %.4G*e**(t/4) + %.4G \n'%(b,a)
#deff('T = T(t)'%('T = b*exp(t/4) + a '
def T(t):
    tt=b*exp(t/4.0)+a
    return tt
                 
T_6 = T(6)
print 'The temperature at t = 6 is',T_6

b =  37.6294062985
a =  20.9234245534
The relationship between T and t is 
T = 37.63*e**(t/4) + 20.92 

The temperature at t = 6 is 189.566723485

Example No. 10_04 Pg NO. 335

In [1]:
from numpy import array,ones,identity
from numpy.linalg import inv

#Curve Fitting

x = array(range(1,5))
y = [6, 11, 18, 27 ]
n = len(x) #Number of data points
m = 2+1       #Number of unknowns
print 'Using CA = B form , we get'
C=identity(m)
B=ones(m)
for j in range(0,m):
    for k in range(0,m):
        C[j,k] = sum(x**(j+k-2))
    
    B[j] = sum( y*( x**( j-1 ) ) )

print 'B',B
print 'C',C
A = inv(C)*B
print 'A = ',A
print 'Therefore the least sqaures polynomial is\n   y = 1J + 1J*x + 1J*x**2 \n',(A[0])
print A[1]
print A[2]

Using CA = B form , we get
B [   6.   62.  190.]
C [[  1.   1.   4.]
 [  1.   4.  10.]
 [  4.  10.  30.]]
A =  [[  20.          103.33333333 -190.        ]
 [  10.          144.66666667 -190.        ]
 [  -6.          -62.           95.        ]]
Therefore the least sqaures polynomial is
   y = 1J + 1J*x + 1J*x**2 
[  20.          103.33333333 -190.        ]
[  10.          144.66666667 -190.        ]
[ -6. -62.  95.]

Example No. 10_05 Pg No. 342

In [1]:
from numpy import array,ones,identity,arange,vstack,transpose
from scipy.sparse.linalg import lsqr
#Plane Fitting

x = range(1,5)
z = range(0,4)
y = arange(12,31,6)
n = len(x) #Number of data points
m = 3         #Number of unknowns
G = vstack([ones(n),x,z])
H = transpose(G)
C = G.dot(H)
B = y.dot(H)
D = lsqr(C,B)
print 'C=\n',C
print 'B=\n',B
print '\n The regression plane is \n y = %d + %.f*x + %d*z \n'%(D[0][0],D[0][1],D[0][2])

C=
[[  4.  10.   6.]
 [ 10.  30.  20.]
 [  6.  20.  14.]]
B=
[  84.  240.  156.]

 The regression plane is 
 y = 5 + 6*x + 0*z