#!/usr/bin/env python
# coding: utf-8

# In[1]:


# To support both python 2 and python 3
from __future__ import division, print_function, unicode_literals

# Common imports
import numpy as np
import os

# to make this notebook's output stable across runs
np.random.seed(42)

# To plot pretty figures
get_ipython().run_line_magic('matplotlib', 'inline')
import matplotlib as mpl
import matplotlib.pyplot as plt
mpl.rc('axes', labelsize=14)
mpl.rc('xtick', labelsize=12)
mpl.rc('ytick', labelsize=12)

# Where to save the figures
PROJECT_ROOT_DIR = "."
CHAPTER_ID = "classification"

def save_fig(fig_id, tight_layout=True):
    path = os.path.join(PROJECT_ROOT_DIR, "images", CHAPTER_ID, fig_id + ".png")
    print("Saving figure", fig_id)
    if tight_layout:
        plt.tight_layout()
    plt.savefig(path, format='png', dpi=300)

# Saving and loading sklearn objects to pickle files in hard disk
import os.path, joblib

def save_pkl(obj, name):
    joblib.dump(obj, name)
    print(name + 'saved!')

def load_pkl(name):
    name = name + '.pkl'
    if os.path.exists(name):
        return joblib.load(name)
    print('there is no' + name)
    return False


# In[3]:


# Some crazy stuff about ordering the data set so I get the same results as in the book
def sort_by_target(mnist):
    reorder_train = np.array(sorted([(target, i) for i, target in enumerate(mnist.target[:60000])]))[:, 1]
    reorder_test = np.array(sorted([(target, i) for i, target in enumerate(mnist.target[60000:])]))[:, 1]
    mnist.data[:60000] = mnist.data[reorder_train]
    mnist.target[:60000] = mnist.target[reorder_train]
    mnist.data[60000:] = mnist.data[reorder_test + 60000]
    mnist.target[60000:] = mnist.target[reorder_test + 60000]


# In[4]:


# Loading MINST dataset
mnist = load_pkl('mnist')
#Fetching MNISTdataset
if not mnist:
    from sklearn.datasets import fetch_openml
    mnist = fetch_openml('mnist_784', version=1, cache=True)
    save_pkl(mnist,'mnist')


# In[5]:


mnist.target = mnist.target.astype(np.int8) # fetch_openml() returns targets as strings
sort_by_target(mnist) # fetch_openml() returns an unsorted dataset


# In[6]:


X, y = mnist['data'], mnist['target']
print(X.shape,y.shape)


# In[7]:


#grab an instance’s feature vector, reshape it to a 28×28 array, and display it
some_digit = X[36000]
some_digit_image = some_digit.reshape(28,28)


# In[8]:


plt.imshow(some_digit_image)


# In[9]:


# The label corresponding to some_digit is indeed a 5
y[36000]


# In[10]:


#training and testing set
X_train, X_test, y_train, y_test = X[:60000], X[60000:], y[:60000], y[60000:]
#shuffling to guarantee cross validation folds will be similar
shuffle_index = np.random.permutation(60000)
#generates an array of indexes permutated
X_train, y_train= X_train[shuffle_index], y_train[shuffle_index]


# # Training a binary classifier: 5-detector

# In[11]:


#creating labels for this classification task
y_train_5 = (y_train == 5)
y_test_5 = (y_test == 5)


# In[15]:


sgd_clf = load_pkl('sgd_clf')

if not sgd_clf:
    from sklearn.linear_model import SGDClassifier
    sgd_clf = SGDClassifier(random_state=42)
    sgd_clf.fit(X_train, y_train_5)
    # Output a pickle file for the model
    save_pkl(sgd_clf,'sgd_clf')


# In[16]:


# Let us verify it works with the digit choosed before
sgd_clf.predict([some_digit]) == 5


# # Model performance

# ## Measuring Accuracy Using Cross-Validation

# In[17]:


from sklearn.model_selection import cross_val_score


# In[18]:


sgd_clf_accuracy_cv = cross_val_score(sgd_clf, X_train, y_train_5, cv=3, scoring="accuracy")


# In[39]:


sgd_clf_accuracy_cv


# In[19]:


np.savetxt('sgd_clf_accuracy_cv.np',sgd_clf_accuracy_cv, fmt='%f')


# Seems impressive that our SGD model has an accuracy over 96%, but with the next example we must check that this result is not to be taken seriously since it has some flaws related to the very same measure of acurracy and science case on a binary classifier for the MNIST data constitution :O

# In[20]:


from sklearn.base import BaseEstimator
# with BaseEstimator we create a classifier from scratch that says any feature from MNIST is not a 5

class Never5Classifier(BaseEstimator):
    
    # All hard work usually must be implemented on the fit method, but since the classifier only has to
    # say any X is labeled with a 5, it is left empty
    
    def fit(self, X, y=None):
        pass
    # predict method just outputs our y vector which is pretty dumb in this case: True, all the time.
    def predict(self, X):
        return np.zeros((len(X), 1), dtype=bool)


# In[21]:


# Now it is time to implement our Mighty Binary classifier XD and check its accuracy
never_5_clf = Never5Classifier()
cross_val_score(never_5_clf, X_train, y_train_5, cv=3, scoring="accuracy")


# This result shows why accuracy is not the ultimate word on how well the ML model performs and it is related to the very nature of the task and the data:
# 
# 1. It is a binary classifier
# 2. The statistical distribution of the data is such that 5 is present 10% of the time
# 
# Regarding item (1), if it were a multi classifier then the accuracy of the Mighty-non-5 classifier wouldn't be over 90%. Regarding item (2), the data has the asymmetry of having 1/10 of its Xs being 5s; if the asymmetry were 1/100 5s, then the accuracy will be higher or if it were 5s all the time, then the accuracy will be 0%.
# 
# These two points are quite important and must be considered when training a model, and a reason why another performance measures must be considered.

# ## Confusion Matrix

# From the book:
# "
# The general idea is to count the number of times instances of class A are classified as class B. To compute the confusion matrix, you first need to have a set of predictions, so they can be compared to the actual targets.
# "

# In[24]:


if os.path.exists('y_train_pred.np'):
    y_train_pred = np.loadtxt('y_train_pred.np', dtype=float)
else:
    from sklearn.model_selection import cross_val_predict
    y_train_pred = cross_val_predict(sgd_clf, X_train, y_train_5, cv=3)
    np.savetxt('y_train_pred.np',y_train_pred, fmt='%f')


# In[25]:


from sklearn.metrics import confusion_matrix
CM = confusion_matrix(y_train_5, y_train_pred)


# In[26]:


print(CM)


# "Each row in a confusion matrix represents an actual class, while each column represents a predicted class."
# 
# In this case, the first row represents the negative class, ie, the true negatives; meaning that CM[0,0] is the number of non-5s classified correctly while CM[0,1] represents the number false positives, ie, the non-5s classified as 5s. The second class represents the possitive class. CM[1,0] is the number of false negatives, ie, the number of 5s classified as non-5s while CM[1,1] represents the number of true positives, ie, the number of 5s classified as 5s.
# 
# "A perfect classifier would have only true positives and true negatives"

# ## Metrics derived from the confusion matrix

# ### Accuracy of the positive predictions: the precision of the classifier

# $$precision = \frac{TP}{TP+FP}$$

# ### Recall: ratio of positive instances that are correctly detected by the classifier

# $$recall = \frac{TP}{TP+FN}$$

# In[27]:


from sklearn.metrics import precision_score, recall_score
precision = precision_score(y_train_5, y_train_pred)
recall = recall_score(y_train_5, y_train_pred)


# In[28]:


print(precision,recall)


# Here we go nuts :O since this means that our classifier is correct only 88% of the time when it says an image is a 5 (precision). On the other hand, it only detects 68% of the 5s (recall).

# ### F$_1$ score: the harmonic mean

# $$F_1 = \frac{TP}{TP+\frac{FN+FP}{2}}= 2 \frac{precision \times recall}{precision + recall}$$

# In[29]:


from sklearn.metrics import f1_score
F_1 = f1_score(y_train_5, y_train_pred)


# In[30]:


print(F_1)


# https://towardsdatascience.com/beyond-accuracy-precision-and-recall-3da06bea9f6c
# 
# "We use the harmonic mean instead of a simple average because it punishes extreme values. A classifier with a precision of 1.0 and a recall of 0.0 has a simple average of 0.5 but an F1 score of 0. The F1 score gives equal weight to both measures and is a specific example of the general Fβ metric where β can be adjusted to give more weight to either recall or precision... If we want to create a balanced classification model with the optimal balance of recall and precision, then we try to maximize the F1 score."

# ## Precision/Recall Tradeoff

# It is clear that there is correlation between precision and recall and that according to our science case, we should chose the combination that gives us the best results for our purposes. With that in mind, lets, examine the the precision and recal from our SGDClassifier as a function the threshold it uses to decide when classifying (an hyperparameter of the model).

# In[38]:


# Let's obtain the decition scares our classifier uses
from sklearn.model_selection import cross_val_predict

if os.path.exists('y_scores.np'):
    y_scores = np.loadtxt('y_scores.np', dtype=float)
else:
    y_scores = cross_val_predict(sgd_clf, X_train, y_train_5, cv=3, method= "decision_function")
    np.savetxt('y_scores.np',y_scores, fmt='%f')


# In[40]:


# Now with the decision scores, we proceed to compute teh precision and recall curve
from sklearn.metrics import precision_recall_curve

precisions, recalls, thresholds = precision_recall_curve(y_train_5, y_scores)


# In[41]:


def plot_precision_recall_vs_threshold(precisions, recalls, thresholds):
    plt.plot(thresholds, precisions[:-1], "b--", label="Precision")
    plt.plot(thresholds, recalls[:-1], "g-", label="Recall")
    plt.xlabel("Threshold")
    plt.legend(loc="upper left")
    plt.ylim([0, 1])
    plt.xlim([-25000, 25000])


# In[42]:


plot_precision_recall_vs_threshold(precisions, recalls, thresholds)
plt.show()


# ### The ROC Curve:  receiver operating characteristic (ROC) curve

# Fom the book:
# "the ROC curve plots the true positive rate (another name for recall) against the false positive rate. The FPR is the ratio of negative instances that are incorrectly classified as positive. It is equal to one minus the true negative rate, which is the ratio of negative instances that are correctly classified as negative. The TNR is also called specificity. Hence the ROC curve plots sensitivity (recall) versus 1 – specificity."

# In[43]:


from sklearn.metrics import roc_curve

fpr, tpr, thresholds = roc_curve(y_train_5, y_scores)


# In[44]:


def plot_roc_curve(fpr, tpr, label=None):
    plt.plot(fpr, tpr, linewidth=2, label=label)
    plt.plot([0, 1], [0, 1], 'k--')
    plt.axis([0, 1, 0, 1])
    plt.xlabel('False Positive Rate')
    plt.ylabel('True Positive Rate')


# In[45]:


plot_roc_curve(fpr, tpr)
plt.show()


# "Once again there is a tradeoff: the higher the recall (TPR), the more false positives (FPR) the classifier produces. The dotted line represents the ROC curve of a purely random classifier; a good classifier stays as far away from that line as possible (toward the top-left corner)."

# ### Multiclass Classification

# There are ML algorithms that are able to handle multiclass classification, but ohers than are able only to do binary classification. Nonetheless binary classifier can be used to generate multiclass classifiers using the one versus all strategy, where as many bianry classifiers as classes to classify we have, then we compare the classification scores and chose the one with the higest score. There are different alternatives too, depending on the dataset and on the problem to handle...

# "Scikit-Learn detects when you try to use a binary classification algorithm for a multiclass classification task, and it automatically runs OvA"

# In[46]:


multi_sgd_clf = load_pkl('multi_sgd_clf')
if not multi_sgd_clf:
    multi_sgd_clf.fit(X_train, y_train)
    save_pkl(multi_sgd_clf,'multi_sgd_clf')


# In[52]:


multi_sgd_clf.predict([some_digit])


# In[53]:


some_digit_scores = multi_sgd_clf.decision_function([some_digit])


# In[54]:


np.argmax(some_digit_scores)


# In[55]:


multi_sgd_clf.classes_[5]


# "When a classifier is trained, it stores the list of target classes in its classes_ attribute, ordered by value. In this case, the index of each class in the classes_ array conveniently matches the class itself (e.g., the class at index 5 happens to be class 5), but in general you won’t be so lucky."

# Now let's evaluate the classifier using cross validation:

# In[56]:


multi_cv_scores = cross_val_score(multi_sgd_clf, X_train, y_train, cv=3, scoring="accuracy")


# In[58]:


np.savetxt('multi_cv_scores.np',multi_cv_scores, fmt='%f')


# In[59]:


multi_cv_scores


# This is a good performance, meaning the data is well behaved on how its classes are distributed

# ## Error Analysis

# The goal here is to tweak the model for improving it by exploring the different erros it makes. To do so, let's build its confussion matrix but scaliting the imput parameters first, because this enhances our classifier.

# In[65]:


# sacating input vectors
from sklearn.preprocessing import StandardScaler
if os.path.exists('X_train_scaled.np'):
    X_train_scaled = np.loadtxt('X_train_scaled.np', dtype=float)
else:
    scaler = StandardScaler()
    X_train_scaled = scaler.fit_transform(X_train.astype(np.float64))
    np.savetxt('X_train_scaled.np',X_train_scaled, fmt='%f')


# In[67]:


# Training the model on the scaled input training vector
if os.path.exists('multi_y_train_pred.np'):
    multi_y_train_pred = np.loadtxt('multi_y_train_pred.np', dtype=float)
else:
    multi_y_train_pred = cross_val_predict(multi_sgd_clf, X_train_scaled, y_train, cv=3)
    np.savetxt('multi_y_train_pred.np',multi_y_train_pred, fmt='%f')


# In[68]:


# Confusion matrix
multi_CM = confusion_matrix(y_train, multi_y_train_pred)


# plt.matshow(multi_CM, cmap=plt.cm.gray)
# plt.show()

# In[71]:


# Computing error rates
row_sums = multi_CM.sum(axis=1, keepdims=True)
norm_multi_CM = multi_CM / row_sums
#diagonal with zeros to keep only the errors
np.fill_diagonal(norm_multi_CM, 0)
plt.matshow(norm_multi_CM, cmap=plt.cm.gray)
plt.show()


# In[ ]:





# ## Multilabel classification

# In[ ]:




