For example, I am trying to perform linear regression on the following set of data

Data examples: $X = [[1, 20], [3, 40], [5, 60]]$ (each row is an example, there are three examples, each with a feature of $2$, arranged in Numpy array)

Targets: $y = [1, 2, 3]$ (whatever you like, it doesn't affect our result.

Fitting a standardscaler gives me,

X = [[1, 20], [3, 40], [5, 60]]
scaler = StandardScaler()
Y = scaler.transform(X)

$Y = [[-1.22474487 -1.22474487] [ 0. 0. ] [ 1.22474487 1.22474487]]$

Now I want to compute the normal equation of a linear regression problem. This inolves calculating the following matrix $Z = (Y^T Y)^{-1} Y^T$

Z = np.linalg.inv(np.dot(np.transpose(Y), Y))*np.transpose(Y)

I get LinAlgError: Singular matrix

Note that this does not seem to be a problem with the original data set $X$

Is this a usual behavior or did I do something wrong?

  • $\begingroup$ Since you have two features, YtY should be a 2×2 matrix and not 3×3 matrix. Try (YYt)inv Y instead. $\endgroup$ Jun 4, 2021 at 5:14

3 Answers 3


I think the root confusion is the nuance between linear and affine relationships, which is not something that becomes a problem in most of data science (we generally allow affine relationships even if we use the word "linear").

The matrix $X$ has full rank: the columns demonstrate an affine relationship ($x_2=10x_1+10$), but not a linear one. So $X^T X$ (which is $2\times2$) is indeed invertible, and everything proceeds normally.

If you add an all-ones column to $X$ (to incorporate an intercept to the OLS), you elevate the affine relationship to a linear one, and you'll find that $X^T X$ is not invertible.

The StandardScaler (in addition to scaling) centers the features, which again rips away the bias/shift, and turns the affine relationship to a linear one (of course, it's the identity relationship).


A linear or affine (with fixed constant term) transformation like scaling should not change the rank of a matrix nor its invertibility. This is one point. However according to documentation of StandardScaler, StandardScaler scales each feature independantly (ie. no longer a uniform linear/affine transformation, constant terms are multiple), so this can in fact alter the rank and invertibility of a matrix of data (for example use a different scaler).

A further point is that the formula for Moore-Penrose pseudo-inverse you use is valid only in certain cases (ie. for linearly independent columns) else a more general formula should be used. Check that as well.


There is a problem with your example data in the first place (not related to scaling). It produces a prefect fit. This will also translate to scaled data. You also need to add an constant term in order to solve $(X'X)^{-1} X'y = \hat{\beta}$.

Example in R:

df = data.frame(y=c(1,2,3),x1=c(1,3,5),x2=c(20,40,60)) summary(lm(y~x1+x2,data=df))

Warning message: In summary.lm(lm(y ~ x1 + x2, data = df)) :
essentially perfect fit: summary may be unreliable

x0 <- c(1,1,1) 
x1 <- c(1,3,5)
x2 <- c(20,40,60)
x <- as.matrix(cbind(x0,x1,x2))
y <- as.matrix(c(1,2,3))

# (X'X)^-1 X'y
beta1 = solve(t(x)%*%x) %*% t(x)%*%y 

Error in solve.default(t(x) %*% x) : Lapack routine dgesv: system is exactly singular: U[3,3] = 0


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.