import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from sklearn import *
import scipy

%matplotlib inline

"""
Explanation: MNIST
End of explanation
"""


training = pd.read_csv("/data/MNIST/mnist_train.csv", header = None)
testing = pd.read_csv("/data/MNIST/mnist_test.csv", header = None)

X_train, y_train = training.iloc[:, 1:].values, training.iloc[:, 0].values 
X_test, y_test = testing.iloc[:, 1:].values, testing.iloc[:, 0].values 

print("Shape of X_train: ", X_train.shape, "shape of y_train: ", y_train.shape)
print("Shape of X_test: ", X_test.shape, "shape of y_test: ", y_test.shape)

"""
Explanation: Dataset description
Datasource: http://yann.lecun.com/exdb/mnist/
The training dataset consists of 60,000 training digits and the test set contains 10,000 samples, respectively. The images in the MNIST dataset consist of  pixels, and each pixel is represented by a gray scale intensity value. Here, we unroll the  pixels into 1D row vectors, which represent the rows in our image array (784 per row or image). The second array (labels) returned by the load_mnist function contains the corresponding target variable, the class labels (integers 0-9) of the handwritten digits.
Csv version of the files are available in the following links.
CSV training set http://www.pjreddie.com/media/files/mnist_train.csv
CSV test set http://www.pjreddie.com/media/files/mnist_test.csv
End of explanation
"""


label_counts = pd.DataFrame({
    "train": pd.Series(y_train).value_counts().sort_index(), 
    "test": pd.Series(y_test).value_counts().sort_index()
})
(label_counts / label_counts.sum()).plot.bar()

plt.xlabel("Class")
plt.ylabel("Frequency (normed)")

"""
Explanation: Distribution of class frequencies
End of explanation
"""


scipy.stats.chisquare(label_counts.train, label_counts.test)

"""
Explanation: Chisquare test on class frequencies
End of explanation
"""


fig, axes = plt.subplots(5, 5, figsize = (15, 9))
for i, ax in enumerate(fig.axes):
    img = X_train[i, :].reshape(28, 28)
    ax.imshow(img, cmap = "Greys", interpolation="nearest")
    ax.set_title("True: %i" % y_train[i])
plt.tight_layout()

"""
Explanation: Display a few sample images
End of explanation
"""


fig, axes = plt.subplots(10, 5, figsize = (15, 20))

for i, ax in enumerate(fig.axes):
    img = X_train[y_train == 7][i, :].reshape(28, 28)
    ax.imshow(img, cmap = "Greys", interpolation="nearest")

plt.tight_layout()

"""
Explanation: View different variations of a digit
End of explanation
"""


scaler = preprocessing.StandardScaler()
X_train_std = scaler.fit_transform(X_train.astype(np.float64))
X_test_std = scaler.transform(X_test.astype(np.float64))

"""
Explanation: Feature scaling
End of explanation
"""


%%time 
lr = linear_model.LogisticRegression()
lr.fit(X_train_std, y_train)
print("accuracy:", lr.score(X_test_std, y_test))

"""
Explanation: Applying logistic regression classifier
End of explanation
"""


y_test_pred = lr.predict(X_test_std)
miss_indices = (y_test != y_test_pred)
misses = X_test[miss_indices]
print("No of miss: ", misses.shape[0])

fig, axes = plt.subplots(10, 5, figsize = (15, 20))
misses_actual = y_test[miss_indices]
misses_pred = y_test_pred[miss_indices]

for i, ax in enumerate(fig.axes):
    img = misses[i].reshape(28, 28)
    ax.imshow(img, cmap = "Greys", interpolation="nearest")
    ax.set_title("A: %s, P: %d" % (misses_actual[i], misses_pred[i]))
plt.tight_layout()

"""
Explanation: Display wrong predictions
End of explanation
"""


inits = np.random.randn(10, 784) 
inits = inits / np.std(inits, axis=1).reshape(10, -1)

%%time
est = linear_model.SGDClassifier(n_jobs=4, tol=1e-5, eta0 = 0.15,  
                                 learning_rate = "invscaling", 
                                 alpha = 0.01, max_iter= 100)
est.fit(X_train_std, y_train, inits)
print("accuracy", est.score(X_test_std, y_test), "iterations:", est.n_iter_)

fig, _ = plt.subplots(3, 4, figsize = (15, 10))
for i, ax in enumerate(fig.axes):
    if i < est.coef_.shape[0]:
        ax.imshow(est.coef_[i, :].reshape(28, 28), cmap = "bwr", interpolation="nearest")
    else:
        ax.remove()

pd.DataFrame(est.coef_[0, :].reshape(28, 28))

"""
Explanation: Applying SGD classifier
End of explanation
"""


