# Regularised regression code

I wrote some code to do regularised linear regression and it works, but I don't like the fact that I've had to double call the functions when plotting, nor the fact that I've sliced those calls to get the parts that I want.

import numpy as np
import matplotlib.pyplot as plt

def fit(phi_fn, xx, yy):
w_fit = np.linalg.lstsq(phi_fn(xx), yy, rcond=None)[0]
grid_size = 0.01
x_grid = np.arange(0,9,grid_size)[:,None]
f_grid = np.matmul(phi_fn(x_grid),w_fit)
return(x_grid, f_grid)

def fitreg(phi_fn, xx, yy,lamb):
yy = np.pad(yy,(0,8),'constant',constant_values=0)
zz = np.concatenate((phi_fn(xx),lamb*np.identity(8)),axis=0)
w_fit = np.linalg.lstsq(zz, yy, rcond=None)[0]
grid_size = 0.01
x_grid = np.arange(0,9,grid_size)[:,None]
f_grid = np.matmul(phi_fn(x_grid),w_fit)
return(x_grid, f_grid)

def phi_poly(xx):
return np.concatenate([np.ones((xx.shape[0],1)), xx,xx**2,xx**3,xx**4,xx**5,xx**6,xx**7], axis=1)

D = 1
N = 10
mu = np.array([0,1,2,3,4,5,6,7,8,9])

xx = np.tile(mu[:,None], (1, D)) + 0.01*np.random.randn(N, D)
yy = 2*xx + 2*np.random.randn(N,D)

plt.clf()
plt.plot(xx,yy,'kx')

plt.plot(fit(phi_poly, xx, yy)[0], fit(phi_poly, xx, yy)[1], 'b-')
plt.plot(fitreg(phi_poly, xx, yy,1)[0], fitreg(phi_poly, xx, yy,1)[1][:,0], 'r-')
plt.plot(fitreg(phi_poly, xx, yy,10)[0], fitreg(phi_poly, xx, yy,10)[1][:,0], 'g-')
plt.plot(fitreg(phi_poly, xx, yy,0.1)[0], fitreg(phi_poly, xx, yy,0.1)[1][:,0], 'y-')

plt.show()

• Did you test the code with the values you posted? When trying to run it, I currently get a TypeError: must be real number, not NoneType in the call to np.linalg.lstsq. Aug 26, 2018 at 11:39
• Yeah, it works fine for me. Maybe my version of numpy is old, lemme go update. Aug 26, 2018 at 12:00
• FYI, I used numpy 1.13.3 in cPython 3.6.3 Aug 26, 2018 at 12:09
• I'm using numpy 1.14.4 and python 3.6.5 Aug 26, 2018 at 12:13

## 1 Answer

You can get rid of some of the repetition using a for loop looping over the lambdas and colors and using tuple assignment to at least get rid of one level of having to use indexing as well as the necessity of having to call fitreg twice for each plot:

if __name__ == "__main__":

D = 1
N = 10
mu = np.arange(10)

xx = mu[:, None] + 0.01 * np.random.randn(N, D)
yy = 2 * xx + 2 * np.random.randn(N, D)

plt.clf()
plt.plot(xx, yy, 'kx')

x, y = fit(phi_poly, xx, yy)
plt.plot(x, y, 'b-')
for lamb, color in zip([1, 10, 0.1], "rgy"):
x, y = fitreg(phi_poly, xx, yy, lamb)
plt.plot(x, y[:, 0], '{}-'.format(color))

plt.show()


I also added a if __name__ == "__main__": guard to allow importing from this script in another script and used numpy.arange for your mu array.

For the construction of xx you don't actually need the numpy.tile, the normal broadcasting rules already do that for you.