Machine-Learning/implement-softmax-in-theano.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 用theano实现softmax分类器\n",
    "\n",
    "上一篇笔记[softmax分类器](./softmax-crossentropy-derivative.ipynb)介绍了softmax分类器的基本原理，涉及交叉熵误差及其对应的梯度更新式的推导过程。这一篇笔记着重于softmax分类器的实现，所用的工具为theano，测试的数据集为MNIST手写数字数据集，该数据集包含60000张训练图片，10000张测试图片。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "载入依赖包"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "DEBUG: nvcc STDOUT mod.cu\r\n",
      "   ���ڴ����� C:/Users/hschen/AppData/Local/Theano/compiledir_Windows-10-10.0.14393-Intel64_Family_6_Model_60_Stepping_3_GenuineIntel-2.7.12-64/tmpgcrtfx/265abc51f7c376c224983485238ff1a5.lib �Ͷ��� C:/Users/hschen/AppData/Local/Theano/compiledir_Windows-10-10.0.14393-Intel64_Family_6_Model_60_Stepping_3_GenuineIntel-2.7.12-64/tmpgcrtfx/265abc51f7c376c224983485238ff1a5.exp\r\n",
      "\n",
      "Using gpu device 0: GeForce GTX 960M (CNMeM is disabled, cuDNN 5103)\n",
      "C:\\Anaconda2\\lib\\site-packages\\theano\\sandbox\\cuda\\__init__.py:600: UserWarning: Your cuDNN version is more recent than the one Theano officially supports. If you see any problems, try updating Theano or downgrading cuDNN to version 5.\n",
      "  warnings.warn(warn)\n",
      "Using Theano backend.\n"
     ]
    }
   ],
   "source": [
    "from __future__ import print_function, division\n",
    "from abc import ABCMeta, abstractmethod\n",
    "from six import add_metaclass\n",
    "import numpy as np\n",
    "import theano\n",
    "import theano.tensor as T\n",
    "from keras.datasets import mnist\n",
    "from collections import OrderedDict"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "随机数生成器和其他的工具函数"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "class NumpyRNG(object):\n",
    "    _rng = None\n",
    "\n",
    "    @classmethod\n",
    "    def get_rng(cls):\n",
    "        if cls._rng == None:\n",
    "            cls._rng = np.random\n",
    "        return cls._rng\n",
    "\n",
    "    @classmethod\n",
    "    def set_rng(cls, seed):\n",
    "        cls._rng = np.random.RandomState(seed)\n",
    "\n",
    "\n",
    "def floatX(arr):\n",
    "    return np.asarray(arr, dtype=theano.config.floatX)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "初始化模块，用于初始化权重"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "@add_metaclass(ABCMeta)\n",
    "class Initializer(object):\n",
    "    @abstractmethod\n",
    "    def create_param(self, size):\n",
    "        \"\"\"\"\"\"\n",
    "\n",
    "    def __call__(self, size, shared = True):\n",
    "        param = self.create_param(size)\n",
    "        if shared:\n",
    "            return theano.shared(param)\n",
    "        else:\n",
    "            return param\n",
    "\n",
    "class Constant(Initializer):\n",
    "    def __init__(self, scale = 0.):\n",
    "        self.scale = scale\n",
    "\n",
    "    def create_param(self, size):\n",
    "        param = floatX(np.ones(size) * self.scale)\n",
    "        return param\n",
    "\n",
    "class Normal(Initializer):\n",
    "    def __init__(self, mean = 0., std = .01):\n",
    "        self.mean = mean\n",
    "        self.std = std\n",
    "\n",
    "    def create_param(self, size):\n",
    "        rng = NumpyRNG.get_rng()\n",
    "        param = floatX(rng.normal(loc=self.mean, scale=self.std, size=size))\n",
    "        return param"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "softmax层的实现"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "class Layer(object):\n",
    "    def __init__(self, **kwargs):\n",
    "        \"\"\"\"\"\"\n",
    "\n",
    "\n",
    "    @abstractmethod\n",
    "    def apply(self, x, mask=None):\n",
    "        \"\"\"\"\"\"\n",
    "\n",
    "    def __call__(self, x, mask=None):\n",
    "        return self.apply(x, mask)\n",
    "class Logistic(Layer):\n",
    "    def __init__(self, n_in, n_out, **kwargs):\n",
    "        self.n_in = n_in\n",
    "        self.n_out = n_out\n",
    "        self.W = Normal(0., .01)((n_in, n_out))\n",
    "        self.b = Constant()(n_out)\n",
    "        self.params = [self.W, self.b]\n",
    "\n",
    "    def apply(self, x, mask=None):\n",
    "        p_y_given_x = T.nnet.softmax(T.dot(x, self.W) + self.b)\n",
    "        return p_y_given_x\n",
    "\n",
    "    def cost(self, x, label):\n",
    "        if label.ndim != 1 or x.ndim != 2:\n",
    "            raise NotImplementedError\n",
    "        prediction = self.apply(x)\n",
    "        cost = T.mean(-T.log(T.clip(prediction[T.arange(x.shape[0]), label], 1e-3, 1.)))\n",
    "        return cost"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "优化算法为adadelta，借用了lasagne的代码。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def adadelta(params, grads, learning_rate=1.0, rho=0.95, epsilon=1e-6):\n",
    "    updates = OrderedDict()\n",
    "\n",
    "    # Using theano constant to prevent upcasting of float32\n",
    "    one = T.constant(1)\n",
    "\n",
    "    for param, grad in zip(params, grads):\n",
    "        value = param.get_value(borrow=True)\n",
    "        # accu: accumulate gradient magnitudes\n",
    "        accu = theano.shared(np.zeros(value.shape, dtype=value.dtype),\n",
    "                             broadcastable=param.broadcastable)\n",
    "        # delta_accu: accumulate update magnitudes (recursively!)\n",
    "        delta_accu = theano.shared(np.zeros(value.shape, dtype=value.dtype),\n",
    "                                   broadcastable=param.broadcastable)\n",
    "\n",
    "        # update accu (as in rmsprop)\n",
    "        accu_new = rho * accu + (one - rho) * grad ** 2\n",
    "        updates[accu] = accu_new\n",
    "\n",
    "        # compute parameter update, using the 'old' delta_accu\n",
    "        update = (grad * T.sqrt(delta_accu + epsilon) /\n",
    "                  T.sqrt(accu_new + epsilon))\n",
    "        updates[param] = param - learning_rate * update\n",
    "\n",
    "        # update delta_accu (as accu, but accumulating updates)\n",
    "        delta_accu_new = rho * delta_accu + (one - rho) * update ** 2\n",
    "        updates[delta_accu] = delta_accu_new\n",
    "\n",
    "    return updates"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "生成batch的index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def get_minibatches_idx(n, minibatch_size, shuffle=False):\n",
    "    \"\"\"\n",
    "    Used to shuffle the dataset at each iteration.\n",
    "    \"\"\"\n",
    "\n",
    "    idx_list = np.arange(n, dtype=\"int32\")\n",
    "    if shuffle:\n",
    "        np.random.shuffle(idx_list)\n",
    "\n",
    "    minibatches = []\n",
    "    minibatch_start = 0\n",
    "    for i in range(n // minibatch_size):\n",
    "        minibatches.append(idx_list[minibatch_start:\n",
    "                                    minibatch_start + minibatch_size])\n",
    "        minibatch_start += minibatch_size\n",
    "    if (minibatch_start != n):\n",
    "        # Make a minibatch out of what is left\n",
    "        minibatches.append(idx_list[minibatch_start:])\n",
    "    return zip(range(len(minibatches)), minibatches)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "训练代码。读取数据的部分借用了keras的模块。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "epoch 0, trn_loss=0.286183655262, trn_acc=0.918216666667, val_loss=0.280671298504, val_acc=0.9188\n",
      "epoch 1, trn_loss=0.267140179873, trn_acc=0.9241, val_loss=0.26878619194, val_acc=0.923\n",
      "epoch 2, trn_loss=0.258325606585, trn_acc=0.926983333333, val_loss=0.264801651239, val_acc=0.925\n",
      "epoch 3, trn_loss=0.252924561501, trn_acc=0.928833333333, val_loss=0.263516664505, val_acc=0.9257\n",
      "epoch 4, trn_loss=0.249107733369, trn_acc=0.93035, val_loss=0.263270288706, val_acc=0.9264\n",
      "epoch 5, trn_loss=0.246360614896, trn_acc=0.931616666667, val_loss=0.263392955065, val_acc=0.9267\n",
      "epoch 6, trn_loss=0.244001120329, trn_acc=0.9324, val_loss=0.263515710831, val_acc=0.9266\n",
      "epoch 7, trn_loss=0.242098540068, trn_acc=0.932816666667, val_loss=0.263630479574, val_acc=0.9274\n",
      "epoch 8, trn_loss=0.240467473865, trn_acc=0.93315, val_loss=0.263916522264, val_acc=0.927\n",
      "epoch 9, trn_loss=0.239036202431, trn_acc=0.93375, val_loss=0.264052331448, val_acc=0.9272\n",
      "epoch 10, trn_loss=0.237787619233, trn_acc=0.93405, val_loss=0.264259874821, val_acc=0.9277\n",
      "epoch 11, trn_loss=0.236647352576, trn_acc=0.934566666667, val_loss=0.264347076416, val_acc=0.9278\n",
      "epoch 12, trn_loss=0.235716700554, trn_acc=0.934883333333, val_loss=0.264674186707, val_acc=0.9279\n",
      "epoch 13, trn_loss=0.234838038683, trn_acc=0.935116666667, val_loss=0.264962106943, val_acc=0.9277\n",
      "epoch 14, trn_loss=0.234002992511, trn_acc=0.9354, val_loss=0.265269070864, val_acc=0.9275\n",
      "epoch 15, trn_loss=0.233251526952, trn_acc=0.9356, val_loss=0.265526205301, val_acc=0.927\n",
      "epoch 16, trn_loss=0.23260357976, trn_acc=0.935733333333, val_loss=0.265771985054, val_acc=0.9266\n",
      "epoch 17, trn_loss=0.231973141432, trn_acc=0.93585, val_loss=0.266015022993, val_acc=0.9265\n",
      "epoch 18, trn_loss=0.231382533908, trn_acc=0.936116666667, val_loss=0.266254454851, val_acc=0.9266\n",
      "epoch 19, trn_loss=0.230891555548, trn_acc=0.93615, val_loss=0.266672164202, val_acc=0.9266\n",
      "epoch 20, trn_loss=0.230336993933, trn_acc=0.93645, val_loss=0.26688337326, val_acc=0.9268\n",
      "epoch 21, trn_loss=0.229838177562, trn_acc=0.936733333333, val_loss=0.267110854387, val_acc=0.9267\n",
      "epoch 22, trn_loss=0.22934114933, trn_acc=0.936883333333, val_loss=0.267413973808, val_acc=0.9268\n",
      "epoch 23, trn_loss=0.228895649314, trn_acc=0.937116666667, val_loss=0.267637491226, val_acc=0.9273\n",
      "epoch 24, trn_loss=0.228484258056, trn_acc=0.937183333333, val_loss=0.26788431406, val_acc=0.9271\n",
      "epoch 25, trn_loss=0.228089660406, trn_acc=0.9373, val_loss=0.26819974184, val_acc=0.927\n",
      "epoch 26, trn_loss=0.22770345211, trn_acc=0.937566666667, val_loss=0.268486142159, val_acc=0.9274\n",
      "epoch 27, trn_loss=0.227417364717, trn_acc=0.937433333333, val_loss=0.268964886665, val_acc=0.927\n",
      "epoch 28, trn_loss=0.226861461997, trn_acc=0.937533333333, val_loss=0.269212245941, val_acc=0.9267\n",
      "epoch 29, trn_loss=0.226534023881, trn_acc=0.93775, val_loss=0.269475758076, val_acc=0.9266\n"
     ]
    }
   ],
   "source": [
    "batch_size = 32\n",
    "n_epoch = 30\n",
    "(X_train, y_train), (X_test, y_test) = mnist.load_data(path=r\"G:\\data\\mnist.pkl.gz\")\n",
    "X_train = X_train.reshape(len(X_train), 28 * 28) / np.float32(255.)\n",
    "X_test = X_test.reshape(len(X_test), 28*28) / np.float32(255.)\n",
    "X_tr = T.matrix(dtype=theano.config.floatX)\n",
    "y_tr = T.vector(dtype='int32')\n",
    "lr = Logistic(784, 10)\n",
    "lr_cost = lr.cost(X_tr, y_tr)\n",
    "grads = T.grad(lr_cost, lr.params)\n",
    "updates = adadelta(lr.params, grads)\n",
    "acc = T.mean(T.eq(T.argmax(lr(X_tr), axis=1), y_tr))\n",
    "fn_train = theano.function([X_tr, y_tr], [lr_cost, acc], updates=updates)\n",
    "\n",
    "mini_batches = get_minibatches_idx(len(X_train), batch_size, shuffle=True)\n",
    "for e in range(n_epoch):\n",
    "    for _, train_idx in mini_batches:\n",
    "        cost = fn_train(X_train[train_idx], y_train[train_idx])\n",
    "    trn_loss,trn_acc = fn_train(X_train, y_train)\n",
    "    val_loss,val_acc = fn_train(X_test, y_test)\n",
    "    print(\"epoch {}, trn_loss={}, trn_acc={}, val_loss={}, val_acc={}\"\\\n",
    "          .format(e, trn_loss, trn_acc, val_loss, val_acc))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 一些总结\n",
    "1. 必须对数据进行预处理，可以简单地将像素值除以255，放缩为[0,1]区间内的值。  \n",
    "```py\n",
    "X_train = X_train.reshape(len(X_train), 28 * 28) / np.float32(255.)\n",
    "X_test = X_test.reshape(len(X_test), 28*28) / np.float32(255.)\n",
    "```\n",
    "这一步很重要，如果不进行预处理，结果准确度大约是0.4~0.6  \n",
    "2. 进行预处理后的结果大约为0.9266，这与UFLDL中给出的0.926一致，说明实现是正确的"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 参考\n",
    "1. [UFLDL exercise-softmax-regression](http://deeplearning.stanford.edu/wiki/index.php/Exercise:Softmax_Regression)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}