

*Images generated by a deep neural network that interprets text to generate images (DALL·E 2)*

# PC lab: intro to Neural networks & PyTorch

Deep learning is the subfield of machine learning that concerns neural networks with representation learning capabilities. As of recent years, it is arguably the most quickly growing field within machine learning, enjoying major breakthroughs every year (Listing a couple ones from last year(s): GPT-3, AlphaFold v2, DALL·E 2, AlphaZero). Although the popularity of neural nets is a recent phenomenon, they were first described by Warren McCulloch and Walter Pitts in 1943. Early progress in training competitive neural networks was stalled by a multitude of reasons, such as the limited computer resources, sub-optimal network architectures and the use of smaller datasets. In this PC-lab we will introduce you to the basics of implementing a neural network using contemporary practices.

# 1.1 Background

### From linear models to neurons to neural networks

The core unit of every (artificial) neural network is considered the neuron. Every neuron can be observed as a linear combination of **one or more inputs** $\mathbf{x}$ with weights $\mathbf{w}$ (and optionally adding a **bias** $b$), outputting a **single output** $z$:

$z = \sum\limits_{i=1}^{n}(w_ix_i) + b$

equivalently written as dot product:

$z = \mathbf{x} \cdot \mathbf{w} + b$

When multiple neurons are placed next to eachother, we get multiple output neurons, this is called a **layer**. Mathematically, our weights and biases (intercepts) now become a matrix and vector respectively, and we obtain as output a vector $\mathbf{z}$:

$\mathbf{z} = \mathbf{x} \cdot \mathbf{W} + \mathbf{b}$

Multiple **layers** can be stacked sequentially to make up a deep neural network. Performing this gives us a big linear model, because multiple matrix multiplications performed consecutively can also be written as just one: (i.e. $\mathbf{x} \cdot \mathbf{W_1} \cdot \mathbf{W_2} ...$ could also be written as just $\mathbf{x} \cdot \mathbf{W_{all}}$). To obtain the *deep learning magic*, we need to make the whole thing **non-linear**. We do this by adding **activation functions** after every (hidden) layer. The most classical activation is the sigmoid activation $\sigma()$, used also in logistic regression. Nowadays, we usually opt for a more simple activation function: the **ReLU** $ReLU(z) = max\{0,z\}$. This function has the favorable property that its derivative is very efficient to compute (1 when $z$ is positive, 0 when it is negative). It also acts as a switch: a neuron will have a "dead" ($0$) activation whenever $z$ is negative.

For the output layer of a neural network: our activation depends on the task at hand. For binary classification, we use sigmoid to constrain our output between 0 and 1. For multi-class, we use a softmax operation so the output of all neurons sums to 1. For regression, we simply do not use an activation (or a custom one depending on your data: if you already know that your outputs can't be negative but can take all positive numbers ($\mathbb{R}^+$), then maybe a ReLU activation in the output nodes makes sense)

In order to build more intuition for neural networks: consider the following figure where we "visually build up" a neural network starting from Linear regression with four input features (**a**), to Logistic regression (**b**), to an archetypical output neuron with ReLU activation (**c**). For multi-output settings, we visualize multi-output regression (**d**), multi-label classification (more than one class can be 1 in a sample) (**e**), and multi-class classification via softmax (**f**). Finally, a simple neural network with two hidden layers for binary classification (sigmoid output head) is shown under **g**.



**This figure makes it crystal clear that the most simple neural network is just a bunch of linear regressions stacked on top of eachother with non-linearities inbetween.** More advanced neural network architectures exist that modify how we make information flow between inputs. In this example, everything is just connected to everything with linear weights. This type of neural network is what we call an **MLP** or a **multi layer perceptron**.

Keep in mind that all of the above methods usually fit a bias/intercept in addition to weights fitted on the input features.
For an MLP, visually it would look like this:



<div class="alert alert-success">

<b>THOUGHT EXERCISE:</b> 
<p> How much weights does the model pictured above have (including biases)? </p>

</div>

Keep in mind that every layer consists of (1) a matrix multiplication of an input $X \in \mathbb{R}^{N, D}$ with a set of weights $W \in \mathbb{R}^{D, M}$ and (2) a non-linearity. With $D$ the number of input features/dimensions/nodes (including one node/feature for the intercept, see the following equation) and $M$ the number of output dimensions/nodes.
\begin{equation}
XW =
\begin{bmatrix}
1 & x_{0,1} & ...  & x_{0,D-1} & x_{0,D} \\
1 & x_{1,1} & ... & x_{1,D-1} & x_{1,D} \\
... & ... & ... & ... & ...\\
1 & x_{N-1,1} & ... & x_{N-1,D-1} & x_{N-1,D} \\
1 & x_{N,1} & ...  & x_{N,D-1} & x_{N,D} \\
\end{bmatrix}
\begin{bmatrix}
W_{0,0} & W_{0,1} & ...  & W_{0,M-1} & W_{0,M} \\
W_{1,0} & W_{1,1} & ... & W_{1,M-1} & W_{1,M} \\
... & ... & ... & ... & ...\\
W_{D-1,0} & W_{D-1,1} & ... & W_{D-1,M-1} & W_{D-1,M} \\
W_{D,0} & W_{D,1} & ...  & W_{D,M-1} & W_{D,M} \\
\end{bmatrix}
\end{equation}

The whole network is trained by first performing a forward pass to get predictions and compute a loss w.r.t. ground truth known labels, then doing backpropagation (essentially applying the chain rule of derivatives) to obtain the gradient of all weights w.r.t. the loss. These gradients are then used by gradient descent (or more modern variants such as Adam) to optimize the neural network.

Training neural networks is (typically) more computationally demanding than more traditional machine learning methods, and we usually use neural networks when we have large datasets. For these two reasons, it is not a good idea to process the whole dataset at once in every training step/iteration. Therefore, to train the network, we process samples in batches, called **(mini-batch) stochastic gradient descent**. This allows for faster training and improved convergence of the loss during gradient descent. Advantages of stochastic gradient descent or other optimization algorithms for loss calculation are not discussed in this PC-lab, but have been [extensively discussed](https://ruder.io/optimizing-gradient-descent/) before.

### Dropout and Normalization layers

Dropout is a popular addition to neural networks. It is a form of regularization during which we stochastically deactivate a percentage of neurons in every training step by putting their activation to zero. This regularization only happens during training, as during testing we (usually) want deterministic outputs. Conceptually, it is similar to other regularization techniques such as ridge regression and subsampling of features in random forest training, in the sense that it will force our model to look at all features, because sometimes one single feature will not be available during training. The difference here is that we do it by stochastically putting nodes to zero during training, and that we can perform it not only on our input features, but also on our hidden nodes.

Mathematically, it can be performed by simply sampling a boolean vector and doing element-wise multiplication.

Visually, it would look a bit like this, where nodes in cyan are dropped out, and the associated cyan weights do not have any influence on training anymore (in that training step):



Normalization is another operation that we usually put between layers in order to more easily train our networks. In essence, this ensure that our hidden outputs are numerically well behaved, which speeds up training and helps with better convergence. The most popular ones are batch normalization and layer normalization.
More info here: [link](https://www.baeldung.com/cs/normalizing-inputs-artificial-neural-network)

With current techniques, most people often put the order of operations like this:
`Layer -> Activation -> Dropout -> Normalization -> and_repeat ...`

## 1.2 PyTorch

To implement neural networks with more ease, a few high-level python libraries are available: (PyTorch, TensorFlow/keras, JAX, ...). These libraries provide functionality in terms of automatic differentiation (backprop), ready-to-use implementations for various layers, loss functions ...

In this lab, we will use [PyTorch](https://pytorch.org). PyTorch is the most popular library for deep learning in academia as of today. For this course it offers the advantage that it has the most 'pythonic' syntax, to the point where almost all NumPy functions have a PyTorch counterpart.

If you want to run this notebook locally, you can find the installation instructions for PyTorch [here](https://pytorch.org/get-started/locally/). Make sure to select the right installation options depending on your system (if you have a GPU or not).



```python
import torch
import numpy as np
```

### Tensors

Tensors are the fundamental data structures in PyTorch. They are analogous to NumPy arrays. The difference is that tensors can also run on GPU hardware. GPU hardware is optimized for many small computations. Matrix multiplications, the building blocks of all deep learning, run orders-of-magnitude faster on GPU than on CPU. Let's see how tensors are constructed and what we can do with them:


```python
x = [[5,8],[9,8]]
print(torch.tensor(x))
print(np.array(x))
```


```python
x_numpy = np.array(x)
print(torch.from_numpy(x_numpy))

x_torch = torch.tensor(x)
print(x_torch.numpy())
```


```python
print(np.random.randn(8).shape)
print(np.random.randn(8,50).shape)


print(torch.randn(8).shape) # an alternative for .shape in PyTorch is .size()
print(torch.randn(8,50).shape)
```


```python
print(np.zeros((8,50)).shape)
print(torch.zeros(8,50).shape) # works with 'ones' as well
```

In PyTorch, the standard data type for floats is `float32`, which is synonymous to `float` within its framework. `float64` is synonymous to `double`.
This is different from the NumPy defaults and naming conventions: NumPy default data type for float is `float64`. Keep this in mind when converting numpy arrays to tensors and back!


```python
print(np.zeros(8).dtype)
print(torch.zeros(8).dtype)
```

Conversion of data types:


```python
x = torch.randn(8)
print(x.dtype)
x = x.to(torch.float64)
print(x.dtype)
```

`torch.long` is synonymous to `torch.int64`. The only difference between int32 and int64 is the amount of bytes with which you will store every integer. If you go up to very high numbers, you will get numerical overflow faster with more compressed data types. We recommend you to always use the defaults: `torch.long` and `torch.float`








```python
x = torch.randint(low=0, high=8, size=(8,), dtype=torch.int32)
print(x)
print(x.dtype)
x = x.to(torch.long)
print(x.dtype)
```

Indexing and other operations work as in NumPy arrays


```python
x = torch.randn(8,50,60)
print(x.shape)
print(x[:4,10:-10].shape)
x[0,0,:10] = 0
print(x[0,0,:16])

print(torch.min(x), torch.max(x), torch.min(torch.abs(x)))
# most of these functions are also tensor methods:
print(x.min(), x.max(), x.abs().min())
```

Joining tensors via concatenation:


```python
print(x.shape)
x_cat0 = torch.cat([x, x], dim=0)
print(x_cat0.shape)
x_cat1 = torch.cat([x, x, x], dim=1)
print(x_cat1.shape)
```

Matrix multiplication: let's say we have an input `x`, consisting of 8 samples with 26 features, that we linearly combine with weights `w` to get a single output for every sample:




```python
x = torch.randn(8,26)
w = torch.randn(26,1)

y_hat = torch.matmul(x, w) # an alternative and equivalent syntax is x @ w
print(y_hat)
print(y_hat.shape)
```

Note that matrix multiplication is different from element-wise multiplication. For element-wise, `*` is used.


```python
x = torch.ones(8)
print(x)
x = x - 1.5
print(x)
x -= 1.5
print(x)
x += torch.randn(1)
print(x)
x += torch.randn(8)
print(x)
```

Broadcasting works as in NumPy: [link](https://pytorch.org/docs/stable/notes/broadcasting.html)

Keep in mind, just like in NumPy, whatever you want to do with a tensor, there's probably an elegant operation for it implemented somewhere, you just have to look for it (on google and in the documentation)

# 2 Building a neural network in PyTorch

## 2.1 The building blocks

Neural networks are initialized through the use of class objects. You have encountered class objects already during this course: sklearn models are all class objects. The difference here is that we will code our own class first, before using it.

Many of the functionalities necessary to create [**all types of neural networks**](http://www.asimovinstitute.org/neural-network-zoo/) have [**already been implemented**](http://pytorch.org/docs/master/nn.html).

Let's inspect the most basic building blocks first: the [linear layer](https://pytorch.org/docs/master/generated/torch.nn.Linear.html#torch.nn.Linear) and the [ReLU](https://pytorch.org/docs/master/generated/torch.nn.ReLU.html#torch.nn.ReLU)

A linear layer is an object that will perform a matrix multiplication once called. Here, we instantiate such a layer with 20 input nodes and 40 output nodes:


```python
import torch.nn as nn
nn.Linear(20, 40)
```

Let's simulate some random data for this layer: A data set (or batch) with 16 samples and 20 features:


```python
x = torch.randn(16, 20)
```

Now let's use our linear layer on this data:


```python
layer = nn.Linear(20, 40)
print(x.shape)
z = layer(x)
print(z.shape)
```

What happens when we try to feed our layer an input with a different number of features?


```python
x = torch.randn(16, 30)
layer(x)
```

Let's see the ReLU in action:


```python
relu = nn.ReLU()

x = torch.randn(2, 4)
print(x)
z = relu(x)
print(z)
```

As you may have noticed, `nn.Module`s are class objects, a bit like scikit-learn models, that you instantiate and then call.

You can chain `nn.Module`s with the use of `nn.Sequential`:


```python
linear_and_relu = nn.Sequential(nn.Linear(20, 40), nn.ReLU())

x = torch.randn(16, 20)
z = linear_and_relu(x)
```

Or even longer constructs:

Always keep in mind what happens with the dimensions of your input and outputs with every layer, if your first layer outputs 40 features/nodes/hidden dimensions, then logically the next will have to take 40 as input.


```python
a_whole_damn_network = nn.Sequential(
    nn.Linear(128, 64),
    nn.ReLU(),
    nn.Linear(64, 16),
    nn.ReLU(),
    nn.Linear(16, 4),
    nn.ReLU(),
    nn.Linear(4, 1)
    )

x = torch.randn(16, 128)
z = a_whole_damn_network(x)
print(z.shape)
print(z)
```

The output $z$ that we now obtain after this whole network are called **logits**. They are the real-numbered $\mathbb{R}$ outputs that we obtain at the end of the network before our last activation function. This last activation function will be a, depending on the task at hand, sigmoid, softmax, or nothing at all for regression

Similar implementations exist for Dropout and Normalization layers in PyTorch, we invite you to look them up.

## 2.2 Class object neural networks and hyperparameters

We've seen how to implement a neural network using PyTorch `nn.Sequential`. It is more flexible however to write our own model class. This allows us to have more control over which operations we use and define our own hyperparameters. To make a PyTorch model, we specify our class object to be a submodule of `nn.Module` and inherit its methods via `super().__init__()`. Further, we specify all necessary attributes (such as layers) in our `__init__` function (executed upon initialization) and implement a `forward` function which will be executed when we call the object after being initialized.

The following code shows two examples, the first one of a very basic construction of a neural network without hyperparameters. The other one shows the same network, but where we set up our `__init__` function to process hyperparameters as input arguments. We can for example specify a hyperparameter whether we want to use dropout or not. (We could go even further to have an extra hyperparameter specifying the probability of dropout, ...).


```python
class BasicModel(nn.Module):
    def __init__(self):
        super().__init__()
        
        self.layer1 = nn.Linear(50, 40)
        self.relu1 = nn.ReLU()
        self.dropout1 = nn.Dropout(0.2)

        self.layer2 = nn.Linear(40, 20)
        self.relu2 = nn.ReLU()
        self.dropout2 = nn.Dropout(0.2)

        self.layer3 = nn.Linear(20, 10)
        self.relu3 = nn.ReLU()
        self.dropout3 = nn.Dropout(0.2)

        self.layer4 = nn.Linear(10, 5)
        # Think again: why do we not want a relu and dropout after our last layer again?

    def forward(self, x):
        # call them in separate lines:
        x = self.layer1(x)
        x = self.relu1(x)
        x = self.dropout1(x)

        # or together:
        x = self.dropout2(self.relu2(self.layer2(x)))
        x = self.dropout3(self.relu3(self.layer3(x)))

        # we could've also wrapped everything in a nn.Sequential ..

        x = self.layer4(x)
        return x

class HyperparameterModel(nn.Module):
    def __init__(self, dimensions_from_input_to_output = [50, 40, 20, 10, 5], dropout = True):
        super().__init__()

        layers = []
        # iterate through all layers:
        for i in range(len(dimensions_from_input_to_output) - 2):
            layer = nn.Linear(dimensions_from_input_to_output[i], dimensions_from_input_to_output[i + 1])
            layers.append(layer)
            layers.append(nn.ReLU())
            if dropout == True:
                layers.append(nn.Dropout(0.2))


        # the last layer separate from the loop because we don't want a ReLU and dropout after the last layer
        layer = nn.Linear(dimensions_from_input_to_output[i+1], dimensions_from_input_to_output[i + 2])
        layers.append(layer)

        # wrap the layers in a sequential
        self.net = nn.Sequential(*layers)

    def forward(self, x):
        return self.net(x)


```

Testing out the basic model:


```python
net = BasicModel()
net
```


```python
x = torch.randn(4, 50)
y = net(x)
print(y)
print(y.shape)

```

Testing out the hyperparameter model, when no arguments are specified the default ones are chosen:


```python
net = HyperparameterModel()
net
```


```python
x = torch.randn(4, 50)
y = net(x)
print(y)
print(y.shape)
```

Or by specifying hyperparameters, the following code shows a bit of a deeper model:


```python
net = HyperparameterModel(dimensions_from_input_to_output = [50, 160, 80, 40, 20, 10, 5], dropout = False)
net
```


```python
x = torch.randn(4, 50)
y = net(x)
print(y)
print(y.shape)
```

<div class="alert alert-success">

<b>EXERCISE:</b> 
<p> Copy over any of the two networks above (preferably the one you understand best), and modify it to take as input 784 features and as output 10 nodes. We will use this network later in the PC lab. You can test your network for bugs using some randomly generated data (as above).</p>

</div>


```python
######## YOUR CODE HERE #########


#################################
```

## 2.3 Data and training

You may have noticed that the outputs of your model have a `grad_fn` attribute. This grad function will be used by PyTorch automatic differentation engine to perform backward passes and compute gradients for every parameter with respect to the loss/cost function.

Now that we have our model and know that PyTorch will automatically compute gradients for us, we can move on to train a model. To do this, we need a couple more things:

The most basic blueprint of PyTorch model training consists of
- Get your data
- Wrap your data splits in a [data loader](https://pytorch.org/docs/stable/data.html)
- Instantiate the model
- Instantiate a [loss function](https://pytorch.org/docs/stable/nn.html#loss-functions)
- Instantiate an [optimizer object](https://pytorch.org/docs/stable/optim.html), to which you pass the parameters you want to optimize
- Iterate through your training data, for every batch:
    - reset the gradients
    - do forward pass
    - compute loss
    - backward pass
    - update parameters

(Optionally):
- After every full iteration through all training data samples (called an epoch), loop through all batches of validation data:
    - forward pass
    - compute loss and validation scores

In this way, we can monitor how good our model is performing on left-out validation data during training, this is used in early stopping. Notably, we do not compute and reset gradients and update parameters during our validation iterations, because we do not want to fit our model on this data.


Let's start with step one: get your data and wrap them in data loaders. We will first illustrate how we can convert our usual pandas or numpy datasets to be compatible with PyTorch with some random data:


```python
X_train = np.random.randn(100, 784)
y_train = np.random.randint(low = 0, high = 10, size = (100))
```


```python
X_train = torch.from_numpy(X_train)
print(X_train.dtype)
```

Remember to look at your data types: by default NumPy is `float64`, but if you instantiate a model, by default it will have weights in `float32`. It is therefore advised to convert your data to `float32`. In PyTorch, simply `float` is shorthand for `float32`.


```python
X_train = X_train.float()
# Equivalent: X_train.to(torch.float) or X_train.to(torch.float32)
print(X_train.dtype)
```


```python
y_train = torch.tensor(y_train)
print(y_train.dtype)
```

Now that we have our X_train and y_train as tensors, we can make them PyTorch-ready by wrapping them in a PyTorch dataset, and then by wrapping them in a DataLoader, which will create batches for us:


```python
train_dataset = torch.utils.data.TensorDataset(X_train, y_train)
train_dataloader = torch.utils.data.DataLoader(train_dataset, batch_size=16, pin_memory=True, shuffle=True)
```

Now we can use our train_dataloader as such:


```python
for batch in train_dataloader:
    X, y = batch
    print(X.shape, y.shape)
```

As you can see, we can iterate through our batches by use of a for loop, and it will spit out a training batch consisting of a list of X and y tensors. We can also test out code by isolating one training batch like this (only necessary for testing out code):


```python
batch = next(iter(train_dataloader))
batch
```

For the remainder of the PC lab, we will work with the MNIST dataset included in `torchvision`. (If you're running this code locally, you may have to pip install torchvision).


```python
from torchvision import datasets
from torchvision.transforms import ToTensor
train_data = datasets.MNIST(
    root = 'data',
    train = True,                         
    transform = ToTensor(), 
    download = True,            
)
test_data = datasets.MNIST(
    root = 'data', 
    train = False, 
    transform = ToTensor()
)

X_train = train_data.data
y_train = train_data.targets

X_test = test_data.data
y_test = test_data.targets
```


```python
print('shapes:')
print(X_train.shape, y_train.shape, X_test.shape, y_test.shape)
print('first training image tensor:')
print(X_train[0])
print('first five labels:')
print(y_train[:5])
```

Our data is images, each data sample has $28 \times 28$ input features, signifying the pixels. In order to feed this data to our model, we will need to flatten these features:


```python
X_train = X_train.reshape(-1, 28 * 28)
X_test = X_test.reshape(-1, 28 * 28)
```

In addition, the grayscale values of our images go from 0 to 255. It is perhaps good practice to min-max standardize these numbers by dividing through 255:


```python
X_train = X_train / 255
X_test = X_test / 255
```

Finally, let's check our datatypes to see if everything is looking good to go:


```python
X_train.dtype, X_test.dtype, y_train.dtype, y_test.dtype
```

Let's split up our training set in a training and validation set and finally wrap our data in a data loader:


```python
np.random.seed(42)
train_indices, val_indices = np.split(np.random.permutation(len(X_train)), [int(len(X_train)*0.8)])
X_val = X_train[val_indices]
y_val = y_train[val_indices]
X_train = X_train[train_indices]
y_train = y_train[train_indices]

train_dataset = torch.utils.data.TensorDataset(X_train, y_train)
train_dataloader = torch.utils.data.DataLoader(train_dataset, batch_size=16, pin_memory=True, shuffle=True)

val_dataset = torch.utils.data.TensorDataset(X_val, y_val)
val_dataloader = torch.utils.data.DataLoader(val_dataset, batch_size=16, pin_memory=True, shuffle=True)

test_dataset = torch.utils.data.TensorDataset(X_test, y_test)
test_dataloader = torch.utils.data.DataLoader(test_dataset, batch_size=16, pin_memory=True, shuffle=True)
```

Let's visualize a random batch:


```python
batch = next(iter(train_dataloader))
X_batch, y_batch = batch
import matplotlib.pyplot as plt


figure = plt.figure(figsize=(10, 8))
cols, rows = 4, 4
for i in range(cols * rows):
    img, label = X_batch[i], y_batch[i]
    figure.add_subplot(rows, cols, i+1)
    plt.title(label.item())
    plt.axis("off")
    plt.imshow(img.reshape(-1, 28, 28).squeeze(), cmap="gray")
plt.show()
```

Now that we have our data ready, let's reiterate our PyTorch model blueprint:

The most basic blueprint of PyTorch model training consists of
- Get your data
- Wrap your data splits in a [data loader](https://pytorch.org/docs/stable/data.html)
- Instantiate the model
- Instantiate a [loss function](https://pytorch.org/docs/stable/nn.html#loss-functions)
- Instantiate an [optimizer object](https://pytorch.org/docs/stable/optim.html), to which you pass the parameters you want to optimize
- Iterate through your training data, for every batch:
    - reset the gradients
    - do forward pass
    - compute loss
    - backward pass
    - update parameters

(Optionally):
- After every full iteration through all training data samples (called an epoch), loop through all batches of validation data:
    - forward pass
    - compute loss and validation scores


In one of the previous exercises, we already implemented a model compatible with MNIST: 784 input features and 10 output nodes (one for each class). Hence, we can move on to loss function and optimizers. For multi-class classification of the digits, we will need to use the [Cross Entropy Loss](https://pytorch.org/docs/stable/generated/torch.nn.CrossEntropyLoss.html#torch.nn.CrossEntropyLoss). Take a look at what kind of inputs this loss function expects. According to the documentation: `The input is expected to contain raw, unnormalized scores for each class.` Meaning that we can pass logits directly to this loss function, and we do not have to apply a softmax operation ourselves.

For optimizer choice, we can choose [vanilla gradient descent](https://pytorch.org/docs/stable/generated/torch.optim.SGD.html#torch.optim.SGD) or the nowadays-industry-standard [Adam optimizer](https://pytorch.org/docs/stable/generated/torch.optim.Adam.html#torch.optim.Adam), which is itself a pimped version of stochastic gradient descent (with momentum).

Take note that we can also specify our desired learning rate in our model. This learning should almost always be tuned as it will influence how fast our model trains but also convergence.


```python
model = # your model from previous exercises here.

loss_function = nn.CrossEntropyLoss()

optimizer = torch.optim.SGD(model.parameters(), lr=0.005) # SGD = stochastic gradient descent
```

Now we're ready to perform training. We'll build up the training loop step-wise in the following codeblocks. Let's first start with passing one batch to the model, computing the loss, performing backpropagation and updating the weights:


```python
batch = next(iter(train_dataloader))
X_batch, y_batch = batch

y_hat_batch = model(X_batch)

loss = loss_function(y_hat_batch, y_batch) # Compute loss

loss.backward()   # Calculate gradients
optimizer.step()   # Update weights using defined optimizer

print(X_batch.shape, y_batch.shape)
print(y_hat_batch.shape)
print("Outputs as logits, first two samples:")
print(y_hat_batch[:2])
print('loss:', loss)
```

Everytime we perform a training step, we should reset the gradients so that the gradients computed on the previous batch do not influence the next. We can do this by calling `.zero_grad()` on our optimizer. In practice, it is most safe to call this before every forward pass. So if we train a complete epoch:


```python
all_losses = []

for batch in train_dataloader:
    optimizer.zero_grad()
    X_batch, y_batch = batch

    y_hat_batch = model(X_batch)

    loss = loss_function(y_hat_batch, y_batch) # Compute loss

    loss.backward()   # Calculate gradients
    optimizer.step()   # Update weights using defined optimizer

    all_losses.append(loss.item())
```

Plotting the loss function of every batch during one epoch:


```python
plt.plot(np.arange(len(all_losses)), all_losses)
smoothed_losses = np.convolve(all_losses, np.ones(50)/50, mode = "valid")
plt.plot(np.arange(len(smoothed_losses)), smoothed_losses)
```

We have evaluated the training progress during one epoch on our training set. What if we want to see the performance of the validation set during training, so that we can see if/when the model starts overfitting? It is common practice to perform a pass through the whole validation set after every training epoch. However, we should not compute gradients now, we can block automatic gradient computation using `torch.no_grad()`. Also, we need to tell PyTorch to put our model in evaluation mode (so it does not perform dropout anymore, etc...). After the validation epoch, we should put the model back to training mode again. The code should look something like this:


```python
predictions = []
true_labels = []
losses = []

model.eval()
with torch.no_grad():
    for batch in val_dataloader:
        X_batch, y_batch = batch

        y_hat_batch = model(X_batch)

        loss = loss_function(y_hat_batch, y_batch)

        losses.append(loss.item())
        predictions.append(y_hat_batch)
        true_labels.append(y_batch)

model.train()

predictions = torch.cat(predictions)
true_labels = torch.cat(true_labels)
accuracy = (true_labels == predictions.argmax(-1)).sum().item() / len(predictions)

print(accuracy)
print(np.mean(losses))
```

<div class="alert alert-success">

<b>EXERCISE:</b> 
<p> Using the code above, put it all together to train a model for multiple epochs. After every epoch, print or save some training and validation statistics.

Monitor how good your model is training. What things could you change? In particular, try the Adam optimizer or try using gradient descent with the momentum argument. How does this influence training speed? What is it doing? Try tweaking the learning rate.

You should be able to obtain an accuracy of +- 96% </p>

</div>


```python
N_EPOCHS = 20

model =  # your model from previous exercises here.

# loss function & optimizer

for i in range(1, N_EPOCHS + 1):
    
    # train loop

    # eval loop

    # record or print some variables that you want to keep track of during training
```

<div class="alert alert-success">

<b>EXERCISE:</b> 
<p> After training, what is your performance on the test dataset? </p>

</div>


```python
### YOUR CODE HERE ###

######################
```

## 2.4 Extra: The autoencoder

An autoencoder tries to reconstruct its own input. By constraining the dimensionality of our hidden layers (the bottleneck), the model is forced to learn a low-dimensional manifold of the data, much like PCA does (in fact, *architecturally* there are only two differences between PCA and an autoencoder: (1) PCA is linear whereas an autoencoder can be arbitrarily complex and (2) the bottleneck in PCA is constrained to be orthogonal, whereas in an autoencoder no such constraints are in place. Other than that, they are also optimized differently.)

Terminology-wise, we can split up our neural network into two parts: the encoder, which encodes the data to a lower-dimensional space, and the decoder, which decodes that encoding back to the original sample space. During training, samples go through both encoder and decoder. After training, we can opt to use different parts of the autoencoder: if we are interested in the lower-dimensional space of samples (for e.g. clustering or visualization), we only use the encoder. If we want to denoise samples (look up denoising autoencoder), we are interested in both encoder and decoder. Finally, if we want to generate samples, we can start from randomly generated vectors in the bottleneck space (e.g. random 10 dimensional vectors) and put them through the decoder to generate new samples (look up variational autoencoders).

It is very easy to implement a simple autoencoder in PyTorch, especially with the code that we have already written. The only things we need to change are:
- The architecture of the network should be changed so that it represents an hourglass with a bottleneck (most easy for visualization is if you have a bottleneck with 2 hidden dimensions, so that you can plot it directly, otherwise if you have e.g. a bottleneck with 10 hidden dimensions, then you need to do an additional PCA or t-SNE first on the samples). A symmetrical network in terms of hidden dimensions is most conventional. 
- The loss function should be adjusted to a regression based loss (such as MSE), since our input/output pixels are continuous as well.
- The loss function should be called on the original input, instead of the target classes
- Accuracy can not be computed on reconstruction, only loss.

The following shows how you can easily make an autoencoder:


```python
model = nn.Sequential(
    HyperparameterModel(dimensions_from_input_to_output= [784, 256, 128, 64]),
    HyperparameterModel(dimensions_from_input_to_output= [64, 128, 256, 784])
)

# Alternatively, we could have made our new class for autoencoders
# or we could have just specified one HyperparameterModel passing [784, 128, 2, 128, 784] as argument
# The construct with the Sequential has the advantage that it allows us to call the encoder and decoder separately
# as model[0](input_data) and model[1](encoded_data)
```

<div class="alert alert-success">

<b>EXTRA EXERCISE:</b> 
<p> Train this autoencoder model by reusing the code used in the supervised learning task and adjusting code where necessary, for a to-do list of necessary changes, consult the list above. </p>

</div>


```python
N_EPOCHS = 20

model =  # your model from previous exercises here.

# loss function & optimizer

for i in range(1, N_EPOCHS + 1):
    
    # train loop

    # eval loop

    # record or print some variables that you want to keep track of during training
```

We can now visualize the manifold our model has learned. In order to do so, we encode samples to our $x$-dimensional bottleneck. After which we can use t-SNE to give us a final 2-dimensional space for visualization. Does our model encode digit identity well in its "latent space"?
What about if we fit t-SNE directly on the pixel values? Does our autoencoder noticeably add value? **Note**: t-SNE may take a couple minutes to run.


```python
inputs = []
encoded_inputs = []
classes = []
with torch.no_grad():
    model.eval()
    for batch in val_dataloader:
        X_batch, y_batch = batch
        inputs.append(X_batch)
        encoded_inputs.append(model[0](X_batch))
        classes.append(y_batch)

classes = torch.cat(classes).numpy()
inputs = torch.cat(inputs).numpy()
encoded_inputs = torch.cat(encoded_inputs).numpy()
```


```python
from sklearn.manifold import TSNE
import seaborn as sns

encoded_tsne = TSNE(verbose = 10, init = "pca", learning_rate = "auto").fit_transform(encoded_inputs)
plt.figure(figsize = (10, 10))
sns.scatterplot(x = encoded_tsne[:,0], y = encoded_tsne[:, 1], hue = classes, alpha = .75, palette="deep")
plt.show()
```


```python
encoded_tsne = TSNE(verbose = 10, init = "pca", learning_rate = "auto").fit_transform(inputs)
plt.figure(figsize = (10, 10))
sns.scatterplot(x = encoded_tsne[:,0], y = encoded_tsne[:, 1], hue = classes, alpha = .75, palette="deep")
```

Our autoencoder may not add much value visually on top of just fitting t-SNE directly on pixel values, but the advantage of an autoencoder do not only come from visualization. For example, we can reconstruct samples from the latent space back to the original space with the decoder. As a last, let's visualize some reconstructions:


```python
with torch.no_grad():
    batch = next(iter(val_dataloader))
    X_batch, y_batch = batch
    model.eval()
    y_hat_batch = model(X_batch)


figure = plt.figure(figsize=(10, 8))
cols, rows = 4, 4
for i in range(cols * rows):
    img, label = X_batch[i], y_batch[i]
    figure.add_subplot(rows, cols, i+1)
    plt.title(label.item())
    plt.axis("off")
    plt.imshow(img.reshape(-1, 28, 28).squeeze(), cmap="gray")
plt.show()


figure = plt.figure(figsize=(10, 8))
cols, rows = 4, 4
for i in range(cols * rows):
    img, label = y_hat_batch[i], y_batch[i]
    figure.add_subplot(rows, cols, i+1)
    plt.title(label.item())
    plt.axis("off")
    plt.imshow(img.reshape(-1, 28, 28).squeeze(), cmap="gray")
plt.show()
```

These are just samples generated from an encoding of a true original sample. What if we want to generate completely new synthetic digits?
For this we need to sample vectors in our latent space.
Naively we could just sample vectors from a random gaussian.


```python
samples = torch.randn(16, 64)

with torch.no_grad():
    model.eval()
    decoded_samples = model[1](samples)


figure = plt.figure(figsize=(10, 8))
cols, rows = 4, 4
for i in range(cols * rows):
    img = decoded_samples[i]
    figure.add_subplot(rows, cols, i+1)
    plt.axis("off")
    plt.imshow(img.reshape(-1, 28, 28).squeeze(), cmap="gray")
plt.show()
```

The images vaguely look like something that could be digits, but are not convincing to say the least.
One way to make these generations better is to not sample random vectors, but try to inform our sampling. For this dataset, we have class labels, so we can try to see what distribution samples from a certain class (e.g. 7) generally follow in the latent space. We can then sample from that distribution, instead of from a standard gaussian.

To do this: we (1) first encode samples in the latent space, (2) take the samples belonging to the cluster we're interested in, (3) fit a distribution to the latent space: in this case we will simply compute the mean and standard deviation of all hidden dimensions, and (4) finally, sample some vectors from that distribution.


```python
class_ = 7

encoded_inputs = []
classes = []
with torch.no_grad():
    model.eval()
    for batch in val_dataloader:
        X_batch, y_batch = batch
        encoded_inputs.append(model[0](X_batch))
        classes.append(y_batch)

classes = torch.cat(classes)
encoded_inputs = torch.cat(encoded_inputs)

encoded_inputs_class = encoded_inputs[classes == class_]

mu = encoded_inputs_class.mean(0)
std = encoded_inputs_class.std(0)
print(mu)
print(std)

figure = plt.figure(figsize=(10, 8))
cols, rows = 4, 4
for i in range(cols * rows):
    sampled_vector = torch.normal(mu, std).unsqueeze(0)
    with torch.no_grad():
        model.eval()
        decoded_sample = model[1](sampled_vector)
    img = decoded_sample
    figure.add_subplot(rows, cols, i+1)
    plt.axis("off")
    plt.imshow(img.reshape(-1, 28, 28).squeeze(), cmap="gray")
plt.show()

```

Looks like conditioning on a class-specific distribution helps the model along a bit, but it is still clear that not every sampled vector taken from the latent space gives us a convincing decoded image. An advanced technique based on autoencoders to give us better image generation quality is the Variational Autoencoder, which is out of scope for this course. Another (non-autoencoder) based technique is GANs.

### Extra: Using GPUs

Matrix multiplication run orders of magnitude faster on GPU hardware. If you have a local GPU and have installed a PyTorch version with GPU, you should be able to run this code locally. Otherwise, In Google Colab, you can request access to a GPU via `Runtime > Change runtime type > Hardware accelerator = GPU`.

Briefly, the steps needed to train on GPUs consist of
1. Putting your model on the GPU
2. During your training loop, putting every batch on the GPU before the forward pass
3. If you have a validation loop, doing the same there for every batch.
4. If you have variables that you will use after training (e.g. predictions on the validation set), remember to return this back to the CPU, as GPUs have limited memory.

In PyTorch, we put variables and models on the GPU by specifying their 'device' to be 'cuda' (the parallel computing platform for nvidia GPUs that PyTorch uses).

The following code illustrates how to train your models on GPU hardware:



```python
X, y = next(iter(train_dataloader))

model = model.to('cuda')
print(X.device)
X = X.to('cuda')
print(X.device)

y_hat = model(X)

y_hat
```

We encourage you to try out training on GPUs during the next PC lab(s)

6a21b4e234b233ba0240d57b35b872efd072163d
ipynb
Jupyter Notebook
predmod/ANN_intro/PClab12_ANN.ipynb
gdewael/teaching
a78155041918422a843f31c863dd11e8afc5646a
MIT
predmod/ANN_intro/PClab12_ANN.ipynb
gdewael/teaching
a78155041918422a843f31c863dd11e8afc5646a
MIT
predmod/ANN_intro/PClab12_ANN.ipynb
gdewael/teaching
a78155041918422a843f31c863dd11e8afc5646a
MIT
Qwen/Qwen-72B
1. YES
2. YES