# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# **These code cells were run on google colab with GPU support.**

# + colab={"base_uri": "https://localhost:8080/"} id="GrQx0ROdcsh1" outputId="6dbb4f19-1b71-4432-b67d-62895f56008f"
import torch
if torch.cuda.is_available():
  device = torch.device("cuda:0")
  print("GPU")
else:
  device = torch.device("cpu")
  print("CPU")
# -

# **Fetching the data**

# + colab={"base_uri": "https://localhost:8080/"} id="hhI29xmxi9T_" outputId="380952c9-4258-4c7c-ab8e-9b20ef729faa"
from sklearn.datasets import fetch_openml
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np

mnist = fetch_openml("mnist_784")
print(mnist.keys())

# -

# **Plotting one example of each class**

# + colab={"base_uri": "https://localhost:8080/", "height": 224} id="3IZkcerWjEsP" outputId="25babe09-1336-42eb-d035-e4d4cb242f26"
X = mnist["data"]
Y = mnist["target"]
Y = Y.astype(int)
X=(X/255 - 0.5)*2

fig, ax = plt.subplots(nrows=2, ncols=5, sharex=True, sharey=True)
ax = ax.flatten()

for i in range(10):
    for x, y in zip(X, Y):
        if y==i:
            img=np.array(x).reshape((28,28))
            ax[i].imshow(img, cmap="Greys")
            break
ax[0].set_yticks([])
ax[0].set_xticks([])
plt.tight_layout()
plt.show()
# -

# **Out of 70k examples, 10k will be used for test set, and remaining will be used for training and validation.**

# + id="6TzuxPwWjJi9"
X_train, X_test, Y_train, Y_test = X[:60000], X[60000:], Y[:60000], Y[60000:]

# + id="5ENyvH4kjKOj"
from sklearn.model_selection import StratifiedShuffleSplit
split = StratifiedShuffleSplit(n_splits=1, test_size=0.1, random_state=42)
for train_index, val_index in split.split(X_train,Y_train):
    X_train_strat = X_train[train_index, :]
    Y_train_strat = Y_train[train_index]
    X_val_strat = X_train[val_index, :]
    Y_val_strat = Y_train[val_index]

# + id="panFa_57kQhw"
import torch.nn as nn
import torch.nn.functional as Func
from torch.autograd import Variable
import torch.optim as optim
import torch.utils.data as data
import random
from scipy.io import savemat
import os
from os import path
from sklearn.preprocessing import normalize
from torch.nn.utils import clip_grad_norm_
import torch.nn.parallel.data_parallel as data_parallel
from sklearn.metrics import confusion_matrix

# -

# **Extending and overriding methods for our own dataset**

# + id="k19sQHq3kVKe"

class mnist_dataset(data.Dataset):
    def label_transformer(self, labels):
        return labels

    def __init__(self, input_data, labels):
        input_data = input_data.reshape((len(input_data),1,28,28))
        self.feats = input_data
        self.labels = self.label_transformer(labels)
    def __len__(self):
        return len(self.labels)
    def __getitem__(self, index):
        x = self.feats[index]
        y = self.labels[index]
        return x,y



# -

# **Creating dataloader for each of the train, validation, test dataset.**

# + id="9_CFo7qQkX9Z"
class hyperparam:
    bs = 100
    lr = 0.05
    num_epochs = 50

params = {
    "batch_size": hyperparam.bs,
    "shuffle": True,
    "num_workers": 2,
    "drop_last": False,
    "pin_memory": True
}

train_set = mnist_dataset(X_train_strat, Y_train_strat)
val_set = mnist_dataset(X_val_strat, Y_val_strat)
test_set = mnist_dataset(X_test, Y_test)

training_gen = data.DataLoader(train_set, **params)
val_gen = data.DataLoader(val_set, **params)
test_gen = data.DataLoader(test_set, **params)
# -

# **Created a DNN with 2 layers of CNN with 12 filters each and adding two fully connected layers of 100 and 10 neurons respectively.
# Used Relu activation function, with initial learning rate = 0.05 with glorot initialization.**

# + id="GgTVa1ckklJ9"
from torch.nn import Conv2d, Linear
from torch import flatten

class optim_cnn(nn.Module):
    
    def glorot_initialize(self, layers):
        for layer in layers:
            torch.nn.init.xavier_normal_(layer.weight)
            torch.nn.init.zeros_(layer.bias)
    
    def __init__(self):
        super(optim_cnn, self).__init__()
        self.conv1 = Conv2d(1,12,kernel_size=(3,3), padding = 1)
        self.conv2 = Conv2d(12,12,kernel_size=(3,3), padding = 1)
        self.fc1 = Linear(588, 100)
        self.fc2 = Linear(100, 10)

        self.glorot_initialize([self.conv1, self.conv2, self.fc1, self.fc2])
        
    def forward(self, sig):
        sig = Func.max_pool2d(Func.relu(self.conv1(sig)), (2, 2))
        sig = Func.max_pool2d(Func.relu(self.conv2(sig)), (2, 2))
        sig = sig.view(-1, 12*7*7)
        sig = Func.relu(self.fc1(sig))
        sig = self.fc2(sig)
        return sig
        # return Func.softmax(sig, dim = 1)



# -

# **We have created three models, each of these will be trained using a different optimizer.**

# + id="qVmZYmAEGHy7"
cnn_models = [optim_cnn().to(device), optim_cnn().to(device), optim_cnn().to(device)]
# -

# **Model is trained for 50 epochs, after each epochs, printing the validation accuracy. and resulting learning rate after adjusting learning rate by 10% each 10 epochs
# Also using early stopping mechanism, which stops the learning if the validation accuracy starts dropping for a consecutive 5 cycles. This is done to prevent overfitting.**

# + colab={"base_uri": "https://localhost:8080/"} id="wov0a1kVkoBq" outputId="70d7479d-7e08-4027-a664-a4221bcfca72"
from tqdm import tqdm
from datetime import datetime
from torch.optim.lr_scheduler import StepLR
from torch.optim import Adam, RMSprop

tr_avg_loss_list = {0: [], 1:[], 2:[]}
tr_accuracy_list = {0: [], 1:[], 2:[]}
val_avg_loss_list = {0: [], 1:[], 2:[]}
val_accuracy_list = {0: [], 1:[], 2:[]}

print(datetime.now())

def get_optimizer(model_num, model):
    if model_num==0:
        return RMSprop(model.parameters(), lr = 0.001, alpha = 0.9)
    elif model_num == 1:
        return torch.optim.SGD(model.parameters(), lr = 0.05, momentum=0.9, nesterov=True) # TODO Nesterov
    elif model_num == 2:
        return Adam(model.parameters(), lr = 0.001, eps = 1e-8, weight_decay=0)

for model_num, cnn_model in enumerate(cnn_models):

    optimizer = get_optimizer(model_num, cnn_model)

    if model_num == 1:
        scheduler = StepLR(optimizer, step_size=10, gamma=0.9)

    loss = nn.CrossEntropyLoss()

    for epoch in range(hyperparam.num_epochs):
        print("Epoch:" + str(epoch) + " model num: " + str(model_num+1))
        tr_num_correct = 0
        tr_num_samples = 0
        tr_total_loss = 0
        val_num_correct = 0
        val_num_samples = 0
        val_total_loss = 0
        
        print("Learning rate: " + str(optimizer.param_groups[0]['lr']))

        with torch.set_grad_enabled(True):
            cnn_model.train(True)
            for ind, (local_batch, local_labels) in enumerate(training_gen):
                
                optimizer.zero_grad()
                local_batch = local_batch
                local_labels = local_labels
                local_batch, local_labels = Variable(local_batch).float(), Variable(local_labels)
                local_batch = local_batch.to(device)
                local_labels = local_labels.to(device)
                out1 = cnn_model(local_batch)
                ploss = loss(out1, local_labels.long())
                tr_total_loss += ploss * hyperparam.bs
                ploss.backward()
                optimizer.step()
                sel_class = torch.argmax(out1, dim=1)
                tr_num_correct += sel_class.eq(local_labels).sum().item()
                tr_num_samples += hyperparam.bs
                    
                                    
            tr_avg_loss = tr_total_loss / len(training_gen.dataset)
            tr_avg_loss_list[model_num].append(tr_avg_loss)
            tr_accuracy = tr_num_correct / len(training_gen.dataset)
            tr_accuracy_list[model_num].append(tr_accuracy)
        
        with torch.set_grad_enabled(False):
            cnn_model.eval()
            for local_batch, local_labels in val_gen:
                local_batch = local_batch.float()
                local_labels = local_labels.float()
                local_batch, local_labels = Variable(local_batch), Variable(local_labels)
                local_batch = local_batch.to(device)
                local_labels = local_labels.to(device)
                out1 = cnn_model(local_batch)
                ploss = loss(out1, local_labels.long())
                val_total_loss += ploss * hyperparam.bs
                sel_class = torch.argmax(out1, dim=1)
                val_num_correct += sel_class.eq(local_labels).sum().item()
                val_num_samples += local_labels.size(0)

            val_avg_loss = val_total_loss / len(val_gen.dataset)
            val_avg_loss_list[model_num].append(val_avg_loss)
            val_accuracy = val_num_correct / len(val_gen.dataset)
            print("Validation accuracy: " + str(val_accuracy))
            val_accuracy_list[model_num].append(val_accuracy)

            if model_num == 1:
                scheduler.step()

            if epoch > 10:
                if sum([val_accuracy_list[model_num][i] < val_accuracy_list[model_num][i-1] for i in range(epoch-5, epoch)]) == 5:
                    break


# -

# **Plotting learning curves for validation and train dataset**

# + colab={"base_uri": "https://localhost:8080/", "height": 1000} id="Wqimq2CM8TOb" outputId="ad693855-7932-4885-e054-ca4428f043ae"
def plot_x_y_vals(x_vals, y_vals, x_label, y_label, label, line_titles):
    for i in range(len(x_vals)):
        plt.plot(x_vals[i], y_vals[i], label=line_titles[i])
    plt.title(label)
    plt.legend()
    plt.xlabel(x_label)
    plt.ylabel(y_label)
    plt.show()

model_names = ["RMSProp","Nesterov", "Adam"]

for i in tr_accuracy_list:
    epocs = [i+1 for i in range(len(tr_accuracy_list[i]))]
    plot_x_y_vals([epocs, epocs], [tr_accuracy_list[i], val_accuracy_list[i]], "Epochs", "Accuracy", "Train & Validation Accuracy, " + model_names[i], ["train", "validation"])
    plot_x_y_vals([epocs, epocs], [tr_avg_loss_list[i], val_avg_loss_list[i]], "Epochs", "Loss", "Train & Validation Loss, " + model_names[i], ["train", "validation"])
# -

# **Overfit or Underfit?**
#
# **RMSPros** , we see that training accuracy is almost 1, while validation is around .985. We see that the validation loss increases with epochs towards the end of training. And we see that the validation accuracy is also not as high as baseline. We can classify this as **mild overfit**.
#
# **Nesterov**, we see that training accuracy is again almost 1, while validation is around 0.99. We also see that validation loss increases towards the end, but not as much. Overall this performs similar to baseline model. However the loss is not increasing as much as baseline model, we can say that it is just right fit. We could have stopped training at around 25 epochs as well. 
#
# **Adam**, We see that validation loss doesnt increase with the epochs towards the end. The accuracy doesnt increase as much, but doesnt increase as well. It is very similar to baseline model. We can say that comparatively this model is similar to baseline. So no relative overfit or underfit compared to baseline.
#

# **Checking the accuracy of test set**

# + id="dR8RFZkglxKX"
total_accurate = [0,0,0]
total_values = [0,0,0]
errors={0:{i:{j:0 for j in range(10)} for i in range(10)}, 1:{i:{j:0 for j in range(10)} for i in range(10)}, 2:{i:{j:0 for j in range(10)} for i in range(10)}}
incorrect_samples = {0:[], 1:[], 2:[]}
correct_samples = {0:[], 1:[], 2:[]}

def calculate_class_wise_errors(local_labels, sel_class, local_batch, model_num):
    true_labels = local_labels[sel_class.not_equal(local_labels)]
    predicted = sel_class[sel_class.not_equal(local_labels)]
    for (i, t), (i,p) in zip(enumerate(true_labels), enumerate(predicted)):
        errors[model_num][t.item()][p.item()] += 1
    
    true_labels = local_labels[sel_class.eq(local_labels)]
    predicted = sel_class[sel_class.eq(local_labels)]
    for (i, t), (i,p) in zip(enumerate(true_labels), enumerate(predicted)):
        errors[model_num][t.item()][p.item()] += 1

    if len(incorrect_samples[model_num]) < 10:
        samples = local_batch[sel_class.not_equal(local_labels)]
        predicted = sel_class[sel_class.not_equal(local_labels)]
        true_labels = local_labels[sel_class.not_equal(local_labels)]
        for (i,s), (i,p), (i, t) in zip(enumerate(samples), enumerate(predicted), enumerate(true_labels)):
            incorrect_samples[model_num].append((s.cpu().numpy(), p.cpu().numpy(), t.cpu().numpy()))

    if len(correct_samples[model_num]) < 10:
        samples = local_batch[sel_class.eq(local_labels)]
        predicted = sel_class[sel_class.eq(local_labels)]
        true_labels = local_labels[sel_class.eq(local_labels)]
        for (i,s), (i,p), (i, t) in zip(enumerate(samples), enumerate(predicted), enumerate(true_labels)):
            correct_samples[model_num].append((s.cpu().numpy(), p.cpu().numpy(), t.cpu().numpy()))

for model_num, cnn_model in enumerate(cnn_models):
    with torch.set_grad_enabled(False):
        cnn_model.eval()
        for local_batch, local_labels in test_gen:
            local_batch = local_batch.float()
            local_labels = local_labels.float()
            local_batch, local_labels = Variable(local_batch), Variable(local_labels)
            local_batch = local_batch.to(device)
            local_labels = local_labels.to(device)
            out1 = cnn_model(local_batch)
            ploss = loss(out1, local_labels.long())
            sel_class = torch.argmax(out1, dim=1)
            calculate_class_wise_errors(local_labels, sel_class, local_batch, model_num)
            total_accurate[model_num] += sel_class.eq(local_labels).sum().item()
            total_values[model_num] += local_labels.size(0)

    

# + colab={"base_uri": "https://localhost:8080/"} id="pshZ1lENmC1V" outputId="89bcb7b2-5247-43b6-9aa2-793de3cb53d8"
print("Predicted " + str(total_accurate) +" correctly out of " + str(total_values) + "for respective models: " + str(model_names))

# -

# **We see that all three models perform well, and ther is very small margin in accuracy. We see that Nesterov has performed the best while Adam is second and RMSProp has performed worst, but not by a lot.**
#
# **Below we are plotting the heatmap, where y-axis represents the actual label and x-axis represents the predicted labels. Pleaase mind that all the diagonal elements have been set to zero. So the heatmap only represents the incorrectly classified label counts.**
#
# **For example row = 4, col = 3 represents the count of images, which were 4, but were actually classified as 3. And the cell (4,4) is left empty, although it should ideally contain count of all the correctly classified images of 4. As this heatmap is generated only to see if there is any pair that is mistaken a lot in the classification OR if there is any bias in our model for any label.** 
#

# + colab={"base_uri": "https://localhost:8080/", "height": 926} id="jr3x3rR1Drwk" outputId="2d2f30a4-e261-4614-cba4-34241449545d"
import seaborn as sns
class_acc = np.zeros((3,10,10))
for d in range(3):
    for i in range(10):
        for j in range(10):
            if i!=j:
                class_acc[d,i,j] = errors[d][i][j]
            else:
                class_acc[d,i,j] = 0
    print("\n\n"+model_names[d])
    sns.heatmap(class_acc[d])
    plt.show()
# -

# **Plotting few correctly classified images by our model.**

# + colab={"base_uri": "https://localhost:8080/", "height": 880} id="oagP_kqzcllW" outputId="31179f2c-9731-4efa-d1db-60059449cae5"
print("Incorrectly classified samples. (True and predicted values)")

for d in range(3):
    print("\n\n"+model_names[d])
    fig, ax = plt.subplots(nrows=2, ncols=5, sharex=True, sharey=True)
    ax = ax.flatten()
    for i in range(10):
        img=np.array(incorrect_samples[d][i][0]).reshape(28,28)
        ax[i].imshow(img, cmap="Greys")
        ax[i].title.set_text(str(int(incorrect_samples[d][i][2])) + "-" + str(incorrect_samples[d][i][1]))
    ax[0].set_yticks([])
    ax[0].set_xticks([])
    plt.tight_layout()
    plt.show()

# -

# **Plotting few correctly classified images by our model.**

# + colab={"base_uri": "https://localhost:8080/", "height": 880} id="GClPJXDjoFDD" outputId="9eb96aba-4405-4730-996e-bd0f9d196bfe"

print("Correctly classified samples. (true and predicted values)")
for d in range(3):
    print("\n\n"+model_names[d])
    fig, ax = plt.subplots(nrows=2, ncols=5, sharex=True, sharey=True)
    ax = ax.flatten()
    for i in range(10):
        img=np.array(correct_samples[d][i][0]).reshape(28,28)
        ax[i].imshow(img, cmap="Greys")
        ax[i].title.set_text(str(int(correct_samples[d][i][2])) + "-" + str(correct_samples[d][i][1]))
    ax[0].set_yticks([])
    ax[0].set_xticks([])
    plt.tight_layout()
    plt.show()


# -

# **Plotting below confusion matrix for each class. By each class, we mean all the correct prediction fot that class = true positive. all the images of a class, that were incorrectly classfied as false negative. All the images, not of that class, but classified as of that class as false positive. And all the images that were not of a class and were also classified as not belonging to that class as true negative.**

# + colab={"base_uri": "https://localhost:8080/", "height": 376} id="9gm8Z_bN5r7q" outputId="9246d32f-86f5-429a-9f16-ee5a23e23a9b"
import pandas as pd
# Confusion matrix
confusion_arr = np.zeros((3, 10, 4))
confusion_dfs = []
for d in range(3):
    for i in range(10):
        confusion_arr[d][i][0] = errors[d][i][i] # tp
        for j in range(10):
            if i!=j:
                confusion_arr[d][i][1]+=errors[d][j][i] # fp
        for j in range(10):
            if i!=j:
                confusion_arr[d][i][2]+= errors[d][i][j] # fn
        confusion_arr[d][i][3] = total_values[d] - sum(confusion_arr[d][i][:3]) # tn
    
    confusion_dfs.append(pd.DataFrame(confusion_arr[d], columns=["tp", "fp", "fn", "tn"]))
    confusion_dfs[d]["precision"] = confusion_dfs[d]["tp"] / (confusion_dfs[d]["tp"] + confusion_dfs[d]["fp"])
    confusion_dfs[d]["recall"] = confusion_dfs[d]["tp"] / (confusion_dfs[d]["tp"] + confusion_dfs[d]["fn"])
    confusion_dfs[d]["accuracy"] = (confusion_dfs[d]["tp"] + confusion_dfs[d]["tn"]) / 10000


print("Overall Accuracy:" + str(total_accurate[0]/total_values[0]) + " For optimizer:" + model_names[0])
confusion_dfs[0]

# + colab={"base_uri": "https://localhost:8080/", "height": 376} id="SjYlJLRurviz" outputId="48225745-1808-47ec-ab8b-078424911b7e"
print("Overall Accuracy:" + str(total_accurate[1]/total_values[1]) + " For optimizer:" + model_names[1])
confusion_dfs[1]

# + colab={"base_uri": "https://localhost:8080/", "height": 376} id="CJcPPMGar7J9" outputId="2271ea25-6ea1-46f1-b8b2-218cbb8aab50"
print("Overall Accuracy:" + str(total_accurate[2]/total_values[2]) + " For optimizer:" + model_names[2])
confusion_dfs[2]
# -

# Conclusion:
# While training all three optimizer, we noticed that all three require different learning rates. While adam when supplied with larger rate, has hard time converging to a solution. Nesterov is able to work with higher learning rates. So the ideal learning rate for the current dataset which have worked are:
# Nesterov: 0.5, Adam:0.001, RMSProd: 0.001
#
# We see that with nesterov, we dont have to worry much about learning rate, as we can use learning rate scheduler, While the adam and rmsprop, use there own learning rates calculations. That is why when training with adam optimizer, we have to be mindfull of learning rates, while nesterov doesnt have those kind of restrictions. In case of RMSProp as well a lower learning rate yielded better results.
#
# Nesterov has also performed the best of all, even though by slight margin, So in all Nesterov would be preferred choice for current dataset. 
#
# Compared to RMSprop, Adam and nesterov optimizers also descent fast, and in few epochs itself they achieve good accuracies on validation set. This is usefull, as in case of large files, we can get away with less epochs.
#
# Overall all three have worked well, and depending on the problem, one could perform better than other. However for the current problem the nesterov seems like the best choice.



StarcoderdataJupyterScriptsDedupFiltered