# ML-Fundamentals - Lineare Regression - Exercise: Multivariate Linear Regression

## Table of Contents
* [Introduction](#Introduction)
* [Requirements](#Requirements) 
  * [Knowledge](#Knowledge) 
  * [Modules](#Python-Modules)
* [Exercises - Multivariate Linear Regression](#Exercises---Multivariate-Linear-Regression)
  * [Create Features](#Create-Features)
  * [Linear Hypothesis](#Linear-Hypothesis)
  * [Generate Target Values](#Generate-Target-Values)
  * [Plot The Data](#Plot-The-Data)
  * [Cost Function](#Cost-Function)
  * [Gradient Descent](#Gradient-Descent)
  * [Training and Evaluation](#Training-and-Evaluation)
  * [Feature Scaling](#Feature-Scaling)
* [Summary and Outlook](#Summary-and-Outlook)
* [Literature](#Literature) 
* [Licenses](#Licenses)

## Introduction

In this exercise you will implement the _multivariate linear regression_, a model with two or more predictors and one response variable (opposed to one predictor using univariate linear regression). The whole exercise consists of the following steps:

1. Generate values for two predictors/features $(x_1, x_2)$
2. Implement a linear function as hypothesis (model) 
3. Generate values for the response (Y / target values)
4. Plot the $((x_1, x_2), y)$ values in a 3D plot.
5. Write a function to quantify your model (cost function)
6. Implement the gradient descent algorithm to train your model (optimizer) 
7. Visualize your training process and results
8. Apply feature scaling (pen & paper)

## Requirements
### Knowledge

You should have a basic knowledge of:
- Univariate linear regression
- Multivariate linear regression
- Squared error
- Gradient descent
- numpy
- matplotlib

Suitable sources for acquiring this knowledge are:
- [Multivariate Linear Regression Notebook](http://christianherta.de/lehre/dataScience/machineLearning/basics/multivariate_linear_regression.php) by Christian Herta and his [lecture slides](http://christianherta.de/lehre/dataScience/machineLearning/multivariateLinearRegression.pdf) (German)
- Chapter 2 of the open classroom [Machine Learning](http://openclassroom.stanford.edu/MainFolder/CoursePage.php?course=MachineLearning) by Andrew Ng
- Chapter 5.1 of [Deep Learning](http://www.deeplearningbook.org/contents/ml.html) by Ian Goodfellow 
- Some parts of chapter 1 and 3 of [Pattern Recognition and Machine Learning](https://www.microsoft.com/en-us/research/people/cmbishop/#!prml-book) by Christopher M. Bishop
- [numpy quickstart](https://docs.scipy.org/doc/numpy-1.15.1/user/quickstart.html)
- [Matplotlib tutorials](https://matplotlib.org/tutorials/index.html)

### Python Modules

By [deep.TEACHING](https://www.deep-teaching.org/) convention, all python modules needed to run the notebook are loaded centrally at the beginning. 



```python
# External Modules
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d.axes3d import Axes3D
from matplotlib import cm

%matplotlib notebook
```

## Exercise - Multivariate Linear Regression

We will only use two features in this notebook, so we are still able to plot them together with the target in a 3D plot. But your implementation should also be capable of handling more (except the plots). 

### Create Features

First we will create some features. The features should be in a 2D numpy array, the rows separating the different feature vectors (training examples), the columns containing the features. Each feature should be **uniformly** distributed in a specifiable range.

**Task:**

Implement the function to generate a feature matrix (numpy array).


```python
def create_feature_matrix(sample_size, n_features, x_min, x_max):
    '''creates random feature vectors based on a lienar function in a given interval
    
    Args:
        sample_size: number feature vectors
        n_features: number of features for each vector
        x_min: lower bound value ranges
        x_max: upper bound value ranges
    
    Returns:
        x: 2D array containing feature vecotrs with shape (sample_size, n_features)
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def create_feature_matrix(sample_size, n_features, x_min, x_max):
    '''creates random feature vectors based on a lienar function in a given interval
    
    Args:
        sample_size: number feature vectors
        n_features: number of features for each vector
        x_min: lower bound value ranges
        x_max: upper bound value ranges
    
    Returns:
        x: 2D array containing feature vecotrs with shape (sample_size, n_features)
    '''
    
    return np.random.uniform(x_min, x_max, (sample_size, n_features))
```


```python
sample_size = 100
n_features = 2
x_min = [1.5, -0.5]
x_max = [11., 5.0]

X = create_feature_matrix(sample_size, n_features, x_min, x_max)
X
```




    array([[ 2.55343578,  4.09882663],
           [ 9.14676162,  1.75409426],
           [ 6.87061128,  3.088036  ],
           [ 3.31754536, -0.10581374],
           [10.04676143,  4.91327452],
           [ 9.53248102, -0.21669548],
           [ 4.74762432,  2.11592889],
           [ 3.04610114,  0.66861665],
           [ 8.99810878,  4.10601558],
           [ 2.03046397,  1.35255697],
           [ 4.6004498 ,  1.01205624],
           [ 5.2053511 ,  4.85018639],
           [ 8.92142801,  3.40039974],
           [ 8.41672928,  3.66919872],
           [ 5.01825985,  2.51420676],
           [ 9.08429314,  1.04931984],
           [ 2.15416736,  0.5268242 ],
           [ 6.2284124 ,  4.54399381],
           [ 8.10072993,  1.84264001],
           [ 3.71944741,  0.94380792],
           [ 4.35883408,  1.26337124],
           [ 4.43690719,  0.99997779],
           [ 2.51908645,  0.34135765],
           [ 4.64194197,  1.37604907],
           [ 5.56066631,  4.24084101],
           [ 3.31660389,  3.4357965 ],
           [ 7.46692545, -0.13872683],
           [ 8.67304278,  1.53969574],
           [ 2.25378335,  2.1360613 ],
           [ 2.6540387 ,  1.03259774],
           [ 5.02662661,  2.87918316],
           [ 4.12726842,  4.23028426],
           [ 5.68731577,  2.70034396],
           [10.74272894, -0.16190299],
           [ 5.57938986,  3.30598652],
           [ 9.26774676,  3.15054548],
           [ 4.23662404,  4.93963246],
           [ 3.84702751,  0.70003255],
           [ 7.23072073,  2.14386805],
           [ 1.57051795,  3.70768626],
           [ 5.40511937,  3.34361318],
           [ 1.86679871,  4.50691358],
           [ 2.64773238,  3.48214973],
           [ 5.29444443,  1.06644367],
           [ 8.15360051,  1.91813851],
           [ 1.58028104,  3.52942058],
           [ 5.57410943,  0.36621873],
           [ 3.51710302,  4.05378562],
           [ 1.9668387 ,  3.08308516],
           [ 5.59368553,  0.48926087],
           [ 4.72382896,  1.27217418],
           [ 4.75903772,  4.02492114],
           [10.32314113,  3.45291417],
           [ 8.3607832 ,  1.41331422],
           [ 5.02277249, -0.17925784],
           [ 4.96284139,  3.16654885],
           [ 5.25042165,  4.12198925],
           [ 3.87229428,  2.14576293],
           [ 9.53374907, -0.05896827],
           [ 8.15028007,  3.9327573 ],
           [ 9.56958557,  2.32766841],
           [ 7.64470758,  1.88364085],
           [ 9.19808575,  2.63258445],
           [ 8.54669469,  2.58779554],
           [ 2.0067687 ,  3.09877182],
           [ 3.62733662,  0.25622579],
           [10.05874243,  0.19438902],
           [ 8.14099486,  4.01139473],
           [ 7.21657038,  2.4667673 ],
           [ 8.1279165 ,  0.24171846],
           [ 1.85488792,  2.11330979],
           [ 6.03449596,  3.593834  ],
           [ 3.28212156, -0.2277349 ],
           [ 7.08842704,  2.76401252],
           [ 8.22070264,  3.87316012],
           [ 7.86751765,  2.55121306],
           [ 3.23420664,  4.64625445],
           [ 2.26114604,  3.59172898],
           [ 9.69223623,  1.5591381 ],
           [ 8.4105376 ,  4.45931284],
           [ 2.22748906,  2.79878555],
           [ 7.32528596,  2.07677038],
           [10.36156224,  3.67408725],
           [ 6.55661923, -0.46412987],
           [ 3.9209253 , -0.47884172],
           [ 8.8067252 ,  4.58892015],
           [ 4.13219722, -0.30789565],
           [10.96105526,  2.83739255],
           [ 6.83442425,  2.91034924],
           [ 9.49990468, -0.17566033],
           [ 6.95625609,  1.48659351],
           [ 4.37699785,  1.08344818],
           [ 2.00884202,  0.23834419],
           [ 3.33201506,  2.12857654],
           [ 1.50306421,  0.71705611],
           [ 9.46405326,  1.78328114],
           [ 9.61235747,  1.24012355],
           [ 1.96521893,  0.76458039],
           [ 6.35091656,  4.27638952],
           [ 6.44543344,  2.78178998]])




```python
assert len(X[:,0]) == sample_size
assert len(X[0,:]) == n_features
for i in range(n_features):
    assert np.max(X[:,i]) <= x_max[i]
    assert np.min(X[:,i]) >= x_min[i]
```

### Linear Hypothesis


A short recap, a hypothesis $h_\theta({\bf x})$ is a certain function that we believe is similar to a target function that we like to model. A hypothesis $h_\theta({\bf x})$ is a function of ${\bf x}$ with fixed parameters $\Theta$. 

Here we have $n$ features ${\bf x} = \{x_1, \ldots, x_n \}$ and $n+1$ parameters $\Theta = \{\theta_0, \theta_1 \ldots, \theta_n \}$:

$$
h_\theta({\bf x}) = \theta_{0} + \theta_{1} x_1 + \ldots \theta_n x_n 
$$

adding an extra element to $\vec x$ for convenience, this could also be rewritten as:

$$
h_\theta({\bf x}) = \theta_{0} x_0 + \theta_{1} x_1 + \ldots \theta_n x_n 
$$

with $x_0 = 1$ for all feature vectors (training examples).

Or treating ${\bf x}$ and $\Theta$ as vectors:

$$
h(\vec x) = \vec x'^T \vec \theta
$$

with:

$$
\vec x = \begin{pmatrix} 
x_1 & x_2 & \ldots & x_n \\
\end{pmatrix}^T
\text{   and   }
\vec x' = \begin{pmatrix} 
1 & x_1 & x_2 & \ldots & x_n \\
\end{pmatrix}^T
$$

and

$$
\vec \theta = \begin{pmatrix} 
\theta_0 & \theta_1 & \ldots & \theta_n \\
\end{pmatrix}^T
$$

Or for the whole data set at once: The rows in $X$ separate the different feature vectors, the columns contain the features. 

$$
\vec h_\Theta(X) = X' \cdot \vec \theta
$$

the vector $\vec h(X) = \left( h(\vec x^{(1)}),h(\vec x^{(2)}), \dots, h(\vec x^{(m)}) )\right)^T$ contains all predictions for the data batch $X$.

with:

$$
\begin{align}
X &= \begin{pmatrix} 
x_1^{(1)} & \ldots & x_n^{(1)} \\
x_1^{(2)} & \ldots & x_n^{(2)} \\
\vdots &\vdots &\vdots \\
x_1^{(m)} & \ldots & x_n^{(m)} \\
\end{pmatrix}
&=
\begin{pmatrix} 
\vec x^{(1)T} \\
\vec x^{(2)T}  \\
\vdots  \\
\vec x^{(m)T}  \\
\end{pmatrix}
\end{align}
$$
respectively
$$
\begin{align}
X' = \begin{pmatrix} 
1 & x_1^{(1)} & \ldots & x_n^{(1)} \\
1 & x_1^{(2)} & \ldots & x_n^{(2)} \\
\vdots &\vdots &\vdots &\vdots \\
1 & x_1^{(m)} & \ldots & x_n^{(m)} \\
\end{pmatrix}
&=
\begin{pmatrix} 
\vec x'^{(1)T} \\
\vec x'^{(2)T}  \\
\vdots  \\
\vec x'^{(m)T}  \\
\end{pmatrix}
\end{align}
$$

**Task:**

Implement hypothesis $\vec h_\Theta(X)$ in the method `linear_hypothesis` and return it as a function. Implement it the computationally efficient (**pythonic**) way by not using any loops and handling all data at once (use $X$ respectively $X'$).

**Hint:**

Of course you are free to implement as many helper functions as you like, e.g. for transforming $X$ to $X'$, though you do not have to. Up to you.


```python
def linear_hypothesis(thetas):
    ''' Combines given list argument in a linear equation and returns it as a function
    
    Args:
        thetas: list of coefficients
        
    Returns:
        lambda that models a linear function based on thetas and x
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def linear_hypothesis(thetas):
    ''' Combines given list argument in a linear equation and returns it as a function
    
    Args:
        thetas: list of coefficients
        
    Returns:
        lambda that models a linear function based on thetas and x
    '''
    return lambda x: np.concatenate((np.ones((len(x[:,0]), 1)), x), axis=1).dot(thetas)
```


```python
assert len(linear_hypothesis([.1,.2,.3])(X)) == sample_size
```

### Generate Target Values

**Task:**

Use your implemented `linear_hypothesis` inside the next function to generate some target values $Y$. Additionally add some Gaussian noise.


```python
def generate_targets(X, theta, sigma):
    ''' Combines given arguments in a linear equation with X, 
    adds some Gaussian noise and returns the result
    
    Args:
        X: 2D numpy feature matrix
        theta: list of coefficients
        sigma: standard deviation of the gaussian noise
        
    Returns:
        target values for X
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def generate_targets(X, theta, sigma):
    ''' Combines given arguments in a linear equation with X, 
    adds some Gaussian noise and returns the result
    
    Args:
        X: 2D numpy feature matrix
        theta: list of coefficients
        sigma: standard deviation of the gaussian noise
        
    Returns:
        target values for X
    '''
    noise = np.random.randn(len(X)) * sigma
    h = linear_hypothesis(theta)
    y = h(X) + noise
    return y
```


```python
theta = (2., 3., -4.)
sigma = 3.
y = generate_targets(X, theta, sigma)
```


```python
assert len(y) == sample_size
```

### Plot The Data

**Task:**

Plot the data $\mathcal D = \{((x^{(1)}_1,x^{(1)}_2)^T,y^{(1)}), \ldots, ((x^{(n)}_1,x^{(n)}_2)^T,y^{(n)})\}$ in a 3D scatter plot. The plot should look like the following:



**Sidenote:**

The command `%matplotlib notebook` (instead of `%matplotlib inline`) creates an interactive (e.g. rotatable) plot.


```python
%matplotlib notebook

def plot_data_scatter(features, targets):
    """ Plots the features and the targets in a 3D scatter plot
    
    Args:
        features: 2D numpy-array features
        targets: ltargets
    """
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

%matplotlib notebook

def plot_data_scatter(features, targets):
    """ Plots the features and the targets in a 3D scatter plot
    
    Args:
        features: 2D numpy-array features
        targets: ltargets
    """
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.scatter(features[:,0], features[:,1], targets, c='r')
    ax.set_xlabel('$x_1$')
    ax.set_ylabel('$x_2$')
    ax.set_zlabel('$y$')
```


```python
plot_data_scatter(X, y)
```


    <IPython.core.display.Javascript object>






### Cost Function
A cost function $J$ depends on the given training data $D$ and hypothesis $h_\theta(\vec x)$. In the context of the linear regression, the cost function measures how "wrong" a model is regarding its ability to estimate the relationship between $\vec x$ and $y$ for specific $\Theta$ values. Later we will treat this as an optimization problem and try to minimize the cost function $J_{\mathcal D}(\Theta)$ to find optimal $\theta$ values for our hypothesis $h_\theta(\vec x)$. The cost function we use in this exercise is the [Mean-Squared-Error](https://en.wikipedia.org/wiki/Mean_squared_error) cost function:

\begin{equation}
    J_{\mathcal D}(\Theta)=\frac{1}{2m}\sum_{i=1}^{m}{(h_\Theta(\vec x^{(i)})-y^{(i)})^2}
\end{equation}

Implement the cost function $J_D(\Theta)$ in the method `mse_cost_function`. The method should return a function that takes the values of $\Theta$ as an argument.

As a sidenote, the terms "loss function" or "error function" are often used interchangeably in the field of Machine Learning.


```python
def mse_cost_function(x, y):
    ''' Implements MSE cost function as a function J(theta) on given traning data 
    
    Args:
        x: vector of x values 
        y: vector of ground truth values y 
        
    Returns:
        lambda J(theta) that models the cost function
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def mse_cost_function(x, y):
    ''' Implements MSE cost function as a function J(theta) on given traning data 
    
    Args:
        x: vector of x values 
        y: vector of ground truth values y 
        
    Returns:
        lambda J(theta) that models the cost function
    '''
    assert(len(x) == len(y))
    m = len(x)
    return lambda theta: 1./(2. * float(m)) * np.sum((linear_hypothesis(theta)(x) - y )**2)
```

Review the cell in which you generate the target values and note the theta values, which were used for it (If you haven't edited the default values, it should be `[2, 3, -4]`)

**Optional:**

Try a few different values for theta to pass to the cost function - Which thetas result in a low error and which produce a great error?


```python
J = mse_cost_function(X, y)
print(J(theta))
```

    4.675479225972847


###  Gradient Descent

A short recap, the gradient descent algorithm is a first-order iterative optimization for finding a minimum of a function. From the current position in a (cost) function, the algorithm steps proportional to the negative of the gradient and repeats this until it reaches a local or global minimum and determines. Stepping proportional means that it does not go entirely in the direction of the negative gradient, but scaled by a fixed value $\alpha$ also called the learning rate. Implementing the following formalized update rule is the core of the optimization process:

\begin{equation}
    \theta_{j}^{new} \leftarrow \theta_{j}^{old} - \alpha * \frac{\partial}{\partial\theta_{j}} J(\vec \theta^{old})
\end{equation}

**Task:**

Implement the function to update all theta values.


```python
def update_theta(x, y, theta, learning_rate):
    ''' Updates learnable parameters theta 
    
    The update is done by calculating the partial derivities of 
    the cost function including the linear hypothesis. The 
    gradients scaled by a scalar are subtracted from the given 
    theta values.
    
    Args:
        x: 2D numpy array of x values
        y: array of y values corresponding to x
        theta: current theta values
        learning_rate: value to scale the negative gradient 
        
    Returns:
        theta: Updated theta vector
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def update_theta(x, y, theta, learning_rate):
    ''' Updates learnable parameters theta 
    
    The update is done by calculating the partial derivities of 
    the cost function including the linear hypothesis. The 
    gradients scaled by a scalar are subtracted from the given 
    theta values.
    
    Args:
        x: 2D numpy array of x values
        y: array of y values corresponding to x
        theta: current theta values
        learning_rate: value to scale the negative gradient 
        
    Returns:
        theta: Updated theta vector
    '''
    m = len(x[:,0])
    # update rules
    x_ = np.concatenate((np.ones([m,1]), x), axis=1)
    theta = theta  - alpha / float(m) * (x_.T.dot(linear_hypothesis(theta)(x) - y) )

    
    
    return theta
```

Using the `update_theta` method, you can now implement the gradient descent algorithm. Iterate over the update rule to find the values for $\vec \theta$ that minimize our cost function $J_D(\vec \theta)$. This process is often called training of a machine learning model. 

**Task:**
- Implement the function for the gradient descent.
- Create a history of all theta and cost values and return them.


```python
def gradient_descent(learning_rate, theta, iterations, x, y):
    ''' Minimize theta values of a linear model based on MSE cost function
    
    Args:
        learning_rate: scalar, scales the negative gradient 
        theta: initial theta values
        x: vector, x values from the data set
        y: vector, y values from the data set
        iterations: scalar, number of theta updates
        
    Returns:
        history_cost: cost after each iteration
        history_theta: Updated theta values after each iteration
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def gradient_descent(learning_rate, theta, iterations, x, y):
    ''' Minimize theta values of a linear model based on MSE cost function
    
    Args:
        learning_rate: scalar, scales the negative gradient 
        theta: initial theta values
        x: vector, x values from the data set
        y: vector, y values from the data set
        iterations: scalar, number of theta updates
        
    Returns:
        history_cost: cost after each iteration
        history_theta: Updated theta values after each iteration
    '''
    history_cost = np.zeros(nb_iterations)
    history_theta = np.zeros([iterations, len(theta)])
    cost = mse_cost_function(x, y);
    for i in range(nb_iterations):
        history_theta[i] = theta
        history_cost[i] = cost(theta)
        theta = update_theta(x, y, theta, learning_rate)
    return history_cost, history_theta
```

### Training and Evaluation

**Task:**

Choose an appropriate learning rate, number of iterations and initial theta values and start the training


```python
# Your implementation:

alpha = 42.42 # assign an appropriate value
nb_iterations = 1337 # assign an appropriate value
start_values_theta = [42., 42., 42.] # assign appropriate values
history_cost, history_theta = gradient_descent(alpha, start_values_theta, nb_iterations, X, y)
```

    <ipython-input-17-bfcd484b6e8e>:15: RuntimeWarning: overflow encountered in square
      return lambda theta: 1./(2. * float(m)) * np.sum((linear_hypothesis(theta)(x) - y )**2)



```python
# Solution

alpha = 0.01
nb_iterations = 100
start_values_theta = [1.,-1.,.5]
history_cost, history_theta = gradient_descent(alpha, start_values_theta, nb_iterations, X, y)
```

Now that the training has finished we can visualize our results.

**Task:**

Plot the costs over the iterations. If you have used `fig = plt.figure()` and `ax = fig.add_subplot(111)` in the last plot, use it again here, else the plot will be added to the last plot instead of a new one.

Your plot should look similar to this one:




```python
def plot_progress(costs):
    """ Plots the costs over the iterations
    
    Args:
        costs: history of costs
    """
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def plot_progress(costs):
    """ Plots the costs over the iterations
    
    Args:
        costs: history of costs
    """
    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.plot(np.array(range(len(costs))), costs)
    ax.set_xlabel('Iterationen')
    ax.set_ylabel('Kosten')
    ax.set_title('Fortschritt')
```


```python
plot_progress(history_cost)
print("costs before the training:\t ", history_cost[0])
print("costs after the training:\t ", history_cost[-1])
```


    <IPython.core.display.Javascript object>






    costs before the training:	  190.45165080692664
    costs after the training:	  4.947880491325002


**Task:**

Finally plot the decision hyperplane (just a plain plane here though) together with the data in a 3D plot.

Your plot should look similar to this one:




```python
def evaluation_plt(x, y, final_theta):
    ''' Plots the data x, y together with the final model
    
    Args:
        cost_hist: vector, history of all cost values from a opitmization
        theta_0: scalar, model parameter for boundary
        theta_1: scalar, model parameter for boundary
        x: vector, x values from the data set
        y: vector, y values from the data set
    '''
    raise NotImplementedError("You should implement this!")
```


```python
# Solution

def evaluation_plt(x, y, final_theta):
    g = 100
    x_1 = np.linspace(np.min(X[:,0]-1), np.max(X[:,0]+1), g)
    x_2 = np.linspace(np.min(X[:,1]-1), np.max(X[:,1]+1), g)
    X1, X2 = np.meshgrid(x_1, x_2)
    
    Y = final_theta[0] + final_theta[1] * X1 + final_theta[2] * X2
    print(Y.shape)

    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')

    ax.plot_surface(X1, X2, Y, cmap=cm.jet, rstride=5, cstride=5, antialiased=True, shade=True, alpha=0.5, linewidth=0. )
    ax.scatter(x[:,0], x[:,1], y, c='r')
    ax.set_xlabel('$x_1$')
    ax.set_ylabel('$x_2$')
    ax.set_zlabel('$y$')
```


```python

evaluation_plt(X, y, history_theta[-1])
print("thetas before the training:\t", history_theta[0])
print("thetas after the training:\t", history_theta[-1])
```

    (100, 100)



    <IPython.core.display.Javascript object>






    thetas before the training:	 [ 1.  -1.   0.5]
    thetas after the training:	 [ 1.24635143  3.02581154 -3.83837415]


### Feature Scaling

Now suppose the following features $X$:


```python
X = np.array([[0.0001, 2000],
       [0.0002, 1800],
       [0.0003, 1600]], dtype=np.float32)

sample_size = len(X[:,0])
print(X)
```

    [[1.0e-04 2.0e+03]
     [2.0e-04 1.8e+03]
     [3.0e-04 1.6e+03]]


**Task:**

This task can be done via **pen & paper** or by inserting some code below. Either way, you should be able to solve both tasks below on paper only using a calculator.

1. Apply feature scaling onto `X` using the *mean* and the *standard deviation*. What values do the scaled features have?
    * *Optional:*

        You can even execute the cell above and start running your notebook again from top (all **except** executing the cell to generate your features, which would overwrite these new features).

        When you start training you should notice that your costs do not decrease, maybe even increase, if you have not adjusted your learning rate (training might also throw an overflow warning).

**Task:**

2. After the training with scaled features your new $\\vec \theta'$ values will something like: $\vec \theta'=\left(-7197,  326, -326\right)^T$ (you can try training with but you do not have to). 

    Suppose $\vec \theta'=\left(-7197,  326, -326\right)^T$. What are the corresponding $\theta_j$ values for the unscaled data?

    * (If you did train your model with the scaled features, the resulting $\theta_j$ should really be $\vec \theta'=\left(-7197,  326, -326\right)^T$


```python
# Solution for 1) + 2)

### Generate new targets for new features X

theta = (2., 3., -4.)
sigma = 3.
y = generate_targets(X, theta, sigma)
```


```python
mu = X.mean(axis=0)
std = X.std(axis=0)
print(mu, std)
```

    [2.0000001e-04 1.8000000e+03] [8.1649661e-05 1.6329932e+02]



```python
X_scaled = (X - mu) / std
X_scaled
```




    array([[-1.2247449e+00,  1.2247449e+00],
           [-1.7822383e-07,  0.0000000e+00],
           [ 1.2247449e+00, -1.2247449e+00]], dtype=float32)




```python
### Calculate mean and standarddeviation and scale X

mu1 = X[:,0].mean()
std1 = np.sqrt(np.var(X[:,0]))
mu2 = X[:,1].mean()
std2 = np.sqrt(np.var(X[:,1]))

print("meanX1:\t", mu1, "\tstdX1:\t", std1)
print("meanX2:\t", mu2, "\tstdX2:\t", std2)

X[:,0] = (X[:,0] - mu1) / std1
X[:,1] = (X[:,1] - mu2) / std2
print("X scaled:\n", X)

```

    meanX1:	 0.00020000001 	stdX1:	 8.164966e-05
    meanX2:	 1800.0 	stdX2:	 163.29932
    X scaled:
     [[-1.2247449e+00  1.2247449e+00]
     [-1.7822383e-07  0.0000000e+00]
     [ 1.2247449e+00 -1.2247449e+00]]



```python
### Train with scaled X

alpha = 0.5
nb_iterations = 100

theta = (2., 3., -4.)
sigma = 3.
y = generate_targets(X, theta, sigma)

start_values_theta = [1.,-1.,.5]
history_cost, history_theta = gradient_descent(alpha, start_values_theta, nb_iterations, X_scaled, y)

print("thetas before the training:\t", history_theta[0])
print("thetas after the training:\t", history_theta[-1])
print("----------------------------")
print("costs before the training:\t", history_cost[0])
print("costs after the training:\t", history_cost[-1])
```

    thetas before the training:	 [ 1.  -1.   0.5]
    thetas after the training:	 [ 1.07576387  2.82165621 -3.32163366]
    ----------------------------
    costs before the training:	 43.92675099402078
    costs after the training:	 14.713938810789191



```python
### Compute theta for unscaled data

theta_hat = np.array([-7197.,  326., -326.])
print("theta_hat:\t", theta_hat)

theta_hat_rescaled = np.array([0.,0.,0.])
theta_hat_rescaled[0] = theta_hat[0] - theta_hat[1]*mu1/std1 - theta_hat[2]*mu2/std2
theta_hat_rescaled[1] = theta_hat[1] / std1
theta_hat_rescaled[2] = theta_hat[2] / std2
print("theta_rescaled:\t", theta_hat_rescaled)

start_values_theta = theta_hat_rescaled
```

    theta_hat:	 [-7197.   326.  -326.]
    theta_rescaled:	 [-4.40213221e+03  3.99266812e+06 -1.99633414e+00]



```python
### Rescale X to original

X = np.array([[0.0001, 2000],
       [0.0002, 1800],
       [0.0003, 1600]], dtype=np.float32)
```


```python
### Test costs for original X and new rescaled theta

alpha = 0.000
nb_iterations = 10
history_cost, history_theta = gradient_descent(alpha, start_values_theta, nb_iterations, X, y)

print("costs with unscaled X and new theta:\t ", history_cost[-1])
```

    costs with unscaled X and new theta:	  26114727.727051158


**Solution 3:**

With polynomila regression features like `1000` would easily skyrocket even more. And worse: `0.00001` would become too small very fast and lose precision. In the worst case values would overflow or become zero.

## Summary and Outlook

During this exercise, the linear regression was extended to multidimensional feature space and feature scaling was practiced. You should be able to answer the following questions:
- How does the implementation of the multivariate regression differ from the univariate one?
- Why do we apply feature scaling?
- Why does feature scaling help?

## Licenses

### Notebook License (CC-BY-SA 4.0)

*The following license applies to the complete notebook, including code cells. It does however not apply to any referenced external media (e.g., images).*

Exercise: Multivariate Linear Regression <br/>
by Christian Herta, Klaus Strohmenger<br/>
is licensed under a [Creative Commons Attribution-ShareAlike 4.0 International License](http://creativecommons.org/licenses/by-sa/4.0/).<br/>
Based on a work at https://gitlab.com/deep.TEACHING.


### Code License (MIT)

*The following license only applies to code cells of the notebook.*

Copyright 2018 Christian Herta, Klaus Strohmenger

Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

0280976f4ec9e9260f474ffff9decb7aedbc6d74
ipynb
Jupyter Notebook
_notebooks/Multivariate_lineare_regression.ipynb
yachty66/openblog
9935b588044bcb9c84853aabcc042bba089916f1
Apache-2.0
_notebooks/Multivariate_lineare_regression.ipynb
yachty66/openblog
9935b588044bcb9c84853aabcc042bba089916f1
Apache-2.0
_notebooks/Multivariate_lineare_regression.ipynb
yachty66/openblog
9935b588044bcb9c84853aabcc042bba089916f1
Apache-2.0
Qwen/Qwen-72B
1. YES
2. YES