# Penalised Regression

## YouTube Videos
1. **Scikit Learn Linear Regression:** https://www.youtube.com/watch?v=EvnpoUTXA0E
2. **Scikit Learn Linear Penalise Regression:** https://www.youtube.com/watch?v=RhsEAyDBkTQ

## Introduction
We often do not want the coefficients/ weights to be too large. Hence we append the loss function with a penalty function to discourage large values of $w$.

\begin{align}
\mathcal{L} & = \sum_{i=1}^N (y_i-f(x_i|w,b))^2 + \alpha \sum_{j=1}^D w_j^2 + \beta \sum_{j=1}^D |w_j|
\end{align}
where, $f(x_i|w,b) = wx_i+b$. The values of $\alpha$ and $\beta$ are positive (or zero), with higher values enforcing the weights to be closer to zero.

## Lesson Structure
1. The task of this lesson is to infer the weights given the data (observations, $y$ and inputs $x$).
2. We will be using the module `sklearn.linear_model`.


```python
import numpy as np
from sklearn.linear_model import LinearRegression

import matplotlib.pyplot as plt
%matplotlib inline

# In order to reproduce the exact same number we need to set the seed for random number generators:
np.random.seed(1)
```

A normally distributed random looks as follows:


```python
e = np.random.randn(10000,1)
plt.hist(e,100) #histogram with 100 bins
plt.ylabel('y')
plt.xlabel('x')
plt.title('Histogram of Normally Distributed Numbers')
plt.show()
```

Generate observations $y$ given feature (design) matrix $X$ according to:
$$
y = Xw + \xi\\
\xi_i \sim \mathcal{N}(0,\sigma^2)
$$

In this particular case, $w$ is a 100 dimensional vector where 90% of the numbers are zero. i.e. only 10 of the numbers are non-zero.


```python
# Generate the data
N = 40 # Number of observations
D = 100 # Dimensionality

x = np.random.randn(N,D) # get random observations of x
w_true = np.zeros((D,1)) # create a weight vector of zeros
idx = np.random.choice(100,10,replace=False) # randomly choose 10 of those weights
w_true[idx] = np.random.randn(10,1) # populate then with 10 random weights

e = np.random.randn(N,1) # have a noise vector
y = np.matmul(x,w_true) + e # generate observations

# create validation set:
N_test = 50
x_test = np.random.randn(50,D)
y_test_true = np.matmul(x_test,w_true)
```


```python
model = LinearRegression()
model.fit(x,y)

# plot the true vs estimated coeffiecients
plt.plot(np.arange(100),np.squeeze(model.coef_))
plt.plot(np.arange(100),w_true)
plt.legend(["Estimated","True"])
plt.title('Estimated Weights')
plt.show()
```

One way of testing how good your model is to look at metrics. In the case of regression Mean Squared Error (MSE) is a common metric which is defined as:
$$ \frac{1}{N}\sum_{i=1}^N \xi_i^2$$ where, $\xi_i = y_i-f(x_i|w,b)$. Furthermore it is best to look at the MSE on a validation set, rather than on the training dataset that we used to train the model.


```python
y_est = model.predict(x_test)
mse = np.mean(np.square(y_test_true-y_est))
print(mse)
```

    6.79499954951


Ridge regression is where you penalise the weights by setting the $\alpha$ parameter right at the top. It penalises it so that the higher **the square of the weights** the higher the loss.


```python
from sklearn.linear_model import Ridge

model = Ridge(alpha=5.0,fit_intercept = False)
model.fit(x,y)

# plot the true vs estimated coeffiecients
plt.plot(np.arange(100),np.squeeze(model.coef_))
plt.plot(np.arange(100),w_true)
plt.legend(["Estimated","True"])
plt.show()
```

This model is slightly better than without any penalty on the weights.


```python
y_est = model.predict(x_test)
mse = np.mean(np.square(y_test_true-y_est))
print(mse)
```

    6.42288072501


Lasso is a model that encourages weights to go to zero exactly, as opposed to Ridge regression which encourages small weights.


```python
from sklearn.linear_model import Lasso

model = Lasso(alpha=0.1,fit_intercept = False)
model.fit(x,y)

# plot the true vs estimated coeffiecients
plt.plot(np.arange(100),np.squeeze(model.coef_))
plt.plot(np.arange(100),w_true)
plt.legend(["Estimated","True"])
plt.title('Lasso regression weight inference')
plt.show()
```

The MSE is significantly better than both the above models.


```python
y_est = model.predict(x_test)[:,None]
mse = np.mean(np.square(y_test_true-y_est))
print(mse)
```

    2.30600134084


Automated Relevance Determination (ARD) regression is similar to lasso in that it encourages zero weights. However, the advantage is that you do not need to set a penalisation parameter, $\alpha$, $\beta$ in this model.


```python
from sklearn.linear_model import ARDRegression

model = ARDRegression(fit_intercept = False)
model.fit(x,y)

# plot the true vs estimated coeffiecients
plt.plot(np.arange(100),np.squeeze(model.coef_))
plt.plot(np.arange(100),w_true)
plt.legend(["Estimated","True"])
plt.show()
```


```python
y_est = model.predict(x_test)[:,None]
mse = np.mean(np.square(y_test_true-y_est))
print(mse)
```

    2.87016741296


### Note:
Rerun the above with setting N=400

## Inverse Problems
The following section is optional and you may skip it. It is not necessary for understanding Deep Learning.

Inverse problems are where given the outputs you are required to infer the inputs. A typical example is X-rays. Given the x-ray sensor readings, the algorithm needs to build an image of an individuals bone structure.

See [here](http://scikit-learn.org/stable/auto_examples/applications/plot_tomography_l1_reconstruction.html#sphx-glr-auto-examples-applications-plot-tomography-l1-reconstruction-py) for an example of l1 reguralisation applied to a compressed sensing problem (has a resemblance to the x-ray problem). 


```python

```

fb5bc37188a96c164c15e44283f685de29d93f14
ipynb
Jupyter Notebook
deepschool.io/Lesson 01 - PenalisedRegression - Solutions.ipynb
gopala-kr/ds-notebooks
bc35430ecdd851f2ceab8f2437eec4d77cb59423
MIT
2018-03-28T09:08:02.000Z
2021-11-17T11:15:38.000Z
Lesson 01 - PenalisedRegression - Solutions.ipynb
OWLYone/deepschool.io
ae6718fc14f3ac499697c97edc97a66dad9d9a6c
Apache-2.0
Lesson 01 - PenalisedRegression - Solutions.ipynb
OWLYone/deepschool.io
ae6718fc14f3ac499697c97edc97a66dad9d9a6c
Apache-2.0
2019-03-24T19:29:08.000Z
2019-07-24T13:38:40.000Z
Qwen/Qwen-72B
1. YES
2. YES