{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 内容概要\n",
    "- 训练集/测试集分割用于模型验证的缺点\n",
    "- K折交叉验证是如何克服之前的不足\n",
    "- 交叉验证如何用于选择调节参数、选择模型、选择特征\n",
    "- 改善交叉验证"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. 模型验证回顾\n",
    "进行模型验证的一个重要目的是要选出一个最合适的模型，对于监督学习而言，我们希望模型对于未知数据的泛化能力强，所以就需要模型验证这一过程来体现不同的模型对于未知数据的表现效果。\n",
    "\n",
    "最先我们用训练准确度（用全部数据进行训练和测试）来衡量模型的表现，这种方法会导致模型过拟合；为了解决这一问题，我们将所有数据分成训练集和测试集两部分，我们用训练集进行模型训练，得到的模型再用测试集来衡量模型的预测表现能力，这种度量方式叫测试准确度，这种方式可以有效避免过拟合。\n",
    "\n",
    "**测试准确度的一个缺点**是其样本准确度是一个高方差估计（high variance estimate），所以该样本准确度会依赖不同的测试集，其表现效果不尽相同。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "### 高方差估计的例子\n",
    "下面我们使用iris数据来说明利用测试准确度来衡量模型表现的方差很高。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.datasets import load_iris\n",
    "from sklearn.cross_validation import train_test_split\n",
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn import metrics"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# read in the iris data\n",
    "iris = load_iris()\n",
    "\n",
    "X = iris.data\n",
    "y = iris.target"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "random_state is  1 , and accuracy score is:\n",
      "1.0\n",
      "random_state is  2 , and accuracy score is:\n",
      "1.0\n",
      "random_state is  3 , and accuracy score is:\n",
      "0.947368421053\n",
      "random_state is  4 , and accuracy score is:\n",
      "0.973684210526\n"
     ]
    }
   ],
   "source": [
    "for i in xrange(1,5):\n",
    "    print \"random_state is \", i,\", and accuracy score is:\"\n",
    "    X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=i)\n",
    "\n",
    "    knn = KNeighborsClassifier(n_neighbors=5)\n",
    "    knn.fit(X_train, y_train)\n",
    "    y_pred = knn.predict(X_test)\n",
    "    print metrics.accuracy_score(y_test, y_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "上面的测试准确率可以看出，不同的训练集、测试集分割的方法导致其准确率不同，而交叉验证的基本思想是：将数据集进行一系列分割，生成一组不同的训练测试集，然后分别训练模型并计算测试准确率，最后对结果进行平均处理。这样来有效降低测试准确率的差异。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. K折交叉验证\n",
    "1. 将数据集平均分割成K个等份\n",
    "2. 使用1份数据作为测试数据，其余作为训练数据\n",
    "3. 计算测试准确率\n",
    "4. 使用不同的测试集，重复2、3步骤\n",
    "5. 对测试准确率做平均，作为对未知数据预测准确率的估计\n",
    "![](Image/cross_validation_diagram.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Iteration                   Training set observations                   Testing set observations\n",
      "    1     [ 5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24]        [0 1 2 3 4]       \n",
      "    2     [ 0  1  2  3  4 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24]        [5 6 7 8 9]       \n",
      "    3     [ 0  1  2  3  4  5  6  7  8  9 15 16 17 18 19 20 21 22 23 24]     [10 11 12 13 14]     \n",
      "    4     [ 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 20 21 22 23 24]     [15 16 17 18 19]     \n",
      "    5     [ 0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19]     [20 21 22 23 24]     \n"
     ]
    }
   ],
   "source": [
    "# 下面代码演示了K-fold交叉验证是如何进行数据分割的\n",
    "# simulate splitting a dataset of 25 observations into 5 folds\n",
    "from sklearn.cross_validation import KFold\n",
    "kf = KFold(25, n_folds=5, shuffle=False)\n",
    "\n",
    "# print the contents of each training and testing set\n",
    "print '{} {:^61} {}'.format('Iteration', 'Training set observations', 'Testing set observations')\n",
    "for iteration, data in enumerate(kf, start=1):\n",
    "    print '{:^9} {} {:^25}'.format(iteration, data[0], data[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. 使用交叉验证的建议\n",
    "1. K=10是一个一般的建议\n",
    "2. 如果对于分类问题，应该使用分层抽样（stratified sampling）来生成数据，保证正负例的比例在训练集和测试集中的比例相同"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. 交叉验证的例子\n",
    "### 4.1 用于调节参数\n",
    "交叉验证的方法可以帮助我们进行调参，最终得到一组最佳的模型参数。下面的例子我们依然使用iris数据和KNN模型，通过调节参数，得到一组最佳的参数使得测试数据的准确率和泛化能力最佳。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from sklearn.cross_validation import cross_val_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 1.          0.93333333  1.          1.          0.86666667  0.93333333\n",
      "  0.93333333  1.          1.          1.        ]\n"
     ]
    }
   ],
   "source": [
    "knn = KNeighborsClassifier(n_neighbors=5)\n",
    "# 这里的cross_val_score将交叉验证的整个过程连接起来，不用再进行手动的分割数据\n",
    "# cv参数用于规定将原始数据分成多少份\n",
    "scores = cross_val_score(knn, X, y, cv=10, scoring='accuracy')\n",
    "print scores"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.966666666667\n"
     ]
    }
   ],
   "source": [
    "# use average accuracy as an estimate of out-of-sample accuracy\n",
    "# 对十次迭代计算平均的测试准确率\n",
    "print scores.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0.95999999999999996, 0.95333333333333337, 0.96666666666666656, 0.96666666666666656, 0.96666666666666679, 0.96666666666666679, 0.96666666666666679, 0.96666666666666679, 0.97333333333333338, 0.96666666666666679, 0.96666666666666679, 0.97333333333333338, 0.98000000000000009, 0.97333333333333338, 0.97333333333333338, 0.97333333333333338, 0.97333333333333338, 0.98000000000000009, 0.97333333333333338, 0.98000000000000009, 0.96666666666666656, 0.96666666666666656, 0.97333333333333338, 0.95999999999999996, 0.96666666666666656, 0.95999999999999996, 0.96666666666666656, 0.95333333333333337, 0.95333333333333337, 0.95333333333333337]\n"
     ]
    }
   ],
   "source": [
    "# search for an optimal value of K for KNN model\n",
    "k_range = range(1,31)\n",
    "k_scores = []\n",
    "for k in k_range:\n",
    "    knn = KNeighborsClassifier(n_neighbors=k)\n",
    "    scores = cross_val_score(knn, X, y, cv=10, scoring='accuracy')\n",
    "    k_scores.append(scores.mean())\n",
    "\n",
    "print k_scores"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.text.Text at 0x6dd0fb0>"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAEPCAYAAACDTflkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xu8lWWd9/HPl60IiAcIUWSjkAKiKCCJbKxpl1rUZKWV\nRr5GK8ecKdMpm9HsqbB6ntEanSyfytLMmScPPSVl05hauX3G4SxbQAURBUNAVAhBUTn9nj+u+5bF\nch3udbjX8fd+vfZrr3UfrnXdLPb6rev63dd1ycxwzjnnytGn3hVwzjnXvDyIOOecK5sHEeecc2Xz\nIOKcc65sHkScc86VzYOIc865sqUaRCRNl7Rc0pOSLs+xf5CkWZIWS5on6biMfV+Q9KikpZJuk7Rf\ntH2wpPslrZB0n6SD07wG55xz+aUWRCR1ADcA04FjgRmSxmUddiWwyMwmAOcB10fnDgc+D0w2s+OB\nDuDj0TlXAPeb2Rjgj9Fz55xzdZBmS2QKsNLMVpvZDuAO4ENZx4wDHgAwsyeAkZIOifbtAwyQtA8w\nAFgbbf8gcGv0+Fbgw+ldgnPOuULSDCLDgTUZz5+NtmVaDJwFIGkKcCTQaWZrgWuBPwPrgJfM7A/R\nOYea2Ybo8Qbg0HSq75xzrpg0g0iS+VSuBg6W1AtcDPQCuyQNIrQ4RgKHA/tLOvdNLxDmbPF5W5xz\nrk72SbHstcCIjOcjCK2RN5jZVuDT8XNJq4CngfcBq8xsY7T9LmAa8HNgg6TDzOw5ScOA53O9uCQP\nLs45VyIzUynHp9kSWQiMljRSUl/gHODuzAMkHRTtQ9KFwINm9jLwDDBVUn9JAk4DHo9Ouxs4P3p8\nPvDrfBUws5b8+frXv173OrTy9f3ud6GBe9FFrXl93/hGuL4f/7i65T76aCi3u7u+15f2T73fvzR/\nypFaEDGznYQuqnsJAeBOM1sm6SJJF0WHHQsslbQceC9waXTufOCXwCJgSXTsj6PfVwOnS1oBvDt6\n7lzVzJ4N73lP+N2K0rq+2bNh333hpZeqW65rbGl2Z2Fm9wD3ZG27MePxHGBsnnNnAjNzbN9EaJk4\nl4o5c+CSS+Dcc2HLFjjwwHrXqHp274a5c+HOO8M1VtOcOfDe98KyZdUt1zU2H7HehLq7u+tdhVTV\n8/p27oT58+Ed74ATT4R586r/GvW8vuXLYfBgOPVUWLcONm6sXtmzZ8PZZ8OuXd3VK7QBtfrfX6k8\niDShVv9PXM/re/RR6OwMH7TTpoVv19VWz+ubMydcV0cHTJkSWiXVsHFjCErTp8OmTd3VKbRBtfrf\nX6k8iDiXYfZs6OoKj7u6Wi8vknl906ZV7/rmzg1BaciQ0JrbsqU65brG50HEuQyzZ4cPVwgftnPn\nhjxCq8i+vmoFkbhcCUaMgDVrip/jWoMHEecyxN09AEOHwiGHtE6ieNMmWLsWxo8Pz6dOhYULQ8uh\nUpn/bh5E2osHEeciGzaED9pjjtmzrZW6tObOhZNOgn2iezIHDQof+EuWFD6vmJ07YcECOPnk8HzE\nCHj22cLnuNbhQcS5yJw54dt5n4y/irSS6/WQ2VqIVeP6liyBI44IQQnCjQneEmkfHkSci2QmnWOt\n1BLJdX3VSK7PmbN3ud6d1V48iDgXyUw6x8aPr/54inqIu5ymTt17ezWCZPa/mweR9uJBxDlg+3bo\n7Q23qWbq6Ah9/dUaT1EvmeNfMo0dG6YpWb++/LKzu8k8iLQXDyLOEQLI0UfnnuKkFbq0cnVlQcj/\ndHWVnxdZvx42b4YxY/ZsixPrZc7n55qMBxHnyJ10jrVCcr3Q9VUSROJ8SObNCAceGJ5v3lxema65\neBBxjvzf1CF0Zy1YUJ3xFPVS6PoqSa5nJ9Vj3qXVPjyIOEfupHps0KBwC2ul4ynqJdf4l0xTpsAj\nj8Drr5dedr5/Nw8i7cODiGt7a9aExPpRR+U/ppm7tHJ1OWUaODDkNHp7Syv39ddh8eI334wAHkTa\niQcR1/birh4VWBS0mZPrhbqyYuV0afX2wujRIQhl81Hr7cODiGt7hZLOsWZviRS7vnKS64W6AH3U\nevvwIOLaXpJv6mPGhLuNKhlPUQ/5xr9ki1sipdyWmy+pDt6d1U48iLi29uqr8Nhj8La3FT6u0vEU\n9RKPfznggMLHjRoFu3Yl/+A3K9wS8SDSPjyIuLa2cCEcdxwMGFD82Gbs0krSlQUhH1RK3mfNmhB0\nRo3Kvd8HHLYPDyKurSXpyoo1Y3K9lOsrJble7GaE/feH/v2bf84xV5wHEdfWkn5Th5BXWLy4vPEU\n9VLK9ZXSXVeoKyvmyfX24EHEta1i/frZyh1PUS/x+Je3vjXZ8ZMnw+OPw7ZtxY8tlFSPeV6kPXgQ\ncW3r6aehb9/wYZdUM3VpJRn/kql/fzj++JAnKmTbthBsJk8ufJwHkfbgQcS1rVJaIbFmSq6X0pUV\nSxIkFy4MwaZ//8LHeRBpDx5EXNsqJekciz9km+Guo3KuL0mQTFquj1pvDx5EXNsq55t6qeMp6iXp\n+JdsSYJk0hacJ9bbgwcR15a2boWVK2HSpNLOk6qzLnnaFi4MS/sW63LK1tkZzlm5Mvd+s2RJdfDu\nrHaRahCRNF3ScklPSro8x/5BkmZJWixpnqTjou1jJfVm/Lwk6ZJo30xJz2bsm57mNbjWNH8+TJwY\nEuulaobkejldWbFCXVorV4Yg09lZvJzOTli7FnbvLq8erjmkFkQkdQA3ANOBY4EZksZlHXYlsMjM\nJgDnAdcDmNkTZjbJzCYBk4FtwKzoHAOui/eb2e/TugbXuspJqseaIbleTlddrFCQLOXfrX//sMrh\n88+XVw/XHNJsiUwBVprZajPbAdwBfCjrmHHAAxACBzBS0iFZx5wGPGVmmQ3jhDctOpdbJd/USxlP\nUQ/x+Jc0WiJJu7JinlxvfWkGkeFA5gf/s9G2TIuBswAkTQGOBLIbyh8Hbsva9vmoC+xmSQdXr8qu\nHezeDXPnlv8h269fsvEU9fL007DffqWNf8k0cSI89RRs2fLmfaW24Dy53vrSDCJJboK8GjhYUi9w\nMdAL7Ip3SuoLnAH834xzfgiMAiYC64Frq1Vh1x6WL4fBg+Gww8ovo5GT65V01QHsuy+ceCLMm7f3\n9i1bQoCaMCF5WZ5cb337pFj2WiDzu9AIQmvkDWa2Ffh0/FzSKuDpjEPeBzxsZi9knPN8xvE3Ab/N\nV4GZM2e+8bi7u5vu7u4SL8G1olK7ZHLp6oJ///fq1KfaKunKisVdWqefvmfbvHkhuJRyM4IHkcbW\n09NDT09PRWWkGUQWAqMljQTWAecAMzIPkHQQ8KqZbZd0IfCgmb2cccgM4Pasc4aZWbw00JnA0nwV\nyAwizsUq/aYO4UP6s58N+Yek04rUypw58KlPVVZGVxf88Id7byvn323ECHjkkcrq4tKT/eX6qquu\nKrmM1LqzzGwnoYvqXuBx4E4zWybpIkkXRYcdCyyVtBx4L3BpfL6k/QlJ9buyir5G0hJJi4F3Al9I\n6xpca6pGS6TYeIp6ice/TJxYWTldXaHlkXl7bjn/bp5Yb31ptkQws3uAe7K23ZjxeA4wNs+5rwBD\ncmw/r8rVdG1k06bQvXL88ZWXFXf5jB5deVnVMn9+GEBZzviXTEOHwpAhsGxZWLQrvhnh3/6ttHI8\nsd76fMS6aytz54Z1QfapwtenRkyuV6OrLpa5vsiyZXDIISG4lGL48LAu/a5dxY91zcmDiGsr1ejK\nijXiyPVqJNVjmUGy3HL32y/cCffcc9Wpk2s8HkRcW6nmN/WJE8Mtr7nGU9RDpeNfsmUHkXL/3fwO\nrdbmQcS1jZ07YcECmDq1OuXlG09RL/H4l0MPrU55xx0XuqI2bqysBefJ9dbmQcS1jUcfDX30gwdX\nr8xGmkerkvmycunoCPmj3/0O1q0LswKXw5Prrc2DiGsb1ezKijVScj2N6+vqgn/9Vzj55BBUyuHd\nWa3Ng4hrG9VMqsemTg15iEaY7ryaSfXYtGlhsGAl5XoQaW2pjhNxDuDFF8OHbKm3h1bb7Nnw5S9X\nt8yhQ8OtrzfdBMOGVbfsUrz+eli7o9wup3zi/FElLZxaB5EdO+CZZ+Doo2v3mu3Mg4hL3Te/GT7k\nfvSj+tVhw4Yw0PCYY6pf9qWXwm/zzuBWO5ddVp3xL5kOPhi+8hU45ZTyy6h1Yv03v4HvfKdxbnho\ndbJCiykDkh4GfgrcZmZ/qUmtqkCSFbs2VxsnnRSCyJIl9avDr38NN94I99xT/FhXXTt2wP77h/VX\nqh3kcvniF+H228OdZa40kjCzkmaDS5IT+ThhHZAFku6Q9F6p0aacc41q27awgNPTT8NLL9WvHmkk\nnV0y++4buvzWravN682eHVqe27fX5vXaXdEgYmZPmtmVwBjC4lA/Bf4s6SpJVbxZ0rWihQtDP/2J\nJ4Z5neoljaS6S65WeZHXXoOlS8O8X2vXpv96LuHdWZImANcB3wF+BXwM2Ar8Kb2quVYQtwDqeSvs\n9u3Q2xvGPLj6qFUQefhhGDcuTIrpAxxro2gPZZQTeQm4CbjczF6Pds2VVEG6zbWDOXPg3HPDHErZ\n61PUyiOPwFFHwYEH1uf1Xe2S63GL84UX/LbiWknSEvmYmb3bzG7LCCAAmNmZKdXLtQCzPS2Rrq76\njafwfEj91WrUevxe+9iU2kkSRP5W0sHxE0mDJH0rxTq5FrFyZVi4qbNzz3iKZctqXw8PIvVXiw91\nsz0tEQ8itZMkiLzfzDbHT6LbfP86vSq5VpGdzK7X1OmeVK+/Wnyor14dlio+8kifr6uWkgSRPpL6\nxU8k9QcqXDfNtYPsFkA9kutr1oQxKkcdVdvXdXurRRCJ/79JPnNwLSUJIj8H/ijpAkl/C/wBKHGR\nTNeOcrVEaj3jbVwHH9lUX4cdFmYMSHPsRub/N+/Oqp0k40SuAb4FHAscA3wj2uZcXlu2wFNPhYWb\nYuPHhwFnGzfWrh6eD2kMHR0hkKQ5diPzvR46NPwffO219F7PBYnGiZjZPWZ2mZl9yczuTbtSrvnN\nmxcGGPbN6Pjs6AhTis+dW7t6eBBpHGm2Dl5+GZ54IvyfA+jTBw4/3Lu0aqFoEJHUJWmBpJcl7ZC0\nW1KDLAjqGlW+ZHYtk+uvvgqPPQZve1ttXs8VlmYQWbAAJkwI45FinlyvjSQtkRuATwBPAv2AC4Af\npFkp1/zytQBqmVxfuDAs8dq/f21ezxWWZhDJ9f/Nk+u1kbQ760mgw8x2mdktwPR0q+Wa2e7doTsr\nV0vk5JPDh/vOnenXw2/tbSxpfqjneq89uV4bSYLIK5L2AxZL+rakLwJ+r4vLa9kyeMtbci9CNWgQ\nHHFEbaaF93xIY0mreylzkGEmDyK1kSSI/E103MXANqAT+EialXLNrdiH97Rp6d/qmznlimsMaX2o\nr1gBBxwQEum1eD23t4JBRNI+wP8ys1fN7CUzm2lmXzSzlTWqn2tCxbqRapFcf/rpcGfYiBHpvo5L\nLq0P9XxfFjyxXhsFg4iZ7QSOjLqznEskSUsk7SDirZDGE4/dePXV6pab7732xHptJFmschXwkKS7\nCd1ZAGZm16VXLdesNm4MAwrHj89/zJgxYZXD9eth2LB06jF7tifVG02fPjB8eBhwePTR1St3zhz4\nu7978/YhQ8LKmtu2wYAB1Xs9t7ckOZGngN9Fxw4EDoh+ipI0XdJySU9KujzH/kGSZklaLGmepOOi\n7WMl9Wb8vCTpkmjfYEn3S1oh6b7MGYZd/c2dGxZ/6ujIf0yfPulPgTJnjrdEGlG1u5g2b4ZnnoET\nTnjzPsm7tGqhaEvEzGaWU7CkDsIYk9OAtYQ12u82s8zJwK8EFpnZmZLGAv8bOM3MngAmReX0ic6f\nFZ1zBXC/mX07CkxXRD+uASTtRoqT62edVf06bN0apqGfNKn6ZbvKVDsvMm8eTJ4c1nEv9Hpjx1bv\nNd3ekoxYfyDHT5JlcacAK81stZntAO4APpR1zDjgAYAocIyUdEjWMacBT5lZ/F/vg8Ct0eNbgQ8n\nqIurkaRjM9JMrs+fH+bs6utzTTecageRYl9avCWSviQ5kX/MeNyPcHtvkqFiw4HMt+9Z4OSsYxYD\nZxFyLlOAIwm3EL+QcczHgdsynh9qZhuixxuAQxPUxdXAzp1h+ompU4sfO2UKLF4cpmnfr8q3bXhS\nvXGNGAFLl1avvDlz4JJLCr+eJ9fTlaQ7a2HWpockLUhQtiU45mrgekm9wFKgF9gV75TUFzgDeFM+\nJaqbScr7OjNnznzjcXd3N93d3Qmq5Mq1ZEn4ox00qPixAwfC6NHQ25ss6JRi9mz4zGeqW6arjhEj\n4J57qlPWrl2hO6vQ/58RI2DRouq8Xivq6emhp6enojKKBhFJgzOe9gHeBhyYoOy1QOZd+iMIrZE3\nmNlW4NMZr7UKeDrjkPcBD5tZZstkg6TDzOw5ScOA5/NVIDOIuPSVmsyOb/WtZhDZvTsk92+5pXpl\nuuqpZvfSY4+F6eWHDMl/zIgR8JvfVOf1WlH2l+urrrqq5DKSdGctYk+rYiewmjAJYzELgdGSRgLr\ngHOAGZkHSDoIeNXMtku6EHjQzF7OOGQGcHtWuXcD5wPXRL9/naAurgZmz4ZTT01+/LRp8Osqv3vL\nl4eW0GGHVbdcVx3VzIkkyb/5qPX0JVmUaqSZjYp+RpvZ6Wb2UILzdhKmSrkXeBy408yWSbpI0kXR\nYccCSyUtB94LXBqfL2l/QlL9rqyirwZOl7QCeHf03DWAUic8jJPrlqTjs4Q6eD6kcQ0ZEgYbvvJK\n5WUlyX15Yj19siJ/wZI+B9xmZn+Jng8CZphZQ08HL8mKXZurnvXrw7TrL74YxoEkYRYGG86fHyZl\nrIYLLgi3fH72s9Upz1Xf6NHw29/CMcdUVs6YMXDXXYUHtpqF/Nv69XBgkk74NicJMytpgt0kf+6f\niQMIQPTY05ZuL3PmhNxG0gACYTBYtW/19ZHqja8ad0y98AI8/zwce2zh4yS/QyttSf7k+0QD/oA3\nBhHmGdrj2lW53UjVnEdr06bwYXH88dUpz6WjGl1Mc+aEtWmSfGnxvEi6kgSRe4E7JJ0q6TTCoMHf\np1st12zKHZtRzWnh586Fk06CfZLcLuLqphof6qXk3zyIpCtJELmcMKr874G/A/4A/FOalXLN5fXX\n4ZFHwgDCUk2eDI8/HibJq5Qn1ZtDNT7US/nS4sn1dCUJIv2AH5vZR83so8BNgE8N797Q2xuSnAMH\nln5uv36h+2lh9pDWMvhI9eZQaRDZsSMMIDw5e/6LAq/nOZH0JAkifwL6ZzwfQGiNOAdUnsyuRnJ9\n585wl1e1R7+76qv0Q33xYhg5Eg46KPnreUskPUmCyH6ZAwCjUeY+O797Q6XdSNVIrj/6aOi2GDy4\n+LGuvirtXiq1xelBJF1JgsgrkibHTyS9Dajy2mSuWVVjLfN4bZFKhvX4rb3NY9Cg0HLcsqW880sd\n1BoHER82lo4kQeQfgF9IekjSQ8CdwOfTrZZrFmvWhA+EUaPKL6OzE/r3D2uAlMuT6s0jHrtRbuug\n1C8tBx4YbgXevLm813OFJZn2ZAFh3Y/47qxjcszs69pU/Aetksa4vlmlt/p6Ur25lBtE1q4NU6aM\nHl3663lyPR1JxxePJcxzNRmYIem89Krkmkm1upEqSa5v2BAGGlY6jYarnXI/1OOurFK/tHheJD1J\nVjacCXwP+D7QDXybsLqgc1XrRqqkJVLOlCuuvspNrpf7/82DSHqS/Nl9lDCb7noz+xQwATg41Vq5\npvDKK2Gg4OTJxY8tZuJEeOqp8pKtnlRvPuV+qJf7XnsQSU+SIPKqme0CdkbrfzzP3otNuTa1cGGY\nQbV//+LHFrPvvnDiiWGlulJ5PqT5lPOh/tprYfXMk04q/fV81Hp6kgSRBdH07z8hLDTVC1Rx3lXX\nrKp9R1Q5XVrbt4cR8+VMueLqp5wgsmgRjBsH++9f3ut5Yj0dSe7O+qyZ/cXMfgS8Bzg/6tZyba7a\n3UjlJNd7e+Hoo32tiGYTf6iXMnajkv9v3p2VnpJSkWa2yswWp1UZ1zzMqt8S6eoK3Vm7dyc/x8eH\nNKdyxm5U8l6XE7RcMn4/iyvLypVh8sTOzuqVOXQovOUtsGxZ8nM8qd68SmkdxDMjlPte779/+P+6\ncWN557v8PIi4sqSVzC51Hi1PqjevUoLI6tVhbMiRR5b/ep5cT0fe5XskFZzKzsw2Vb86raOnB669\ntt61SM/y5XDxxdUvd9o0uPpquPvu4sfu2hUS60cdVf16uPQdcQR85Svwox8VP/aFFyqfGSHu0po0\nqfwyrrgCvv716tyRmOm228Io/HLuPKu3QmvALQIMEHAEEK+zPgh4BqhgtqTW9/vfh66Zj3yk3jVJ\nT3d39cv8m78J3xiT9l13dlY+5Yqrj699LdwYkdSECZW9XqXJ9XXr4Jpr4P3vh7/6q8rqku3734cz\nzmixIGJmIwEk/QSYZWb/GT1/H3BmTWrXxNasgenTw38Ml9z++8MHPlDvWrha6Oysbk6tmEqDSHz7\n+ezZ1Q0ir70GDz8cBtw2oyQ5ka44gACY2T2A90IX8eyztf0Dcc4VVo0gMnFiZROF5rJoUVitsVnH\nsSQJIusk/Q9JIyWNkvQVYG3aFWt2a9aE/7TOucZQjcWwLrss/K7mrcKzZ8MppzRv0j9JEJkBDAVm\nAXdFj2ekWalmt3t3mLLaWyLONY5KRq2//npYlvfDH6587Ztss2fD2We3cBAxs41mdgnwdjObZGaX\n+p1ZhT3/fFj/uV+/etfEORfr7Axf7koZzBpbtAjGjoWBA/esxFkN8aDdM86AbdvCT7NJMhX8NEmP\nA8uj5xMk/SD1mjUx78pyrvH07w8HHBBuFy5V5nikUscyFRKPfxk5EoYPb87WSJLurO8C04EXAaJp\nT96ZZqWanSfVnWtM5SbXM6dcqXQVzlzlxksGN2NyPdGIdTP7c9amnUnOkzRd0nJJT0q6PMf+QZJm\nSVosaZ6k4zL2HSzpl5KWSXpc0snR9pmSnpXUG/1MT1KXWvKWiHONqZzkuhn893/vmXJlwoTy177J\nljmVS7NOEpkkiPxZ0ikAkvpK+hJQdHYjSR3ADYRWzLGEZXXHZR12JbDIzCYA5wHXZ+y7HvhPMxsH\nnEDUnUYYAHldlJ+ZZGa/T3ANNeVBxLnGVM63/WeeCb9Hjgy/+/Ytf+2bbJndZK0cRP4e+BwwnHBr\n76ToeTFTgJVmttrMdgB3AB/KOmYc8ACAmT0BjJR0SLT41TvM7KfRvp1m9lLGeQ09RtmDiHONqZwP\n6lzrulcjuf7yy/DEEyEglVu3RpAkiIwxs0+Y2VAzO8TMzgWOSXDecCDzn+TZaFumxcBZAJKmAEcC\nnYQpVV6QdIukRZJ+ImlAxnmfj7rAbpbUcEv1ehBxrjGV80Gda5LPaiTXFywIXWP77Reed3Y2Z06k\n0NxZsRsIrY9i27IlGY5zNXC9pF5gKWHVxF1AX+BE4GIzWyDpu8AVwNeAHwLfiM7/JnAtcEGuwmfO\nnPnG4+7ubrrTmOwpB0+sO9eYyg0iM7JGxnV1wSc/GW4X7lPmXOjZ66PUoyXS09NDT09PRWXI8gy9\nlNRFmN7kC8B17OlCOgA4M8pj5C9YmgrMNLPp0fMvA7vN7JoC56wCjgcGAnPMbFS0/e3AFWb2gazj\nRwK/NbPjc5Rl+a4tTbt2wYABsHVr6Dt1zjWOp5+Gd71rT56jmFdeCevcbNz45nFfRx8Nv/kNHHdc\n7nOL+cAH4FOf2jNJ66ZN8Na3lrZQV7VJwsxKShcUiqF9CQGjI/o9MPrZAnw0QdkLgdHRdCl9gXOA\nvSb4lnRQtA9JFwIPmtnLZvYcsEbSmOjQ04DHouOGZRRxJqEF0zCeew4GD/YA4lwjGj4c1q8PX/aS\nWLAATjgh98DhSrq04kGGmYtsDRoU5tDaurW8Muul0Cy+DwIPSvqZma0utWAz2ynpYuBeQiC62cyW\nSboo2n8j4a6tn0ky4FH27pb6PPDzKMg8BcTrul8jaSKhu2wVcFGpdUuT50Oca1z77Re+5G3YAIcf\nXvz47A/6THFy/cILS6/HihVh4GNmHeKxImvWwLHHll5mvSTJiWyT9C+ED/x4KRYzs3cXOzGa8fee\nrG03ZjyeA4zNc+5i4E2z65vZeQnqXDceRJxrbPEHdZIgMns2nH9+7n3TpsH11+fel6TcXCtyxsn1\nZgoiSVJCPyeM0XgrMBNYTeiqcjl4Ut25xpY0gR13OeVbfnn8+NA1Vs667fnKbcbbfJMEkbeY2U3A\ndjN70Mw+BRRthbQrb4k419iSflA/+WRYJC1fi6WjI6xEOHdu6XXIHKleTt0aSZIgsj36/ZykD0g6\nkbBErsvBg4hzjS3peIx8XU6Zykmub94c7g474YQ372vVIPI/owF9lwFfAm4i3PbrcvAg4lxjS/pB\nXSipHitn5Pq8eTB5Muy7b/l1ayRFE+tm9tvo4WagO9XatAAPIs41tqQf1LNnF7/zaurUcBvwzp2w\nT5LblCjcwmnGUet5L1vS9zOeGnsGGxpAtFCVy7BjR1irYNiw4sc65+ojSRB56SVYtSpMS1LIoEFw\nxBGwdClMKjaHR2TOHLgkz6dnXDezvefqamSFurMejn72I0xBsgJ4kjDdiQ+ly2H9+jC6Nek3Eudc\n7Q0bFlYf3VlgQYtCXU7ZurqS50V27QplT52ae/9BB4Xg8dJLufc3orxBxMx+ZmY/AyYA7zKz75vZ\n9wh3ZiWMue3Fu7Kca3z77guHHBK+9OWTJKkeKyW5/thjcNhhMGRI/mOaLS+SJLF+MHBgxvMDom0u\niwcR55pDsQ/qJEn1WCnJ9STltmIQuRpYJOlWSbcCi4B/TrdazcmDiHPNodAH9a5dYexH0iAydmy4\nbbdQyyaWpIXTbMn1okHEzG4BpgKzgLuAqVE3l8vio9Wdaw6Fgsjjj8Ohh4YuryT69EneGkkSRFqm\nJRIvZStpMjCMsMDUs8Dh0YBDl8VbIs41h0If1KV0ZcWSBJEXXgg/xebFarYgUug+oi8CFxIWfcq1\nMMe7UqmD8BarAAAWFklEQVRRE/Mg4lxz6OyEhx7Kva+UpHps2jT46lcLHzNnDpx8cvFFrFomiJjZ\nhdHv7prVpsl5EHGuORRriXzxi6WVN2UKPPIIvP76nuVuc5WbpIXTMkFE0kcosMStmd2VSo2a1Ouv\nh5XJDj203jVxzhWT74P6xRfDwnKlrlY4cCCMGQO9vfnHgMyeDV/5SvGy4sR6sww4LNSddQaF10n3\nIJJh3bowiKmjo941cc4Vc9hh4Uvf9u17r0IadzmV83ccjxfJFUR27ICHHw5lFzNwYGjNbNoEb3lL\n6fWotULdWZ+sYT2anndlOdc8OjpCIFm7FkaN2rO9nKR6rKsrrLmey+LF4XUOOihZWXFLqRmCSJJx\nIkRTwP+TpK/FP2lXrNl4EHGuueQaj1FOUj0Wt0QsR/9NqeU2U16kaBCRdCNwNnAJYRLGs4EjU65X\n0/Eg4lxzyf6gLqXLKZdRo8J8XLk+/Ett4bRUEAGmReuabzKzqwgDD3Oui97OPIg411yyP6iXLIEj\nj4SDy5zUSco/j1apLZFmGrWeJIi8Gv3eJmk4sBM4LL0qNScfre5cc8kOIpV0ZcVyBZG1a+GVV2D0\n6PLr1siSBJH/kDQI+A5havjVwO1pVqoZeUvEueaS/UFdSVI9lmvkelxuKbfrtlQQMbNvmNlfzOxX\nwEjgGDMrMjaz/XgQca65ZHcZVaMlMnlymHtr27bKym2pICJpiaQrJR1lZq+Z2eZaVKyZvPoqbN2a\nfMI251z9ZX5Qr1sX/obHjKmszP79Yfx4WLhwz7ZyWjidnaEbbPfuyupTC0m6sz4I7AJ+IWmhpC9J\nOiLlejWVZ5+F4cOLz4njnGscQ4fCli3w2mvldTnlk5kXee21kLA/6aTSyujfPww6fPHFyuuTtiTd\nWavN7BozmwzMAE4AVqVesybiSXXnmk+fPnD44eHvtxpdWbHMIPLwwzBuHOy/f+nlNEuXVtLBhiMl\nXQ7cARwD/FOqtWoyng9xrjnFH9TVSKrH4uS6WWXlNksQKTR3FgCS5gF9gV8AHzOzp1OvVZPxIOJc\nc+rshKeeCtOSlNrlVKjM/v1h5crQIvnYx8orp1mCSJKWyPlmNsnM/rnUACJpuqTlkp6MWjLZ+wdJ\nmiVpsaR5ko7L2HewpF9KWibpcUlTo+2DJd0vaYWk+yTVfb13DyLONacRI8J8V2PHhhxEtXR1hQDS\nDi2RJDmR5eUULKkDuAGYDhwLzIhXS8xwJbDIzCYA5wHXZ+y7HvhPMxtHyMMsi7ZfAdxvZmOAP0bP\n68qDiHPNacQIuPfe6uVDYtOmwW23hUT9kWVOEtUso9bTvJ9oCrAySszvIORTPpR1zDjgAQAzewIY\nKekQSQcB7zCzn0b7dprZS9E5HwRujR7fCnw4xWtIxBPrzjWnESPCnFlpBJH77gu/y73jq2VaIhUY\nTliXPfZstC3TYuAsAElTCBM7dgKjgBck3SJpkaSfSBoQnXOomW2IHm8AKloGatOmsKBUJbwl4lxz\niv9uq5VUj02YAP36VVZuswSRJIn1s4Hfm9kWSV8FTgS+aWaLipxaaEGr2NXA9ZJ6gaVAL2FMSt/o\ndS42swWSvkvottprCnozM0l5X2fmzJlvPO7u7qa7u/tNx3z4w2Ft5NNPT1DbHF55JdwL3gzz/jvn\n9vbWt8J73gMjR1a33L594ZxzQtnlGj4c1q+HXbvSW+yup6eHnp6eisqQ5Zr8PvMAaamZHS/p7cC3\ngH8BvmpmBSdMjhLhM81sevT8y8BuM7umwDmrgOOBgcAcMxsVbX8HcLmZfUDScqDbzJ6TNAx4wMyO\nyVGWFbs2gMsvhwED4OtfL3poTsuXwwc/CCtWlHe+c87lc9hhsGhRGM9SC5Iws5I64JJ0Z+2Kfn8A\n+ImZ/QehpVDMQmB0NMakL3AOcHdWhQ+K9iHpQuBBM3vZzJ4D1kiKJyE4FXgsenw3cH70+Hzg1wnq\nkte0aW+eMK0U3pXlnEtLMyTXkwSRtZJ+TAgCv5PUL8l5ZrYTuBi4F3gcuNPMlkm6SNJF0WHHAkuj\n1sV7gUszivg88HNJiwl3Z/2vaPvVwOmSVgDvjp6XrasL5s4tf44aT6o759LSDHmRojkRwkqG04Hv\nmNnmqAvpH5MUbmb3APdkbbsx4/Ec8ixwZWaLgTcN/zGzTcBpSV4/iaFDYcgQWLYMjjuu+PHZvCXi\nnEtLMwSRJC2Rw4DfmdmTkt5FCCrz061WbeVbjSwJDyLOubS0ShC5C9gp6WjgRsItuLelWqsai0eX\nlsODiHMuLa0SRHZH+Y2zgO+b2T8Cw9KtVm1Vklz3IOKcS0urJNa3S/oEYVqS/4i27ZtelWpv/Piw\nKM3GjaWf64l151xaWqUl8mmgC/ifZrZK0luBf0+3WrXV0QFTpoS7tEqxZUsYCHRw3aeAdM61osMP\nhw0bYOfOetckvyS36j4GfAl4VNJ4YE2hAYPNqpzketyVVY3V0JxzLtu++4Zlt9evr3dN8kuyxno3\nsAL438APgCclvTPletVcOcl1z4c459LW6F1aScaJXAe8J5pll2gU+R2Eua1axtSpsHBhaDbuk+Rf\nBQ8izrn0NXpyPUlOZJ84gACY2QqSBZ+mMmgQHHEELFmS/BxPqjvn0tboLZEkQeRhSTdJ6pb0Lkk3\nEebFajnx2shJeUvEOZe2Vggif0dYVfASwnxWjwF/n2al6qXU5LoHEedc2ho9iBTslpK0D7A4mmr9\n2tpUqX66uuBb30p+vAcR51zaGj2IFGyJRCPVn5BU5irBzWXsWNi8OdntdGYeRJxz6WuFxPpg4DFJ\nf5L02+jn7qJnNaE+fZLnRTZvDndxHXBA+vVyzrWvYcPgxRdh+/Z61yS3JHdZfTXreZJlb5tWHETO\nOqvwcd4Kcc7VQkdHWOFw3brqL+NbDXlbIpJGS3q7mfVk/hBWOmzgxlVlkibXPYg452qlkfMihbqz\nvgtsybF9S7SvJU2ZAo88Aq+/Xvg4DyLOuVpp1iByqJm9aehdtG1UelWqr4EDYcwY6O0tfJwHEedc\nrTRycr1QECk0N22/alekkSTp0vLR6s65WmnWlshCSZ/J3ijpQuDh9KpUf0nu0PKWiHOuVho5iBS6\nO+sfgFmSzmVP0JgM7AecmXbF6mnaNLj88jAWJN807x5EnHO10shBRGb579iVJOBdwHjCrb2Pmdmf\nalS3ikiyQtdWiFm4pW7BgjApY679AwaElRAHDKiwos45V8Rzz8EJJ8Dzz6f7OpIws5JWSCo4TiT6\nFP5T9NM2pD15kVxB5MUXQ/DwAOKcq4WhQ+Gll+C116Bfg2Wkk4xYb0uFkuueVHfO1VKfPmGp3LVr\n612TN/Mgkkeh5LrnQ5xztdaoeREPInlMngyPPw7btr15nwcR51yteRBpMv37w/jxYcncbB5EnHO1\n5kGkCeXLi3gQcc7V2ogRjTlqPdUgImm6pOWSnpR0eY79gyTNkrRY0jxJx2XsWy1piaReSfMzts+U\n9Gy0vVfS9LTqny+IeGLdOVdrnZ1t1hKR1AHcAEwHjgVmSBqXddiVwCIzmwCcB1yfsc+AbjObZGZT\nsrZfF22fZGa/T+sa4uR69nATb4k452qtHbuzpgArzWy1me0A7gA+lHXMOOABADN7Ahgp6ZCM/fkG\nvZQ0GKZcnZ0hN7Jy5Z5tu3eH2+y8JeKcq6V2DCLDgcxLfjbalmkxcBaApCnAkUD88WzAHyQtjObr\nyvT5qAvsZkmFJoqsWPatvs8/Dwcd1HgDfpxzrW3IkHC3aK47RuspzSCSZM6Rq4GDJfUCFwO9hEWv\nAN5uZpOA9wGfk/SOaPsPCVPRTwTWA9dWtdZZsvMi3pXlnKsHqTGnhE+yPG651gKZH7cjyFoR0cy2\nAp+On0taBTwd7VsX/X5B0ixC99h/mdnzGcffBPw2XwVmzpz5xuPu7m66u7tLvohp0+Dmm/c896S6\nc65e4uT6mDHVKa+np4eenp6Kyig4AWNFBUv7AE8ApwLrgPnADDNblnHMQcCrZrY96rI6xcw+KWkA\n0GFmWyXtD9wHXGVm90kaZmbro/O/AJxkZp/I8fplT8CYaft2GDw4rG984IHwve/BihVwww0VF+2c\ncyU57zx497vhk59Mp/yqT8BYCTPbKeli4F6gA7jZzJZJuijafyPhrq2fSTLgUeCC6PRDCdPQx3X8\nuZndF+27RtJEQnfZKuCitK4BoG9fOPFEmDcPTj/du7Occ/XTiMn1NLuzMLN7gHuytt2Y8XgOMDbH\neasIOY9cZZ5X5WoWFSfX4yAyaVKta+CccyGILFpU71rszUesJ5CZXPeWiHOuXhpx1LoHkQS6umDu\n3DBGxBPrzrl6acRR6x5EEhg6NNyj/eijYYWx4dmjXZxzrgYaMSfiQSShadNg1qxwp1bfvvWujXOu\nHQ0aBDt2wNat9a7JHh5EEurqgjvv9HyIc65+pMZrjXgQSWjaNFi2zIOIc66+Gi257kEkofHjYeBA\nT6o75+qr0ZLrqY4TaSUdHXDyyd4Scc7V14gRcNttjRNIPIiU4JvfhGHD6l0L51w7+8QnoE8D9SGl\nNndWvVVr7iznnGsX5cyd1UDxzDnnXLPxIOKcc65sHkScc86VzYOIc865snkQcc45VzYPIs4558rm\nQcQ551zZPIg455wrmwcR55xzZfMg4pxzrmweRJxzzpXNg4hzzrmyeRBxzjlXNg8izjnnyuZBxDnn\nXNk8iDjnnCubBxHnnHNl8yDinHOubKkGEUnTJS2X9KSky3PsHyRplqTFkuZJOi5j32pJSyT1Spqf\nsX2wpPslrZB0n6SD07wG55xz+aUWRCR1ADcA04FjgRmSxmUddiWwyMwmAOcB12fsM6DbzCaZ2ZSM\n7VcA95vZGOCP0fO20tPTU+8qpMqvr7n59bWXNFsiU4CVZrbazHYAdwAfyjpmHPAAgJk9AYyUdEjG\n/lwLxn8QuDV6fCvw4arWugm0+n9iv77m5tfXXtIMIsOBNRnPn422ZVoMnAUgaQpwJNAZ7TPgD5IW\nSrow45xDzWxD9HgDcGi1K+6ccy6ZfVIs2xIcczVwvaReYCnQC+yK9r3dzNZFLZP7JS03s//a6wXM\nTFKS13HOOZcCmaXzGSxpKjDTzKZHz78M7Dazawqcswo43sxeztr+dWCrmV0naTkhV/KcpGHAA2Z2\nTI6yPLg451yJzCxXGiGvNFsiC4HRkkYC64BzgBmZB0g6CHjVzLZHXVYPmtnLkgYAHWa2VdL+wHuA\nq6LT7gbOB66Jfv8614uX+g/hnHOudKkFETPbKeli4F6gA7jZzJZJuijafyPhrq2fRa2GR4ELotMP\nBWZJiuv4czO7L9p3NfALSRcAq4Gz07oG55xzhaXWneWcc671tdyI9WIDHJtdvkGYzUrSTyVtkLQ0\nY1vLDCjNc30zJT0bvYe9kqbXs47lkjRC0gOSHpP0qKRLou0t8f4VuL5Wef/6RYO8H4mub2a0vaT3\nr6VaItEAxyeA04C1wAJghpktq2vFqii6+WCymW2qd12qQdI7gJeBfzOz46Nt3wZeNLNvR18EBplZ\nUw4qzXN9b9woUtfKVUjSYcBhZvaIpIHAw4RxW5+iBd6/Atd3Ni3w/gFIGmBm2yTtAzwEXAp8hBLe\nv1ZriSQZ4NgKWuamgei27b9kbW6ZAaV5rg9a4D00s+fM7JHo8cvAMsJYsJZ4/wpcH7TA+wdgZtui\nh32BfQlDM0p6/1otiCQZ4Njs8g3CbCXtMKD089GccTc3a3dPpuguzEnAPFrw/cu4vrnRppZ4/yT1\nkfQI4X26z8zmU+L712pBpHX65vI7xcwmAe8DPhd1l7QsC/2trfa+/hAYBUwE1gPX1rc6lYm6en4F\nXGpmWzP3tcL7F13fLwnX9zIt9P6Z2W4zm0iYKeRkSeOz9hd9/1otiKwFRmQ8H0FojbQMM1sf/X4B\nmEXowms1G6L+aKIBpc/XuT5VZWbPWwS4iSZ+DyXtSwgg/25m8Zitlnn/Mq7v/8TX10rvX8zMXiLM\nY/heSnz/Wi2IvDHAUVJfwgDHu+tcp6qRNEDSAdHjeBDm0sJnNaV4QCkUGFDarKI/zNiZNOl7qDCQ\n62bgcTP7bsaulnj/8l1fC71/Q+KuOEn9gdMJeZ+S3r+WujsLQNL7gO+yZ4DjP9e5SlUjaRSh9QF7\nBmE29fVJuh14JzCE0P/6NeA3wC+AI4gGlJrZ5nrVsRI5ru/rQDehK8SAVcBFGX3QTUPS24H/Byxh\nT5fHl4H5tMD7l+f6riTMvNEK79/xhMR5B6FBcaeZfUvSYEp4/1ouiDjnnKudVuvOcs45V0MeRJxz\nzpXNg4hzzrmyeRBxzjlXNg8izjnnyuZBxDnnXNk8iLimJulPkt6Tte0fJP2gwDk9kianXK/bo7mV\nLs3aPlPSZdHjftGU21/Lcf7HJD0u6Y8V1OHljMfvl/SEpCOiOrwi6ZA8x+6W9C8Zz78UzTzs3Jt4\nEHHN7nbg41nbzgFuK3BOqvM5RVNGvM3MJpjZ9bleO5pR4VfAAjP7Ro5iLgD+1sxOTfiauVYptWjf\nqcD1wHQz+3O070XgsuxjI9uBMyW9Jcc+5/biQcQ1u18Bfx1/iEazrR5uZg9J+qGkBZkL7mTL+gb+\nUUm3RI8PkfRLSfOjn2k5zu0n6RaFRcIWSeqOdt0HDI8WLHp7jpfdl7BMwRNmdmWOcr8GnAL8VNI1\nkvbL9TqSPinp7qi1cn+e6/sr4MfAX5vZqmizAT8FzskzA+2O6Jwv5CrTuUweRFxTixbnmg+8P9r0\nceDO6PGVZnYSMAF4ZzTNw5uKyPP4euBfzWwK8FHCRHvZPgfsMrMTCFNh3Bq1MM4AnjKzSWb2UNY5\nAv4JeN3Mvpjnmr5BmAfuE2Z2OXBxjtfZLzp8EvARM3tXjqL6EabJ+ZCZrcja9zIhkPxDrjoAPwDO\nlXRgnv3OAR5EXGvI7NI6J3oO4Zv2w8Ai4DhgXAllngbcIKmXMJfXAZIGZB1zCvB/AMzsCeAZYAyF\nFywywgpy0ySNTliXfK9jwP0F5jXaDvw38Ld56vE94PxoqvO9d4Yp3f8NuCRhHV2b8iDiWsHdwKmS\nJgEDzKw3mqzyMuDdZjYB+B3hm3m2zNZH/4zHAk6OWhOTzGxExipwZB1Xqv9H6Cq6J55yO4F8r/NK\ngXN2E5ZynSLpy9nlRdN/30Zo6eTyXUJuZv+EdXRtyIOIa3rRQkEPALewJ6F+IOEDdoukQwmLeOWy\nQdIxkvoQpvWOg8p9ZHwLlzQxx7n/BZwb7R9DmPX0iYR1vgv4F+D3kg4qcniu11lOggBmZq8Bf03o\nmvp0jkOuAy4izAqdfe5fCLO5XoAn110eHkRcq7gdOD76jZktBnoJH7Y/J3Qh5XIF8B+Ebp91Gdsv\nAd4W3ab7GPCZHOf+AOgjaQkhUX6+me2I9hX60LWojj8i5Czuzshx5JLvdYrdZRa/zl+A6cD/kHRG\n1r6NwF2ENbb3Oi9yLWEae+dy8qngnXPOlc1bIs4558rmQcQ551zZPIg455wrmwcR55xzZfMg4pxz\nrmweRJxzzpXNg4hzzrmyeRBxzjlXtv8PoX6zsOaCb/4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x6dc5310>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(k_range, k_scores)\n",
    "plt.xlabel(\"Value of K for KNN\")\n",
    "plt.ylabel(\"Cross validated accuracy\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "上面的例子显示了偏置-方差的折中，K较小的情况时偏置较低，方差较高；K较高的情况时，偏置较高，方差较低；最佳的模型参数取在中间位置，该情况下，使得偏置和方差得以平衡，模型针对于非样本数据的泛化能力是最佳的。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.2 用于模型选择\n",
    "交叉验证也可以帮助我们进行模型选择，以下是一组例子，分别使用iris数据，KNN和logistic回归模型进行模型的比较和选择。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.98\n"
     ]
    }
   ],
   "source": [
    "# 10-fold cross-validation with the best KNN model\n",
    "knn = KNeighborsClassifier(n_neighbors=20)\n",
    "print cross_val_score(knn, X, y, cv=10, scoring='accuracy').mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.953333333333\n"
     ]
    }
   ],
   "source": [
    "# 10-fold cross-validation with logistic regression\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "logreg = LogisticRegression()\n",
    "print cross_val_score(logreg, X, y, cv=10, scoring='accuracy').mean()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.3 用于特征选择\n",
    "下面我们使用advertising数据，通过交叉验证来进行特征的选择，对比不同的特征组合对于模型的预测效果。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "from sklearn.linear_model import LinearRegression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# read in the advertising dataset\n",
    "data = pd.read_csv('http://www-bcf.usc.edu/~gareth/ISL/Advertising.csv', index_col=0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# create a Python list of three feature names\n",
    "feature_cols = ['TV', 'Radio', 'Newspaper']\n",
    "\n",
    "# use the list to select a subset of the DataFrame (X)\n",
    "X = data[feature_cols]\n",
    "\n",
    "# select the Sales column as the response (y)\n",
    "y = data.Sales"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[-3.56038438 -3.29767522 -2.08943356 -2.82474283 -1.3027754  -1.74163618\n",
      " -8.17338214 -2.11409746 -3.04273109 -2.45281793]\n"
     ]
    }
   ],
   "source": [
    "# 10-fold cv with all features\n",
    "lm = LinearRegression()\n",
    "scores = cross_val_score(lm, X, y, cv=10, scoring='mean_squared_error')\n",
    "print scores"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这里要注意的是，上面的scores都是负数，为什么均方误差会出现负数的情况呢？因为这里的mean_squared_error是一种损失函数，优化的目标的使其最小化，而分类准确率是一种奖励函数，优化的目标是使其最大化。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 3.56038438  3.29767522  2.08943356  2.82474283  1.3027754   1.74163618\n",
      "  8.17338214  2.11409746  3.04273109  2.45281793]\n"
     ]
    }
   ],
   "source": [
    "# fix the sign of MSE scores\n",
    "mse_scores = -scores\n",
    "print mse_scores"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 1.88689808  1.81595022  1.44548731  1.68069713  1.14139187  1.31971064\n",
      "  2.85891276  1.45399362  1.7443426   1.56614748]\n"
     ]
    }
   ],
   "source": [
    "# convert from MSE to RMSE\n",
    "rmse_scores = np.sqrt(mse_scores)\n",
    "print rmse_scores"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.69135317081\n"
     ]
    }
   ],
   "source": [
    "# calculate the average RMSE\n",
    "print rmse_scores.mean()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.67967484191\n"
     ]
    }
   ],
   "source": [
    "# 10-fold cross-validation with two features (excluding Newspaper)\n",
    "feature_cols = ['TV', 'Radio']\n",
    "X = data[feature_cols]\n",
    "print np.sqrt(-cross_val_score(lm, X, y, cv=10, scoring='mean_squared_error')).mean()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "由于不加入Newspaper这一个特征得到的分数较小（1.68 < 1.69），所以，使用所有特征得到的模型是一个更好的模型。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 参考资料\n",
    "- scikit-learn documentation: [Cross-validation, Model evaluation](http://scikit-learn.org/stable/modules/model_evaluation.html)\n",
    "- scikit-learn issue on GitHub: [MSE is negative when returned by cross_val_score](https://github.com/scikit-learn/scikit-learn/issues/2439)\n",
    "- Scott Fortmann-Roe: [Accurately Measuring Model Prediction Error](http://scott.fortmann-roe.com/docs/MeasuringError.html)\n",
    "- Harvard CS109: [Cross-Validation: The Right and Wrong Way](http://nbviewer.ipython.org/github/cs109/content/blob/master/lec_10_cross_val.ipynb)\n",
    "- Journal of Cheminformatics: [Cross-validation pitfalls when selecting and assessing regression and classification models](http://www.jcheminf.com/content/pdf/1758-2946-6-10.pdf)"
   ]
  }
 ],
 "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.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}