{"cells":[{"metadata":{"trusted":true},"cell_type":"code","source":"import numpy as np\nnp.random.seed(1)\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Loading the data"},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"train_data = pd.read_csv('../input/digit-recognizer/train.csv')\nfig, ax = plt.subplots(4, 4, figsize=(8,8))\nfor i in range(4):\n    for j in range(4):\n        ax[i, j].imshow(train_data.loc[i*37 + j*23, 'pixel0':].to_numpy().reshape(28, 28), cmap = 'gray')\n        ax[i, j].set_title('label = %s' % (train_data.loc[i*37 + j*23, 'label']))\n        ax[i, j].set_xticks([])\n        ax[i, j].set_yticks([])\ntrain_data","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Training data consists of 42000 training examples.  \nEach training example is a 28x28 gray-scale images with a number from 0 to 9.  \nFor convenience each image was reshaped from 28x28 matrix into 28^2 = 784-element row.  \nEach pixel has a single pixel-value associated with it,  \nindicating the lightness or darkness of that pixel, with higher numbers meaning darker.  \n  \n\nData should be normalized"},{"metadata":{"trusted":true},"cell_type":"code","source":"train_data.loc[:, 'pixel0':'pixel783'] = (train_data.loc[:, 'pixel0':'pixel783'] - train_data.loc[:, 'pixel0':'pixel783'].mean(axis = 0)) / 255\ntrain_data","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Building a network\n\nWe will build an L-layer network, final *number of layers L* will be defined according to the perfomance.  \n\nBuilding steps:  \n\n- Initialize parameters  \n- Forward propagation\n- Compute cost with L2 regularization\n- Backward propagation\n- Update parameters\n- Train model\n- Make a prediction with trained parameters"},{"metadata":{"trusted":true},"cell_type":"code","source":"def sigmoid(Z):\n    \"\"\"\n    Arguments:\n    Z - a scalar or numpy array of any size\n    \n    Return:\n    sigmoid(Z)\n    \"\"\"\n    \n    return 1/(1 + np.exp(-Z))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def softmax(Z):\n    \"\"\"\n    Arguments:\n    Z - a numpy array of any size\n    \n    Return:\n    softmax(Zi) = exp(Zi) / sum(exp(Zi))\n    \"\"\"\n    \n    A = np.exp(Z) / np.sum(np.exp(Z), axis = 0)\n    \n    return A","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def random_mini_batches(X, Y, mini_batch_size = 256):\n    \"\"\"\n    Create a list of random mini batches from (X, Y)\n    \n    Arguments:\n    X - dataset, size = (n, m)\n    Y - label vecor, size = (number of labels, m)\n    \n    Return:\n    mini_batches - list of synchronous (mini_batch_X, mini_batch_Y)\n    \"\"\"\n    \n    m = X.shape[1]\n    mini_batches = list()\n    \n    permutation = list(np.random.permutation(m)) # create a list of permutations\n    shuffled_X = X.iloc[:, permutation]\n    shuffled_Y = Y[:, permutation].reshape((Y.shape[0], m))\n    \n    # Partition\n    num_complete_minibatches = np.floor(m/mini_batch_size).astype(int)\n    for k in range(num_complete_minibatches):\n        \n        mini_batch_X = shuffled_X.iloc[:, k*mini_batch_size : (k+1)*mini_batch_size]\n        mini_batch_Y = shuffled_Y[:, k*mini_batch_size : (k+1)*mini_batch_size]\n        \n        mini_batch = (mini_batch_X, mini_batch_Y)\n        mini_batches.append(mini_batch)\n        \n    # Handling the end casse (last mini-batch < mini_batch_size)\n    if m % mini_batch_size != 0:\n        \n        mini_batch_X = shuffled_X.iloc[:, mini_batch_size*(m // mini_batch_size) : ]\n        mini_batch_Y = shuffled_Y[:, mini_batch_size*(m // mini_batch_size) : ]\n        \n        mini_batch = (mini_batch_X, mini_batch_Y)\n        mini_batches.append(mini_batch)\n        \n    return mini_batches","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def initialize_parameters(dims):\n    \"\"\"\n    Arguments:\n    dims - list of the dimensions of layers, the first element is number of features,\n        the last is a number of labels, type = list\n    \n    Return:\n    parameters - Wi - weights, bi - bias term for i-th layer, type = dict\n    parameters = {..., Wi : , bi : , ...}\n    \"\"\"\n    \n    L = len(dims) - 1 # The input layer doesn't count in layer number\n    parameters = dict()\n    \n    for l in range(L):\n        \n        var_l = 1 / np.sqrt(dims[l])\n        parameters['W'+str(l+1)] = np.random.randn(dims[l+1], dims[l]) * var_l\n        parameters['b'+str(l+1)] = np.zeros((dims[l+1], 1))\n        \n    return parameters","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def forward_propagation(X, parameters):\n    \"\"\"\n    Arguments:\n    X - dataset\n    parameters - dictionary with W - weights, b - bias term\n    L - number of layers\n    \n    Return:\n    forward_cache - dictionary with Zi and Ai values\n    \"\"\"\n    L = len(parameters) // 2\n    forward_cache = {'A0' : X, 'Z0' : X}\n    \n    for l in range(L - 1):\n        \n        forward_cache['Z'+str(l+1)] = np.dot(parameters['W'+str(l+1)], forward_cache['A'+str(l)]) + parameters['b'+str(l+1)]\n        forward_cache['A'+str(l+1)] = sigmoid(forward_cache['Z'+str(l+1)])\n    \n    forward_cache['Z' + str(L)] = np.dot(parameters['W'+str(L)], forward_cache['A'+str(L-1)]) + parameters['b'+str(L)]\n    forward_cache['A' + str(L)] = softmax(forward_cache['Z' + str(L)])\n    \n    return forward_cache","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def compute_cost(A, Y, parameters, lambd):\n    \"\"\"\n    Arguments:\n    A - an output layer value\n    Y - true label\n    parameters - weights and bias terms\n    lambd - regularization hyperparameter\n    \n    Return:\n    cost - cross-entropy cost\n    \"\"\"\n    reg = 0\n    \n    for l in range(len(parameters)//2):\n        reg += np.sum(parameters['W'+str(l+1)]**2)\n    \n    cost = -1 * np.sum(Y*np.log(A)) + lambd/2 * reg\n    \n    return cost","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def backward_propagation(Y, forward_cache, parameters, lambd):\n    \"\"\"\n    Arguments:\n    forward_cache - dictionary with Zi and Ai\n    parameters - weights and bias terms\n    lambd - regularization hyperparameter\n    \n    Return:\n    grads - dictionary with dWi and dbi\n        dWi - derivative of loss with the respect to Wi\n        dbi - derivative of loss with the respect to bi\n    \"\"\"\n    \n    m = Y.shape[1]\n    \n    L = len(parameters) // 2\n    grads = {'dZ'+str(L) : forward_cache['A'+str(L)] - Y}\n    dZL = grads['dZ' + str(L)]\n    grads['dW'+str(L)] = 1/m * np.dot(dZL, forward_cache['A' + str(L-1)].T) + lambd/m * parameters['W'+str(L)]\n    grads['db'+str(L)] = 1/m * np.sum(dZL, axis = 1, keepdims = True)\n    \n    for l in range(L - 1, 0, -1):\n        \n        grads['dZ'+str(l)] = 1/m * np.dot(parameters['W' + str(l+1)].T, grads['dZ'+str(l+1)])*forward_cache['A'+str(l)]*(1-forward_cache['A'+str(l)])\n        grads['dW'+str(l)] = 1/m * np.dot(grads['dZ'+str(l)], forward_cache['A'+str(l-1)].T)+lambd/m*parameters['W'+str(l)]\n        grads['db'+str(l)] = 1/m * np.sum(grads['dZ'+str(l)], axis = 1, keepdims = True)\n        \n    return grads","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def initialize_adam(parameters):\n    \"\"\"\n    Arguments:\n    parameters - dictionary with W - weights, b - bias terms\n    \n    Return:\n    v -- python dictionary that will contain the exponentially weighted average of the gradient.\n                    v[\"dW\" + str(l)] = ...\n                    v[\"db\" + str(l)] = ...\n    s -- python dictionary that will contain the exponentially weighted average of the squared gradient.\n                    s[\"dW\" + str(l)] = ...\n                    s[\"db\" + str(l)] = ...\n    \"\"\"\n    \n    L = len(parameters) // 2\n    v = dict()\n    s = dict()\n    \n    \n    for l in range(L):\n        \n        v[\"dW\" + str(l+1)] = np.zeros((parameters['W' + str(l+1)].shape[0], parameters['W' + str(l+1)].shape[1]))\n        v[\"db\" + str(l+1)] = np.zeros((parameters['b' + str(l+1)].shape[0], parameters['b' + str(l+1)].shape[1]))\n        s[\"dW\" + str(l+1)] = np.zeros((parameters['W' + str(l+1)].shape[0], parameters['W' + str(l+1)].shape[1]))\n        s[\"db\" + str(l+1)] = np.zeros((parameters['b' + str(l+1)].shape[0], parameters['b' + str(l+1)].shape[1]))\n    \n    return v, s","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def update_parameters(parameters, grads, learning_rate):\n    \"\"\"\n    Arguments:\n    parameters - weigths and bias terms, dictionary\n    grads - dictionary with derivatives dW and db\n    learning_rate - learning rate hyperparameter, scalar\n    \n    Return:\n    updated parameters\n    \"\"\"\n    \n    L = len(parameters) // 2\n    \n    for l in range(L):\n        \n        parameters['W'+str(l+1)] = parameters['W'+str(l+1)] - learning_rate * grads['dW'+str(l+1)]\n        parameters['b'+str(l+1)] = parameters['b'+str(l+1)] - learning_rate * grads['db'+str(l+1)]\n        \n    return parameters","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def update_parameters_with_adam(parameters, grads, v, s, t, learning_rate, beta1 = 0.9, beta2 = 0.999,  epsilon = 10**(-8)):\n    \"\"\"\n    Arguments:\n    parameters - python dictionary containing your parameters:\n                    parameters['W' + str(l)] = Wl\n                    parameters['b' + str(l)] = bl\n    grads - python dictionary containing your gradients for each parameters:\n                    grads['dW' + str(l)] = dWl\n                    grads['db' + str(l)] = dbl\n    v - Adam variable, moving average of the first gradient, python dictionary\n    s - Adam variable, moving average of the squared gradient, python dictionary\n    learning_rate -- the learning rate, scalar.\n    beta1 -- Exponential decay hyperparameter for the first moment estimates \n    beta2 -- Exponential decay hyperparameter for the second moment estimates \n    epsilon -- hyperparameter preventing division by zero in Adam updates\n\n    Returns:\n    parameters - python dictionary containing your updated parameters \n    v - Adam variable, moving average of the first gradient, python dictionary\n    s - Adam variable, moving average of the squared gradient, python dictionary\n    \"\"\"\n    \n    L = len(parameters) // 2\n    v_corrected = dict()\n    s_corrected = dict()\n    \n    for l in range(L):\n        \n        v[\"dW\" + str(l+1)] = beta1 * v[\"dW\" + str(l+1)] + (1 - beta1) * grads['dW' + str(l+1)]\n        v[\"db\" + str(l+1)] = beta1 * v[\"db\" + str(l+1)] + (1 - beta1) * grads['db' + str(l+1)]\n        \n        v_corrected[\"dW\" + str(l+1)] = v['dW' + str(l+1)] / (1 - beta1**t)\n        v_corrected[\"db\" + str(l+1)] = v['db' + str(l+1)] / (1 - beta1**t)\n        \n        s[\"dW\" + str(l+1)] = beta2 * s[\"dW\" + str(l+1)] + (1 - beta2) * grads['dW' + str(l+1)]**2\n        s[\"db\" + str(l+1)] = beta2 * s[\"db\" + str(l+1)] + (1 - beta2) * grads['db' + str(l+1)]**2\n\n        s_corrected[\"dW\" + str(l+1)] = s['dW' + str(l+1)] / (1 - beta2**t)\n        s_corrected[\"db\" + str(l+1)] = s['db' + str(l+1)] / (1 - beta2**t)\n        \n        parameters[\"W\" + str(l+1)] = parameters['W' + str(l+1)] - learning_rate * v_corrected['dW' + str(l+1)] / (np.sqrt(s_corrected['dW' + str(l+1)]) + epsilon) \n        parameters[\"b\" + str(l+1)] = parameters['b' + str(l+1)] - learning_rate * v_corrected['db' + str(l+1)] / (np.sqrt(s_corrected['db' + str(l+1)]) + epsilon)\n\n    return parameters, v, s","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def model(X, Y, X_val, Y_val, dims, algos, learning_rate = 0.001, num_epochs = 15000, mini_batch_size = 0, lambd = 0, beta1 = 0.9, beta2 = 0.999, epsilon = 1e-8, print_cost = False):\n    \"\"\"\n    Arguments:\n    X - dataset\n    Y - true label vector\n    dims - list of the dimensions of layers\n    learning_rate -  learning rate hyperparameter, default -- 0.001\n    num_iterations - number of iterations for gradient descent, default -- 15000\n    lambd - regularization hyperparameter, default -- 0\n    print_cost - bool, if True - print cost each 1000 iterations, default -- False\n    \n    Return:\n    parameters - learned parameters, dict\n    \"\"\"\n    \n    parameters = dict()\n    forward_cache = dict()\n    grads = dict()\n    costs = list() # list of costs after each # of epochs for training set\n    costs_val = list() # list of costs after # of epochs for cross-validation set\n    m = X.shape[1]\n    m_val = X_val.shape[1]\n    \n    parameters = initialize_parameters(dims)\n    v, s = initialize_adam(parameters)\n    L = len(parameters) // 2\n    \n    for i in range(num_epochs):\n        \n        if mini_batch_size == 0:\n            forward_cache = forward_propagation(X, parameters)\n            forward_cache_val = forward_propagation(X_val, parameters)\n        \n            if i % 100 == 0:\n                cost = compute_cost(forward_cache['A'+str(L)], Y, parameters, lambd)\n                costs.append(cost/m)\n                cost = compute_cost(forward_cache_val['A'+str(L)], Y_val, parameters, lambd)\n                costs_val.append(cost/m_val)\n                if i % 1000 == 0 and print_cost is True:\n                    print('Train cost after {} epochs: {}'.format(i, cost))\n        \n            grads = backward_propagation(Y, forward_cache, parameters, lambd)\n        \n            if algos == 'grad':\n                parameters = update_parameters(parameters, grads, learning_rate)\n            elif algos == 'adam':\n                parameters, v, s = update_parameters_with_adam(parameters, grads, v, s, i+1, learning_rate,\n                                beta1, beta2,  epsilon)\n                \n        elif mini_batch_size != 0:\n            cost_total = 0\n            minibatches = random_mini_batches(X, Y, mini_batch_size)\n            \n            for minibatch in minibatches:\n                \n                (minibatch_X, minibatch_Y) = minibatch\n                \n                forward_cache = forward_propagation(minibatch_X, parameters)\n                cost_total += compute_cost(forward_cache['A'+str(L)], minibatch_Y, parameters, lambd)\n                \n                grads = backward_propagation(minibatch_Y, forward_cache, parameters, lambd)\n                \n                if algos == 'grad':\n                    parameters = update_parameters(parameters, grads, learning_rate)\n                elif algos == 'adam':\n                    parameters, v, s = update_parameters_with_adam(parameters, grads, v, s, i+1, learning_rate,\n                                beta1, beta2,  epsilon)\n            \n            cost_avg = cost_total / m\n            \n            if i % 10 == 0:\n                cost = cost_avg\n                forward_cache_val = forward_propagation(X_val, parameters)\n                costs_val.append(compute_cost(forward_cache_val['A'+str(L)], Y_val, parameters, lambd)/m_val)\n                costs.append(cost)        \n                if i % 25 == 0 and print_cost is True:\n                    print('Train cost after {} epochs: {}'.format(i, cost))\n        \n    ## Plot learning curve:\n    plt.plot(costs, label = 'Training set')\n    plt.plot(costs_val, label = 'Dev set')\n    plt.xlabel('Iterations')\n    plt.ylabel('Cost')\n    plt.title('Cost with learning rate = {}, lambda = {}, dims = {}'.format(learning_rate, lambd, dims))\n    plt.legend()\n    plt.show()\n    ##    \n        \n    return parameters","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def predict(X, Y, parameters, comp_accur = True):\n    \"\"\"\n    Arguments:\n    X - dataset\n    Y - true label vector\n    parameters - weights and bias terms\n    \n    Return:\n    Y_pred - predicted label vector\n    accur - accuracy on predicted\n    \"\"\"\n    \n    L = len(parameters) // 2 \n    \n    forward_cache = forward_propagation(X, parameters)\n    A = forward_cache['A'+str(L)]\n    \n    Y_pred = np.argmax(A, axis = 0)\n    Y = np.argmax(Y, axis = 0)\n    if comp_accur is False:\n        return Y_pred\n    accur = np.mean(Y == Y_pred) * 100\n    \n    return Y_pred, accur","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We need to prepare data as follows:  \n- Label vector for training should have the shape (n, m), where m is a number of training examples, and n is a number of labels. n equals to 10 in this case (10 numbers from 0 to 9).\n- Each value of label vector will be an index where 1 will be put, and in other values will be zeros.  \n\nFor example:"},{"metadata":{"trusted":true},"cell_type":"code","source":"Y_examp_raw = np.array([7, 3, 4]) # unprepared vector\nY_examp = np.zeros((10, 3)) # n = 10, m = 3\nfor i in range(3):\n    for j in range(10):\n        Y_examp[j, i] = 1 if Y_examp_raw[i] == j else 0\n\nprint(Y_examp)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Let's apply this on our training set"},{"metadata":{"trusted":true},"cell_type":"code","source":"Y_raw = train_data.loc[:, 'label'].T\nY = np.zeros((10, Y_raw.shape[0]))\n\nfor i in range(Y.shape[1]):\n    for j in range(Y.shape[0]):\n        Y[j, i] = 1 if Y_raw[i] == j else 0\n    \nprint(Y.shape)\nprint(Y)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"X = train_data.loc[:, 'pixel0':].T\n# Training set\nX_train = X.loc[:, :37999]\nY_train = Y[:, :38000] \n# Cross-validation set\nX_val = X.loc[:, 38000:]\nY_val = Y[:, 38000:]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Let's train\n\nWe will try different combinations of hyperparameters and choose the best according to the perfomance"},{"metadata":{"trusted":true},"cell_type":"code","source":"dims = [X_train.shape[0], 100, Y_train.shape[0]]\nparameters = model(X_train, Y_train, X_val, Y_val, dims = dims, algos = 'adam', num_epochs = 500, print_cost = True, lambd = 0.001, mini_batch_size = 128)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Prediction  "},{"metadata":{"trusted":true},"cell_type":"code","source":"Y_pred_train, accur_train = predict(X_train, Y_train, parameters)\nY_pred_val, accur_val = predict(X_val, Y_val, parameters)\nprint('Prediction accuracy on training set: {} %'.format(accur_train))\nprint('Prediction accuracy on cross-validation set: {} %'.format(accur_val))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# Return training and cross-validation label vectors to original shape\nY_train = np.argmax(Y_train, axis = 0)\nY_val = np.argmax(Y_val, axis = 0)\nprint(Y_pred_train.shape, Y_train.shape, Y_pred_val.shape, Y_val.shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"df_train = pd.DataFrame(data = {'Y_Actual_train':Y_train, 'Y_Pred_train':Y_pred_train})\ndf_val = pd.DataFrame(data = {'Y_Actual_val':Y_val, 'Y_Pred_val':Y_pred_val})\nconf_matrix_train = pd.crosstab(df_train['Y_Actual_train'], df_train['Y_Pred_train'], rownames = ['Actual'], colnames = ['Predicted'])\nconf_matrix_val = pd.crosstab(df_val['Y_Actual_val'], df_val['Y_Pred_val'], rownames = ['Actual'], colnames = ['Predicted'])\n# Confusion matrix for training set\nplt.figure(figsize = (10, 8))\nsns.heatmap(data = conf_matrix_train, annot = True)\nplt.title('Training set confusion matrix')\nplt.show()\n# Confusion matrix for training set\nplt.figure(figsize = (10, 8))\nsns.heatmap(data = conf_matrix_val, annot = True)\nplt.title('Cross-validation set confusion matrix')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Make submission\n\nMaking submission for Kaggle"},{"metadata":{"trusted":true},"cell_type":"code","source":"test_data = pd.read_csv('../input/digit-recognizer/test.csv')\ntest_data = (test_data - test_data.mean(axis = 0)) / 255\nX_test = test_data.T\nY_subm = predict(X_test, Y, parameters, comp_accur = False)\nsubmission = pd.DataFrame(data = {'ImageId':list(i for i in range(1, len(Y_subm) + 1)), 'Label':Y_subm.T}).set_index('ImageId')\nsubmission.to_csv('./submission.csv')","execution_count":null,"outputs":[]}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}