# Is this a mistake in the feature standardization of ridge regression in book Machine Learning in Action?

I am new to ML. I am currently reading the classic book Machine Learning in Action by Peter Harrington. In its implementation of ridge regression in P165 Listing 8.3, the book standardizes feature matrix $X$ by subtracting means of attributes and then dividing it by variances of attributes not by standard deviations of attributes! As follows:

def ridgeRegres(xMat,yMat,lam=0.2):
xTx = xMat.T*xMat
denom = xTx + eye(shape(xMat)[1])*lam
if linalg.det(denom) == 0.0:
print "This matrix is singular, cannot do inverse"
return
ws = denom.I * (xMat.T*yMat)
return ws

def ridgeTest(xArr,yArr):
xMat = mat(xArr); yMat=mat(yArr).T
yMean = mean(yMat,0)
yMat = yMat - yMean     #to eliminate X0 take mean off of Y
#regularize X's
xMeans = mean(xMat,0)   #calc mean then subtract it off
xVar = var(xMat,0)      #calc variance of Xi then divide by it
xMat = (xMat - xMeans)/xVar
numTestPts = 30
wMat = zeros((numTestPts,shape(xMat)[1]))
for i in range(numTestPts):
ws = ridgeRegres(xMat,yMat,exp(i-10))
wMat[i,:]=ws.T
return wMat


Then calls ridgeTest by

>>> abX,abY=regression.loadDataSet('abalone.txt')
>>> ridgeWeights=regression.ridgeTest(abX,abY)


I do not think this makes any sense. According to standardization, To create a unit variance, we should divide $X-\mu$ by standard deviation $\sigma$, because: $$D(\frac{X-\mu}{\sigma})=\frac{1}{\sigma^2}D(X)=\frac{1}{\sigma^2}\sigma^2=1,$$ where $D(X)=\sigma^2$ is the variance of $X$.

If according to the book, divide $X-\mu$ by variance $\sigma^2$, we only get: $$D(\frac{X-\mu}{\sigma^2})=\frac{1}{\sigma^4}D(X)=\frac{1}{\sigma^4}\sigma^2=\frac{1}{\sigma^2},$$ which is not a unit variance.

The reason why I am not so confident that this is a mistake is that, first of all, the book is handling $Y$ carefully by subtracting its mean and not include an all 1 column in $X$ passed, both of which I find to be proper (see L2-normalization does not punish intercept). Hence I do not believe the book will make an obvious mistake like this. Secondly, nor the errata (see errata) or the book forum discussed it before. It is not a new book.

Although I posted divide by std not variance in the book forum, no reply has been received yet. So I turn to stack exchange. Does the feature standardization of the book mistakenly divide $X-\mu$ by variance or it is a meaningful action?

Thanks a lot in advance! Any clue will help me.

• There are different ways to scale your data. The one you describe (subtract mean, divide by std dev) is arguably the most common one, called z-scaling. If one divides by variance, as you point out they do not achieve unit variance across all features, so you end up having mean-centred scaled data where each feature has a different variance. – Zhubarb Jan 9 '18 at 11:08

Yes, it seems to be a mistake and may have an imporant impact on the way some variables are penalized.

If you have two variables with different magnitude, this rescaling with basically give a low standard deviation to the value with the (initial) largest standard deviation. Therefore, this variable will suffer a smaller penalization than it should. And this will make your regression depend on the unit in which the values are expressed (which is not a desirable feature).

Side note

The algorithm scales the input, but returns the outputs evaluated on the scaled inputs. This makes it quite hard to implement a "predict" function...

Numerical experiment

I fixed the issue (in the case where there is not intercept) so that the algorithm returns wMat[i,:]=ws.T / xVar.

import numpy as np
import matplotlib.pyplot as plt

def ridgeRegres(xMat,yMat,lam=0.2):
xTx = xMat.T*xMat
denom = xTx + np.eye(np.shape(xMat)[1])*lam
if np.linalg.det(denom) == 0.0:
print "This matrix is singular, cannot do inverse"
return
ws = denom.I * (xMat.T*yMat)
return ws

# Assuming there is not intercept : y = x_0 + 2 * x_1

def ridgeTest(xArr,yArr):
xMat = np.mat(xArr); yMat=np.mat(yArr).T
yMean = np.mean(yMat,0)
yMat = yMat - yMean     #to eliminate X0 take mean off of Y
#regularize X's
xMeans = np.mean(xMat,0)   #calc mean then subtract it off
xVar = np.var(xMat,0)      #calc variance of Xi then divide by it
xMat = (xMat - xMeans)/xVar
numTestPts = 30
wMat = np.zeros((numTestPts,np.shape(xMat)[1]))
for i in range(numTestPts):
ws = ridgeRegres(xMat,yMat,np.exp(i-10))
wMat[i,:]=ws.T / xVar
return wMat

X = np.mat(np.random.rand(50,2))
y = np.mat([1,2]) * X.T

res_orig = ridgeTest(X,y)

plt.plot(res_orig)
plt.title('Original Scaling')
plt.ylabel('Values of the (estimated) coefficient.')
plt.xlabel('Value of the penalization.')
plt.show()

X[:,1] *= 10
y = np.mat([1,2]) * X.T
res = ridgeTest(X,y)

plt.plot(res)
plt.title('x_1 has been multiplied by 10')
plt.ylabel('Values of the (estimated) coefficient.')
plt.xlabel('Value of the penalization.')
plt.show()


Result

The result with the original scaling provide the actual values ($y=x_0+2x_1$) whereas the second case exhibits a strange behavior.

Edit

As stated by @Berkmeister in the comments, this scaling can be defended : in many models, you can decide, arbitrarly, to change the weights of your features.

But, the fact that the coefficients returned by the model are not rescaled, and that this choice of (uncommon) scaling is proposed without prior explanations let me think that this code can be improved.

• Thanks a lot! You not only tell me this is a mistake, but also prove it experimentally. I learn a lot about how to test a model from this answer. – Naomi Jan 9 '18 at 13:54
• Machine Learning in Action should not have mistakes like this. Since it has many readers and most of them (like me) are new to ML and can not tell right from wrong. The influence can be very bad. I just updated my post about this problem in the book forum, hoping the author will notice and add it to the errata. – Naomi Jan 9 '18 at 14:00
• @Naomi thank you, I fixed the index issue. In case of doubt, you can also compare it to other widespread implementations (sk-learn has a ridge regression implementation and R has the famous glmnet package, per example), but the other caveat is that you need to make sure that parameters are expressed in a similar way... – RUser4512 Jan 9 '18 at 14:12