{ "cells": [ { "cell_type": "code", "execution_count": null, "metadata": { "id": "DA42OJmJ1PUv" }, "outputs": [], "source": [ "import pandas as pd\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import sys" ] }, { "cell_type": "code", "source": [ "from google.colab import drive\n", "drive.mount('/content/drive')" ], "metadata": { "id": "c9PIpG12xXEw", "outputId": "542b36b1-9566-4746-a52f-364177efe649", "colab": { "base_uri": "https://localhost:8080/" } }, "execution_count": null, "outputs": [ { "output_type": "stream", "name": "stdout", "text": [ "Mounted at /content/drive\n" ] } ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "_xMpZeS71PUz" }, "outputs": [], "source": [ "\n", "def sigmoid(z):\n", " return 1.0/(1 + np.exp(-z))\n", "\n", "def onehot(y, n_classes):\n", " # Encode labels into one-hot representation\n", " # y: labels [n_samples]\n", " # n_classes: number of classes\n", " # return onehot: [n_samples, n_classlabels]\n", " n_samples = y.shape[0]\n", " y_onehot = np.zeros((n_samples, n_classes))\n", " for i in range(n_samples):\n", " y_onehot[i, y[i]]=1\n", "\n", "\n", "\n", " return y_onehot\n", "\n", "def forward(X, w_h, b_h, w_out, b_out):\n", " # step 1: net input of hidden layer, z_h [n_samples, n_hidden]\n", " # step 2: activation of hidden layer, a_h [n_samples, n_hidden]\n", " # step 3: net input of output layer, z_out [n_samples, n_classlabels]\n", " # step 4: activation output layer, a_out [n_samples, n_classlabels]\n", "\n", " z_h = np.matmul(X, w_h) + b_h\n", " a_h = sigmoid(z_h)\n", " z_out = np.matmul(a_h, w_out) + b_out\n", " a_out = sigmoid(z_out)\n", " return z_h, a_h, z_out, a_out\n", "\n", "def compute_cost(y_onehot, output, l2, w_h, w_out):\n", " # Compute cost function\n", " # y_onehot: one-hot encoded output [n_samples, n_classlabels]\n", " # output: activation output, a_out [n_samples, n_classlabels]\n", " # l2: regularization coefficient lambda\n", " # w_h, w_out: current weight for hidden and output layers\n", " # return cost: cost value with l2 reguralization\n", " k = len(y_onehot)\n", " m = len(output)\n", " ylog1 = np.matmul(y_onehot, np.log(output).T)\n", " ylog2 = np.matmul((1-y_onehot), np.log(1 - output).T)\n", " hsquare = w_h**2\n", " outsquare = w_out ** 2\n", "\n", " reg = (l2 / (2 * k)) * (np.sum(hsquare) + np.sum(outsquare))\n", " cost = -(np.sum(np.sum(ylog1 + ylog2))) + (l2) * (np.sum(hsquare) + np.sum(outsquare))\n", "\n", " return cost\n", "\n", "def predict(X, w_h, b_h, w_out, b_out):\n", " # predict class labels\n", " # X: input array [n_samples, n_features]\n", " # y_pred: prediced class labels [n_samples]\n", "\n", " _, _, _, a_out = forward(X, w_h, b_h, w_out, b_out)\n", " y_pred = np.argmax(a_out, axis=1)\n", "\n", " return y_pred\n", "\n", "def train(X_train, Y_train, n_hidden=30, eta=0.0001, n_iter=100, minibatch_size=1, l2=0):\n", " # Learn weights from training data\n", " # X_train: input array [n_samples, n_features]\n", " # Y_train: target class labels [n_samples]\n", " # n_hidden: number of units in hidden layer\n", " # eta: learning rate\n", " # n_iter: max number of iterations\n", " # minibatch_size: number of minibatch\n", " # l2: regularization coefficient lambda\n", " # return learned weights\n", "\n", " # 1. Initialize the wegiths\n", "\n", "\n", " n_output = len(np.unique(Y_train))\n", " w_h = np.random.randn(X_train.shape[1] ,n_hidden)\n", " w_out = np.random.randn(n_hidden, n_output)\n", " b_h = np.zeros(n_hidden)\n", " b_out = np.zeros(n_output)\n", "\n", "\n", " # 2. One-hot encoding of Y_train\n", " y_onehot = onehot(Y_train, n_output)\n", "\n", "\n", "\n", " # 3. Do iterations\n", " for i in range(n_iter):\n", "\n", " idx = np.arange(X_train.shape[0])\n", "\n", " # 4. Do iterations for minibatch\n", " for start_idx in range(0, X_train.shape[0]-minibatch_size+1, minibatch_size):\n", "\n", " batch_idx = idx[start_idx:start_idx+minibatch_size]\n", " X_train_batch = X_train[batch_idx]\n", "\n", " # 5. Forward X_train_batch\n", " z_h, a_h, z_out, a_out = forward(X_train_batch, w_h, b_h, w_out, b_out)\n", "\n", " # 6. Set delta_out and grad_w_out, grad_b_out\n", " delta_out = (a_out - y_onehot[batch_idx])\n", " grad_w_out = np.matmul(a_h.T, delta_out) + l2 * w_out\n", " grad_b_out = np.sum(delta_out, axis=0)\n", "\n", " # 7. Set delta_h and grad_w_h, grad_b_h\n", "\n", " delta_h = np.matmul(delta_out, w_out.T) * (a_h * (1 - a_h))\n", " grad_w_h = np.matmul(X_train_batch.T, delta_h) + l2 * w_h\n", " grad_b_h = np.sum(delta_h, axis=0)\n", "\n", " # 8. Update weights: w_h, b_h, w_out, b_out\n", " w_out -= eta *(grad_w_out + l2 * w_out)\n", " b_out -= eta * grad_b_out\n", " w_h -= eta * (grad_w_h + l2 * w_h)\n", " b_h -= eta * grad_b_h\n", "\n", " # 9. Get a_out\n", " a_out = forward(X_train_batch, w_h, b_h, w_out, b_out)[3]\n", "\n", "\n", "\n", " # 10. Compute cost\n", " cost = compute_cost(y_onehot, a_out, l2, w_h, w_out)\n", "\n", "\n", "\n", " # 11. Predict y_train_pred\n", " y_train_pred = predict(X_train, w_h, b_h, w_out, b_out)\n", "\n", " # 12. Check the training results\n", " train_acc = ((np.sum(Y_train == y_train_pred)).astype(float) / X_train.shape[0])\n", " sys.stderr.write('%d/%d | Cost: %.2f '\n", " '| Train Acc.: %.2f%% \\n' %\n", " (i+1, n_iter, cost, train_acc*100))\n", " sys.stderr.flush()\n", "\n", " return w_h, b_h, w_out, b_out" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 146 }, "id": "9UpXaLX91PUz", "outputId": "03520328-f8b4-44a0-c329-287e74265977" }, "outputs": [ { "output_type": "display_data", "data": { "text/plain": [ "
" ], "image/png": "iVBORw0KGgoAAAANSUhEUgAAAnUAAACBCAYAAACma0xyAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAASi0lEQVR4nO3dd3BU5ffH8ZsgLZBQBClDLEEYLDQHEQyEJihFEAwzCAjSlKLCgKJIGZHi0AQnooJikCIKERBmkOLgADaKRBGMjKAQIDQdY5ZAaMnvD8f93XP0e3fX7N3sPvt+/XU/c3ezz3iXzePdk3NiioqKiiwAAABEtNiSXgAAAACKj00dAACAAdjUAQAAGIBNHQAAgAHY1AEAABiATR0AAIAB2NQBAAAY4AZ/HlRYWGjl5ORY8fHxVkxMjNtrQpAUFRVZHo/Hql27thUbG9j+nWsembjm0YdrHn245tHH32vu16YuJyfHSkxMDNriEFonTpyw6tSpE9BzuOaRjWsefbjm0YdrHn18XXO/NnXx8fHeH5aQkBCclcF1eXl5VmJiovf6BYJrHpm45tGHax59uObRx99r7tem7u9btAkJCbwJItB/ucXONY9sXPPowzWPPlzz6OPrmvOHEgAAAAZgUwcAAGAANnUAAAAGYFMHAABgADZ1AAAABmBTBwAAYAA2dQAAAAZgUwcAAGAANnUAAAAGYFMHAABgADZ1AAAABmBTBwAAYAA2dQAAAAa4oaQXAABAOLh48aLIDz30kMhXrlwR+ZtvvnF9TZAKCgpEPnv2rMirVq0See/evSKvXbtW5NhYeW9r1KhRIk+aNEnkm266yf/FlgDu1AEAABiATR0AAIAB2NQBAAAYwMiaury8PJHPnz8vcpkyZUROTEwUOT8/X+QKFSo4vl5ubq7I6enpImdmZnqPV6xYIc61atVK5O3bt4t8ww1GXiIACLmcnByRJ0yYIPLq1atFjouLE3n27NnuLAx+a9GihcgHDx4M6Pm6hi4mJkbkN9980zHba/SaNm0a0GuHAnfqAAAADMCmDgAAwABs6gAAAAwQEQVbhYWFIn/++ecijxs3TuRz586JfObMGZHLli0rcoMGDUTWNXjVq1d3XJ+u09DPt9fF9evXT5zT39dTQxeedN3kvn37RF6+fPn/zLfffrs4l5GRIXKjRo2CsEIg8n377bciT58+XeQvvvjC8fl33HGHyFlZWSLreuurV6+KXKVKFZFnzZol8pAhQxxfH8Wnr9FTTz0l8g8//CByjRo1RO7Zs6fIuu9c+fLlRda/r1955RWRN2/eLPKIESO8x1u3bhXnEhISrJLGnToAAAADsKkDAAAwQER817d+/XqRU1NTHR9frlw5kQcNGhTQ6xUVFYms/+RZ07dcu3XrJrL9z+b1V3EoGZcvXxb5119/FXnq1Kki66/89S17p/fM0aNHxbnx48eLrG/v4y8XLlwQeePGjSJv2rRJ5P3794v8448/ity1a1eR9b/T0qVLi9ynTx+RfbU2QuBOnTolcq9evUQ+ceKE4/P1yCZfX88OHjxY5A4dOois3xPx8fGOPw/Bt2TJEpHXrFkjsv59u2vXLpHr1q0b0OvddtttIq9bt07kevXqiWxvaaI/Y3S7lZLAnToAAAADsKkDAAAwAJs6AAAAA0RETZ2v76n1KJdt27aJ3LJly6Cvye7atWsi05bEfb///rvIHo9H5FtvvVXkn3/+WeTOnTuLrGvqfOnevbvIkyZNErl3797e4+PHj4tzv/32W0CvFS2+++47kYcNGyZyQUGByIcOHXL8eUlJSSJnZ2eL/Pbbb4t85MgRkfV7asyYMY6vh8ClpaWJrGvoFixYIHL//v1F1u2pdK2sVrVq1QBXCLfp358bNmwQWdcrHzhwQGQ95rO4dG2tHium1xNuuFMHAABgADZ1AAAABmBTBwAAYICIKP7So7Q0+9gOy3K/hk6jhs59upfb8OHDRdaj4UaOHCnywoULRda1N3rUTPv27UWeMmWKyLrfYKlSpURu06aN93jZsmXiXEpKihWt7CP19EioFStWiKz7zjVt2lRkXY/VpUsXkWvVqiWy7l+pffXVVyL36NFD5GbNmoncqlUrx5+Hf9JjufR4PX3N9FguX70C6SUYeXSN8ffffy/ygAEDRA52DZ0vuk+tr761JY07dQAAAAZgUwcAAGAANnUAAAAGCMtisEuXLomsa+p0n5gXXnjB8efpeitdTxXq7+jhW2Zmpsi6XspXXcP8+fMdz+s+c/o9Vrt2bV9LdKTrQux0bZjJdu/eLbK9Tk3/u9SmTZsm8sSJE4O3sH9x7733iqx7IZ48edLV148Gujfg6dOnRX711VdF1jVyV65cEfnYsWMily9fXuTKlSuLzCzX8FOzZk2R9Xzetm3bhnA1lpWVlSXymTNnRLbXfTZu3DgkawoEd+oAAAAMwKYOAADAAGzqAAAADBCWNXW6P1Vubq7Iup7queeeE1nXM+m5nrqmrkGDBiLPnTtX5Pvvv19kXbeB4NM1dYHO23v55ZdF7tu3r8i6z1xx6fdUXl6e91ivPTk5OaivHc4ee+wxkZ3q6PRMxzvvvNOVNf0vut+kntmrZ83m5+eLbK/F0Z859lnA0czX3ONKlSqJrGc2z5kzR+R33303oNcfNWqUyHqWbJMmTUTWs2XhPl1XqXuAuk33RNWf7Q8++KD3OBz3AtypAwAAMACbOgAAAAOwqQMAADBAWNbULVmyJKDH69mad999t8iPPvqoyBcuXBB5zZo1Infs2FFkXVexfft2kUM9a9ZE+/btE3n06NEi6zrKihUriqxnSD788MMix8a6+/8vAwcOFNleU9WoUSNxjr6I/85X3aTuaXb48GGRfV1j3YeuTJkyIhcUFIis6wF99anLzs72Huv3Z/369UUOx/5W4WDmzJki6xnOuu5Rf/a2aNFC5EWLFoms++TpmdA33nijyBMmTBBZz5yOi4uzEFylS5cO6evpmv0FCxaIrD+XUlNTXV5R8XCnDgAAwABs6gAAAAzApg4AAMAAYVlTp+f1lStXTmTdVy4jI0PkOnXqiKxrZwoLC0VeunSpyLqO48MPPxS5ffv2Iuu+dk8++aTIoa4RiETvv/++yBcvXnR8vO5DZ58pWhJWr14tsr0GUK81mt4Pd911l8h6VqedrjPT83n37t0rsq6xC5S935RlWdaWLVuK9fPsdN3uzp07Raam7t/pukXdt279+vUit2nTxvHnzZs3T2R73aNlWdbatWtF1jV0ugeqrsHbtWuX97i486IRGteuXRP5tddeE/nUqVMi33zzzSJ37tzZnYUFCXfqAAAADMCmDgAAwABs6gAAAAwQljV1M2bMEHns2LEi65q5QOl+Vnp+W3p6usiTJ08WuVOnTiI/88wzIusZlrpXEv5i7xml+0dpui5tzJgxLqzIf5988onj+fj4eO9xSkqK28sJW7rW0D57c9OmTeKcngu7YcMG9xZmBbeGzpeaNWuG7LUiie4Bpuup9+zZI3K9evWK9Xq6Pkp/jgwbNkxk3TdPzyW1zyfeuHGjONe6dev/uky4aNWqVSLra6ylpaW5uZyg404dAACAAdjUAQAAGIBNHQAAgAHCsqZO17gVt4auuJKSkkQ+cuSIyM2aNRN58eLFIjdv3tx7PGTIkCCvLnLo/kD22gbdO1DPSx03bpzIbs9y1f7880+RBw8e7Ph4e++jqlWrurKmSKD/Lb/33nveY93L7YMPPhD58uXLIv/yyy8i65q7hIQEkfU10o/XffCysrJEzszMFHn37t2Wv/T7t2vXrn4/12S675uuNdQzm0P92V+hQgWRJ06c6Ph4e43dgAEDxDn7/GeEj2effdbxvP63+8ADD7i5nKDjTh0AAIAB2NQBAAAYgE0dAACAAcKypi7S6PmC7dq1E3ndunXe42iuqRs6dKjI9rmJupZF963T59126dIlkVNTU0X+448/RO7SpYvIffv2dWdhBqlYsaLIemayL/Pnzw/o8bqfpC96BrSvmjp7D7UdO3aIc3FxcQG9tqns/Rsty7JGjBhRQivxj75ukyZNEnnlypXe4xMnTohzugZU12ZHC4/HI7Ken6v99NNPIi9fvlxk3ffzlltuEVnPhtczo3V9tO5daP99bVn/nB0f7rhTBwAAYAA2dQAAAAZgUwcAAGAAauqCwD7/D/8vLy9P5G3btokcExPjPbb3dbMsy7rvvvvcW9i/0DMo7f3ULMuytm/fLrJ97ZZlWRkZGSKXK1cuiKtDKFy/fl3kzz77LKDn9+/f33tcqVKloKwJ4UX3XaxRo4b3ODs7W5zT7yeT6Z6Sc+bM8R7Pnj1bnMvPzxdZf5Zq+vzOnTsDer6vn3fy5EmR69atK7L937Vlyf6WTZo0Eed0r8xq1aoFtLZg4E4dAACAAdjUAQAAGIBNHQAAgAGoqQsC3VMNfyldurTIus6oc+fO3mM9NzHUvvzyS5F9zQccP368yNTQRb6DBw+KrGc8ay+99JLIL774YtDXBISjq1eviqx7QNprkvUs1Y4dO4o8fPjwgF77wIEDIuuZ0R9//LHj85OTk0Vu2LChyMePHxd5xYoVItv75un6vOrVq4usP0NC0W+VO3UAAAAGYFMHAABgADZ1AAAABgjLmjo9m23YsGEiv/766yLXqlXL9TXZbdmyReTp06c7Pr5Xr15uLids6Z5O+/fvF7kk69B0TUh6errj4/V7TNdTIfLo3oTPP/+84+NjY+X/A+v+VbqGFObRPc3sdZi6ZljXV5nknXfeEVn39Xz88ce9x2+99ZY4V9zP/WnTpon86aefOj5+w4YNInfo0EHksmXLiqx/N+g54HZ6BvhHH30ksp4HrOfSuoE7dQAAAAZgUwcAAGAANnUAAAAGCMuaukWLFom8adMmkXVNndtmzZol8oQJExwfr+eYDh48OOhrikTh1MvtjTfeEHnp0qWOjx89erTI8fHxwV4SQmzz5s0i+5r1OnXqVJFDUR9jutzcXJF1fXJKSorI9rmboXD27FmRu3TpIrK93io1NVWcq1y5smvrKml6rreuT50yZYr3ONDPfV3TNnToUJHtfeIs65+94rp16yayrxo6TdfGOtXK6lmvun9pSeBOHQAAgAHY1AEAABggbL5+LSgo8B7PnDlTnNMjpIrbwqSwsFDkY8eOiaxbkBw+fNjx540cOVLkUaNG/ffFwRW7d+8WeezYsSLrW/i6ZYmvdheIPNWqVQvo8S1btnRpJdHD4/GI3Lp1a5EPHTok8sCBA11dj25XoUdOPf300yJfvnxZ5MaNG3uPFy9eHOTVRQ79+WkvZ5k8ebI4l52dLfL58+dF1mUOW7dudXwt/ft62bJlIvv6utU03KkDAAAwAJs6AAAAA7CpAwAAMEDY1NTZvyfXfyasv4Pft2+fyDVr1hT5zJkzjq+lW6KsXLnS8fH6T7IzMjJE7tGjh8ilSpVy/Hlwn66V0XWZui6jZ8+eItv/JB9mWrNmjeP5e+65R2Rd/4XA2WunLeufY7e04cOHi6zrrcqUKSNy/fr1RV67dq3Iu3btEll/ll+/fl1kPRquT58+Ittrx/RaTPbEE0+InJmZKfKMGTO8x/q/+Y4dO0TWn8W+6GvaqVMnkfV4ymjDnToAAAADsKkDAAAwAJs6AAAAA4RNTZ29l8yePXvEOT32o3nz5iLrmjddtxEoXSOXnp4ussnjX0yxcOFCkY8cOeL4+OTkZJGjqT4mWuTn54us+1lpaWlpIjuNC4J/qlevLvK6detEfuSRR0T++uuvRdb1U8HWrl07kRcsWCByw4YNXX39SKHr2PVoz7i4OO+xHuulJSYmijxixAiRBw0aJLJ+D0HiTh0AAIAB2NQBAAAYgE0dAACAAcKmps5Of1+v6yp0vdTp06dF1rNaGzRoIHKVKlVE7t27t8hJSUn+LxZhad68eQE9Pisry6WVIFzoXlrnzp0TWfeX1J8TCL62bduKnJubWyLrQPHoOnN7HbquSYe7uFMHAABgADZ1AAAABmBTBwAAYICwrKnTdH+oMWPGlMxCEDH69esn8vz580Xu3r27yHqmJMyje6Jp48ePF1nX4gJAuONOHQAAgAHY1AEAABiATR0AAIABIqKmDgjU3LlzHTOij57pvH79epEnT54cwtUAQPBxpw4AAMAAbOoAAAAMwKYOAADAANTUAYgKKSkpIh89erSEVgIA7uBOHQAAgAHY1AEAABjAr69fi4qKLMuyrLy8PFcXg+D6+3r9ff0CwTWPTFzz6MM1jz5c8+jj7zX3a1Pn8Xgsy7KsxMTEYi4LJcHj8ViVKlUK+DmWxTWPVFzz6MM1jz5c8+jj65rHFPmx1S8sLLRycnKs+Ph4KyYmJqgLhHuKioosj8dj1a5d24qNDeybdq55ZOKaRx+uefThmkcff6+5X5s6AAAAhDf+UAIAAMAAbOoAAAAMwKYOAADAAGzqAAAADMCmDgAAwABs6gAAAAzApg4AAMAA/wd1o8bC8k3H8gAAAABJRU5ErkJggg==\n" }, "metadata": {} } ], "source": [ "# load the training data\n", "data = pd.read_csv('/content/sample_data/mnist_train_small.csv')\n", "data = np.array(data)\n", "Y_train = data[:,0]\n", "X_train = data[:,1:]/255\n", "\n", "# Visualize the input data\n", "fig, ax = plt.subplots(nrows=1, ncols=5, sharex=True, sharey=True)\n", "ax = ax.flatten()\n", "for i in range(5):\n", " img = X_train[i,:].reshape(28,28)\n", " ax[i].imshow(img,cmap='Greys',)\n", "ax[0].set_xticks([])\n", "ax[0].set_yticks([])\n", "plt.tight_layout()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "-lxUAyS81PU0", "colab": { "base_uri": "https://localhost:8080/" }, "outputId": "34c51da7-279b-4055-ae10-c98784220ae3" }, "outputs": [ { "output_type": "stream", "name": "stderr", "text": [ "1/20 | Cost: 9120690.70 | Train Acc.: 40.38% \n", "2/20 | Cost: 9409831.18 | Train Acc.: 57.23% \n", "3/20 | Cost: 9805421.97 | Train Acc.: 66.06% \n" ] } ], "source": [ "# Train the ANN\n", "w_h, b_h, w_out, b_out = train(X_train, Y_train, n_hidden=30, eta=0.001, n_iter=20, minibatch_size=100, l2=0.1)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "id": "T2THY7sf1PU1" }, "outputs": [], "source": [ "# load the test data\n", "data = pd.read_csv('/content/sample_data/mnist_test.csv')\n", "data = np.array(data)\n", "Y_test = data[:,0]\n", "X_test = data[:,1:]/255\n", "\n", "\n", "y_test_pred = predict(X_test, w_h, b_h, w_out, b_out)\n", "acc = (np.sum(Y_test == y_test_pred).astype(float) / X_test.shape[0])\n", "print('Test accuracy: %.2f%%' % (acc * 100))" ] } ], "metadata": { "kernelspec": { "display_name": "ML4ME", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.8.18" }, "colab": { "provenance": [] } }, "nbformat": 4, "nbformat_minor": 0 }