Normal Equation for linear regression

Question

I have the following X and y matrices:

for which I want to calculate the best value for theta for a linear regression equation using the normal equation approach with:

theta = inv(X^T * X) * X^T * y

the results for theta should be : [188.400,0.3866,-56.128,-92.967,-3.737]

I implement the steps with:

X=np.matrix([[1,1,1,1],[2104,1416,1534,852],[5,3,3,2],[1,2,2,1],[45,41,30,36]])
y=np.matrix([460,232,315,178])

XT=np.transpose(X)

XTX=XT.dot(X)

inv=np.linalg.inv(XTX)

inv_XT=inv.dot(XT)

theta=inv_XT.dot(y)

print(theta)

But I dont't get the desired results. Instead it throws an error with:

Traceback (most recent call last): File "C:/", line 19, in theta=inv_XT.dot(y) ValueError: shapes (4,5) and (1,4) not aligned: 5 (dim 1) != 1 (dim 0)

What am I doing wrong?

Answer 1

I have solved the problem by using numpy.linalg.pinv() that is the "pseudo-inverse" instead of numpy.linalg.inv() for the inverting of the matrix since the ducumentation sais:

"The pseudo-inverse of a matrix A, denoted A^+, is defined as: “the matrix that ‘solves’ [the least-squares problem] Ax = b,” i.e., if \bar{x} is said solution, then A^+ is that matrix such that \bar{x} = A^+b."

and solving the least-squares problem is exactly what I want to achieve in the context of linear regression.

Consequently the code is:

X=np.matrix([[1,2104,5,1,45],[1,1416,3,2,40],[1,1534,3,2,30],[1,852,2,1,36]])
y=np.matrix([[460],[232],[315],[178]])

XT=X.T
XTX=XT@X

inv=np.linalg.pinv(XTX)

theta=(inv@XT)@y
print(theta)

[[188.40031946]
 [  0.3866255 ]
 [-56.13824955]
 [-92.9672536 ]
 [ -3.73781915]]

Edit: There is also the possibility of regularization to get rid of the problem of NoN-invertibility by changing the normal-equation to:

theta = (XT@X + lambda*matrix)^(-1)@XT@y where lambda is a real number and called regularization parameter and matrix is a (n+1 x n+1) dimensional matrix of the shape:

 0 0 0 0 ... 0 0 
 0 1 0 0 ... 0 0 
 0 0 1 0 ... 0 0
 0 0 0 1 ... 0 0
 .
 .
 .
 0 0 0 0 0 0 0 1

That is a eye() matrix with the element [0,0] set to 0

More about the concept of regularization can be read here