\n", "\t\t井出 剛著の「入門機械学習による異常検出」(以降、井出本と記す)の例題をSageを使ってお復習いします。\n", "\t
\n", "\n", "\n", "\n", "\t\n", "\t\tこの章でのポイントは、線形回帰モデルを使って、観測量(y', x')の負の対数尤度を異常度とする部分です。\n", "\t
\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\t\n", "\t\t線形回帰では、関数fを以下のような1次関数で表現します。\n", "$$\n", "\tf(x) = \\alpha_0 + \\alpha^T x = \\alpha_0 + \\alpha_1 x_1 + \\cdots + \\alpha_M x_M\n", "$$\t\t\n", "\t
\n", "\t\n", "\t\t係数αは、M+1個をもち、これをデータ$\\mathcal{D}$から以下の確立モデルを使って求めます。\n", "$$\n", "\\begin{eqnarray}\n", "\n", "\t\tp(y | x) & = & \\mathcal{N} (y | \\alpha_0 + \\alpha^T x, \\sigma^2) \\\\\n", "\t\t & = & \\frac{1}{\\sqrt{2 \\pi \\sigma^2}} exp \\left [ - \\frac{1}{2\\sigma^2}(y - \\alpha_0 - \\alpha^T x)^2 \\right ]\n", "\\end{eqnarray}\t\t\n", "$$\t\t\n", "
\n", "\n", "未知のパラメータ$\\alpha_0, \\alpha, \\sigma^2$ に対する尤度は、以下の様に定義されます。\n", "\n", "$$\n", "\tp(\\mathcal{D} | \\alpha_0, \\alpha, \\sigma^2) = \\prod_{n=1}^N \\mathcal{N}(y^{(n)} | \\alpha_0, \\alpha x^{(n)}, \\sigma^2)\n", "$$\n", "\n", "これの対数尤度は、以下の様になります。\n", "$$\n", "\tL(\\alpha_0, \\alpha, \\sigma^2 | \\mathcal{D}) = - \\frac{N}{2} ln(2\\pi\\sigma^2) - \\frac{1}{2\\sigma^2} \\sum_{n-=1}^N \\left [ y^{(n)} - \\alpha_0 - \\alpha^T x^{(n)} \\right ]^2\n", "$$\n", "
\n", "\n", "最初に$\\alpha_0$で微分して、その値が0となる$\\hat{\\alpha_0}$を求めると、以下の様になります。\n", "$$\n", "\t\\hat{\\alpha_0} = \\frac{1}{N} \\left [ y^{(n)} - \\alpha^T x^{(n)} \\right ] = \\bar{y} - \\alpha^T \\bar{x}\n", "$$\n", "\n", "ここで、$\\bar{x}, \\bar{y}$は、変数xと観測値yの平均です。\n", "$$\n", "\t\\bar{y} = \\frac{1}{N} \\sum_{n=1}^N y^{(n)}, \\bar{x} = \\frac{1}{N} \\sum_{n=1}^N x^{(n)}\n", "$$\n", "
\n", "\n", "求まった$\\alpha_0$を定数尤度の式に入れ、$\\alpha$で微分し0（勾配0）とすると、\n", "\n", "$$\n", "\t0 = \\frac{1}{2} \\frac{\\partial}{\\partial \\alpha} \\sum_{n=1}^N \\left [ y^{(n)} - \\bar{y} - \\alpha (x^{(n)} - \\bar{x}) \\right]^2\n", "$$\n", "\n", "ここで、変数xと観測値yと平均の差（中心化）を$\\tilde{X}, \\tilde{y}$で表すと、\n", "$$\n", "\t0 = \\frac{1}{2} \\frac{\\partial}{\\partial \\alpha} \\sum_{n=1}^N \\left [ \\tilde{y} - \\alpha \\tilde{X} \\right]^2 = \\sum_{n=1}^N \\left \\{ \\tilde{y} - \\alpha^T \\tilde{X} \\right \\} \\tilde{X}^T\n", "$$\n", "\t整理すると、\t\n", "$$\n", "\t0 = \\sum_{n=1}^N \\tilde{y} \\tilde{x}^T - \\alpha^T \\left ( \\sum_{n=1}^N \\tilde{X} \\tilde{X}^T \\right )\n", "$$\n", "\t\tこれを$\\alpha$について解くと\t、以下の様に求まります。\n", "$$\n", "\t\\hat{\\alpha} = \\left ( \\tilde{X}^T \\tilde{X} \\right )^{-1} \\tilde{X}^T \\tilde{y}\n", "$$\n", "
\n", "\n", "\n", "同様に$\\sigma^{-2}$で微分すると、\n", "$$\n", "\t\\sigma^2 = \\frac{1}{N} \\sum_{n=1}^N \\left \\{ \\tilde{y} - \\hat{\\alpha}^T \\tilde{X} \\right \\} ^2\n", "$$\n", "
\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\t\n", "\t\t観測量$(y', x')$についての異常度は、この観測量に対する負の対数尤度から、以下の様になります。\n", "$$\n", "\ta(y', x') = \\frac{1}{\\hat{\\sigma^2}} \\left [ y' - \\bar{y} - \\hat{\\alpha}^T (x' - \\bar{x}) \\right ] ^2\n", "$$\t\t\n", "\t
\n", "\n", "\n", "\n", "\t\n", "\t\t実際の観測量では、すべて変数xが独立ではなく互いに関連があるため、行列のランクが不足して逆行列が求まりません。\n", "\t
\n", "\t\n", "\t\tそこで、井出本ではリッジ回帰を使って$\\alpha, \\sigma^2$を求めています。リッジ回帰での$\\hat{\\alpha}_{ridge}$は、\n", "\t\t以下の様に求まります。\n", "$$\n", "\t\\hat{\\alpha}_{ridge} = \\left [ \\tilde{X}^T \\tilde{X} + \\lambda I_M \\right ]^{-1} \\tilde{X}^T \\tilde{y}\n", "$$\t\t\n", "\t
\n", "\t\n", "\t\t分散$\\hat{\\sigma}^2_{ridge}$は、以下の様になります。\n", "$$\n", "\t\\hat{\\sigma}^2_{ridge} = \\frac{1}{N} \\left \\{ \\hat{\\lambda} \\hat{\\alpha}^T_{ridge} \\hat{\\alpha}_{ridge} + \\sum_{n=1}^N \\left [ \\tilde{y}^{(n)} - \\hat{\\alpha}^T_{ridge} \\tilde{x}^{(n)} \\right ] \\right \\}\n", "$$\t\t\n", "\t
\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\t\n", "\t\tいつものように必要なライブラリを読み込み、テストデータとしてRのMASSパッケージに含まれているUScribeを使用します。\n", "\t
\n", "\t\n", "\t\tUScribeは、都市の犯罪率yで、以下の表6.1のような変数を使って回帰分析します。\n", "\t
\n", "" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "" ], "text/plain": [ "| 記号 | 定義 |
|---|---|
| M | 14-24歳の男性の割合 |
| Ed | 平均就学期間 |
| Po1 | 1960年における警察予算 |
| Po2 | 1959年における警察予算 |
| LF | 労働力率 |
| M.F | 女性1000人当たりの男性の数 |
| Pop | 州の人口 |
| NW | 1000人当たりの非白人数 |
| U1 | 都市部男性(14-24歳)の失業率 |
| U2 | 都市部男性(35-39)の失業率 |
| GDP | 州の1人当たりGDP |
| Ineq | 経済的不平等の度合い |
| Prob | 収監率 |
| Time | 刑務所での平均収監期間 |
\n", "\t\tリッジ回帰は、sklearnのlinear_modelパッケージのRidgeを使用しました。\n", "\t\tsklearnのパッケージはRと比べて使い方が統一されています。\n", "\t
\n", "\t\n", "\t\t回帰分析の結果求まった予測値(pred)と観測値(y)がどの程度分散しているか\n", "\t\tプロットしてみます。\n", "\t
\n", "" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": true }, "outputs": [], "source": [ "# データを変数にセット\n", "X = uscrime[vars].values\n", "y = uscrime['y'].values\n", "N = len(y)\n", "M = len(vars)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEDCAYAAADayhiNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3X14U+XdB/DvOUlfktKmKaRitcLqSztZXwlWnDKUCQ6h\n+IBUCzLoFC2IwkQeHO7hEiYXImwgVOeAFVRgOmWIrUBRQOXdksIcojLxQpBRoSk20PQt59zPH5HD\nAm3pS04T2u/nurxKcifNnZ9pvufc9zn3kYQQAkRERADkQHeAiIiCB0OBiIg0DAUiItIwFIiISMNQ\nICIiDUOBiIg0DAUiItIwFIiISMNQICIiDUOBiIg0LQ6FY8eOYfjw4ejWrRuuvvpq5ObmwuVyAQC2\nbt2KzMxMWCwWJCcnY82aNT7PXbx4MZKSkhAdHY1+/fqhtLRUa6utrUVeXh7i4+MRGxuL7OxsVFRU\nNLtfQgi4XC5w1Q4iotaTWrr2UWpqKvr06YMlS5bgzJkzuO+++5CWlobZs2fjxhtvRH5+PnJycrB9\n+3ZkZWVh+/btyMjIQGFhIcaNG4fi4mIkJyfjpZdewqJFi3DkyBGYTCZMnToVO3bswLp162A2m/HI\nI4+gvr4e69evb1a/XC4XLBYLKisrERUV1aIinD59tkWP9zdZlhATE4GKiiqoKkNND6yx/lhjfTVW\nX5st0r+v05IHV1ZWok+fPpg7dy5MJhPi4uIwduxYfPLJJ1i9ejUSExMxduxYhIaGYsCAAcjKysLy\n5csBAEuXLkVubi7sdjvCwsIwbdo0SJKEwsJCKIqCgoICzJw5E3FxcYiOjsacOXNQVFSEsrIyv77h\nYCTLEiRJgixLge5Kh8Ua64811kd9PbBxoxEffGCEqupf3xaFgsViwfLly2Gz2bT7jh8/jmuuuQYO\nhwMZGRk+j8/IyEBJSQkAXNIuSRLS0tJQUlKCI0eOoLKyEunp6Vp7YmIiTCYTHA5Hq94YEdGVTlGA\nUaNMGDvWhJyccIwcCeg9Qt6mieZ9+/YhPz8fzz77LJxOJ6xWq097TEwMysvLAaDJdqfTCUmSLmm3\nWq3a84mIOpuvvpLx8cdG7fa6dcCJE/ruKRgv/5CG7dy5E1lZWZg3bx7uuusuzJs377KTvG1t14ss\nB3aX12CQfX6S/7HG+mON/a9rVwmyLKCq3u+n0FAgKkqC0ahfjVsVCoWFhRgzZgxefvlljB49GgBg\ns9ngdDp9Hud0OhEbG9tke3JyMmw2G4QQcDqdMJvNWntFRYX2fD3FxERAkgI/DhoVZQp0Fzo81lh/\nrLH/WK3Aq68CU6cCBgPwyitAjx761rfFobBr1y6MGzcOa9euxYABA7T77XY7Vq5c6fPYkpISZGZm\nau0OhwNjxowBAKiqitLSUowfPx4JCQmwWq1wOByIj48HABw8eBB1dXWw2+2tfW/NVlFRFfA9hago\nE1yuaiiKGrB+dGSssf5YY33cf7/3v8bqa7VG+PX1WhQKiqJg/PjxmDdvnk8gAMDo0aPx3HPPoaCg\nAKNHj8aWLVuwceNG7N27FwAwYcIE5OTkICcnBykpKZg/fz7Cw8MxePBgyLKMRx99FHPmzIHdbofJ\nZMKMGTMwYsQIn0ltvaiqCIpD6BRFhcfDPyY9scb6Y431pXd9WzQwtXv3bnz55Zd48sknYTKZYDab\ntZ81NTUoKirCkiVLEB0djalTp2L16tXo1asXAGDQoEGYO3cusrOz0bVrV2zZsgUbNmxAWFgYAGD2\n7Nm49dZbkZqaiuuvvx4WiwXLli3z/zsmIqJGtfjktWB1JZ+8ZjTKsFojcOZMFbewdMIa64811ldj\n9Q3oyWtERNSxMRSIiEjDUCAiIg1DgYiINAwFIiLSMBSIiEjDUCAiIg1DgYiINAwFIiLSMBSIiEjD\nUCAiIg1DgYiINAwFIiLSMBSIiEjDUCAiIg1DgYiINC2+RjMREenn889lFBUZce21Ajk59ZDbedOd\noUBEFCQOH5Zx771muN0SAOCLL2Q8/3xtu/aBw0dEREFi2zaDFggAsGFD+2+3MxSIiILEDTeoTd5u\nDwwFIqIgMWCAgtmza5CSouBXv6rHkiU17d4HzikQEQWRvLx65OXVB+z1uadAREQahgIREWkYCkRE\npGEoEBGRhqFAREQahgIREWkYCkREpGEoEBGRhqFAREQahgIREWkYCkREpGEoEBGRhqFAREQahgIR\nEWkYCkREpGEoEBGRhqFAREQahgIREWkYCkREpGEoEBGRhqFAREQahgIREWlaHArFxcXo3r07Ro0a\n5XP/xx9/DFmWYTabYTabYTKZYDabsXbtWu0xixcvRlJSEqKjo9GvXz+UlpZqbbW1tcjLy0N8fDxi\nY2ORnZ2NioqKNrw1IiJqqRaFwvz58zFlyhTcdNNNDbb37NkTbrcbbrcb1dXVcLvdGDFiBACgsLAQ\ns2bNwqpVq/D9999jyJAhGDJkCKqrqwEAM2bMwP79+7F3714cPnwYqqoiNze3jW+PiIhaokWhYDKZ\n8Omnn+L6669v8QstXboUubm5sNvtCAsLw7Rp0yBJEgoLC6EoCgoKCjBz5kzExcUhOjoac+bMQVFR\nEcrKylr8WkRE1DotCoVJkyYhMjKy0XaXy4Xhw4fDZrMhPj4eCxcu1NocDgcyMjK025IkIS0tDSUl\nJThy5AgqKyuRnp6utScmJsJkMsHhcLSki0RE1AZ+m2iOiopCSkoKnnrqKZw8eRIFBQWYNWsWVq5c\nCQBwOp2wWq0+z4mJiUF5eTmcTickSbqk3Wq1ory83F9dJCKiyzD66xelp6dj69at2u27774beXl5\nWLFiBcaNGwcAEEI0+Tsu164XWZYgy1JAXhsADAbZ5yf5H2usP9ZYX+1VX7+FQkN69uypHX1ks9ng\ndDp92p1OJ5KTk2Gz2SCEgNPphNls1torKioQGxurZxcBADExEZCkwIXCeVFRpkB3ocNjjfXHGutL\n7/r6LRTeeecdlJeXIy8vT7vv0KFDSEhIAADY7XY4HA6MGTMGAKCqKkpLSzF+/HgkJCTAarXC4XAg\nPj4eAHDw4EHU1dXBbrf7q4uNqqioCvieQlSUCS5XNRRFDVg/OjLWWH+ssb4aq6/VGuHX1/FbKISG\nhuLpp5/GDTfcgP79+2Pbtm1YuXIl3njjDQDAhAkTkJOTg5ycHKSkpGD+/PkIDw/H4MGDIcsyHn30\nUcyZMwd2ux0mkwkzZszAiBEjYLPZ/NXFRqmqgKoGZujqvymKCo+Hf0x6Yo31xxrrS+/6tigUTCYT\nJElCfX09AGDdunWQJAlutxtZWVlYtGgRJk2ahOPHj6N79+5YvHgxhg0bBgAYNGgQ5s6di+zsbJw+\nfRp9+vTBhg0bEBYWBgCYPXs2zp07h9TUVCiKgqFDh+KVV17x89slIqKmSCJQs7t+5nK5YLFYUFlZ\niaioqBY99/Tpszr1qnmMRhlWawTOnKniFpZOWGP9scb6aqy+Nlvjpwm0Bg8TICIiDUOBiIg0DAUi\nItLoep4CEZE/lJbK2LPHgNRUFT//uRLo7nRoDAUiCmoffWRATo4JiiJBkgT+/OcaDB/uCXS3OiwO\nHxFRUHv3XSMUxXtyqRAS/vGPkAD3qGNjKBBRULvuOt+j5uPjebirnjh8RERB7fHH63DsmIQdO4xI\nS1MwY0ZtoLvUoTEUiCiohYUBixbVAmAYtAcOHxERkYahQJ3emjVGjB5twv/9Xxjc7kD3hiiwOHxE\nndrmzQZMmXJhffqzZ88PVRB1TtxToE7twAFDk7eJOhuGAnVqt92mQJIuHPJ4++08W5Y6Nw4fUad2\n++0K3nijGhs2GJGQIDBxYl2gu0QUUAwF6vQGDlQwcCD3EIgADh8RUTtTeUJyUGMoEFG7+OYbCbff\nbkZcXBeMGWNCTU2ge0QNYSgQUbv4/e/DcfiwAaoqobjYiIICLmwXjBgKRNQufvhB8rldWSk18kgK\nJIYCEbWLxx6rg8HgPfy3a1cVDzxQH+AeUUN49BFRAO3fL+PYMRm33abAZhOXf8IVbNgwD266yY0j\nR2T06aPgqqs69vu9UjEUiALktddC8L//GwYhJMTGqiguduOaazr2F+VPf6ripz/l4UfBjMNHRAGy\nbFkIhPCOq586JePdd7mNRoHHUCAKEKvVd6/giy8MqOMJ1RRgDAWiAFmwoBYJCQoAbzj8/e8hGDfO\n1PST6BI1NcC330qo57y1XzAUiAIkMVHFggW1AC4cmvnhh0a4XIHr05Xm3/+WkZkZgT59uuCOOyJw\n8iQPc20rhgJRAMXFqTAaLwwjde2qIiIigB26wrz4YihOnvR+jX3zjYz8/NAA9+jKx1AgCqCEBIH8\n/BrceKOC1FTviq0GXtKh2ZSL1jH0eALTj46EhztQp1VWJuG77yQkJano0iVw/Rg+3IPhw/lt1hq/\n/W0ddu0yoKJCRlycigkTOFPfVgwF6pS2bjUgN9eE6moJPXuqKCx0B/xkqgMHZDgcBvTurSAtjcfy\nN0dysoq9e6tw7JiMn/wksOHeUXD4iDql+fPDUF3tnZQ8elTGqlWBXZztww8N+NWvzPjd78Lxq1+Z\nsWULx5Cay2LxhgMDwT8YCtQufvjB+8V3+HBwfOTCwnz3CkIDPD/59tshUBRvSCmKhL//nSuIUmAE\nx18odWhlZRLuvDMCo0aZ8YtfmPHOO4EftXzuuVp06+YdoundW0FubmDHort39w2puLiOvdwFBa/A\n/3VSh/f22yE4ccK7/aEoEpYsCcX99wd2YjUtTcWBA1U4c0ZCbKyAFODD26dNq8WxYxI+/dSAzEwF\nU6fWBrZD1GkxFEh3kZGiydst9c9/yvjsMwPsdqVNi6uFhiLgk8vndekCrFjR+KXIvvlGQnm5hJQU\nFeHh7dgx6nQ4fES6GzWqHvfc412DIC5OxQsvtH4reNMmA+65x4ypU8Nx991m7N7d8SdkV60KwW23\nRWDIkAgMHWpGVVWge0QdGUOBWuXik4aaEhoKvP56DY4fP4sDB6rws5+1fut+zZoLE7J1dRLeeqvj\n7+y+8EIoVNX7nv/5TwOKijr+e6bAYShQi6gqMHlyOK65pgtSUiLgcDT/IxQW1vbXv3i4J1iGf85z\nOiW8+aYRH37ovz2Yi4eLTFwzj3TEUKAWKSoy4m9/C4GqSigrk/Hb37bvAPeMGbW4804PrFaBwYPr\n8eSTwXMGq9MpYeBAM5580oRRo8x47jk/pCCAF1+sQUSEN/wGD67Hvffy7GfSD/dDqUUCffF1qxV4\n663qdn3N5tq2zYDjxy9sZ73xRgiee67tRxHddZeCL788h3PnJHTtGlx7RtTxcE+BWmTIkHr07Hlh\nTuDxx4NnSz3QYmNFk7fbIiwMDARqF9xToBaJiQE++KAKu3cbcPXVAqmprZ803rHDgO++k9C/v3LJ\nyVtXon79vOcXrFgRgm7dBF5+ufFDTImCFUOBWsxiAe65pwWHHzUgPz8Es2d75yO6dVOxebMb1157\n5QfD9Ol1mD6de0905Wrx8FFxcTG6d++OUaNGXdK2detWZGZmwmKxIDk5GWvWrPFpX7x4MZKSkhAd\nHY1+/fqhtLRUa6utrUVeXh7i4+MRGxuL7OxsVFRUtOIt0ZXgtdcuLDZUXi7j/fe5fUIUDFoUCvPn\nz8eUKVNw0003XdJWVlaGYcOGYeLEiTh9+jQWLVqE8ePHa1/8hYWFmDVrFlatWoXvv/8eQ4YMwZAh\nQ1Bd7Z00nDFjBvbv34+9e/fi8OHDUFUVubm5fniLFIz0HH8notZrUSiYTCZ8+umnuP766y9pW716\nNRITEzF27FiEhoZiwIAByMrKwvLlywEAS5cuRW5uLux2O8LCwjBt2jRIkoTCwkIoioKCggLMnDkT\ncXFxiI6Oxpw5c1BUVISysjL/vFMKKosW1eBnP1NgsQg8/HAd7ruPh1kSBYMWhcKkSZMQGRnZYJvD\n4UBGRobPfRkZGSgpKWmwXZIkpKWloaSkBEeOHEFlZSXS09O19sTERJhMJjgcjpZ0ka4QN96oYutW\nN/7973OYO7c24AvSEZGX3w5JdTqdsFqtPvfFxMSgvLz8su1OpxOSJF3SbrVatecTEZH+/Dq7J0TT\n48JtbdeLLEuQ5cBtqhoMss9P8j/WWH+ssb7aq75+CwWbzQan0+lzn9PpRGxsbJPtycnJsNlsEELA\n6XTCbDZr7RUVFdrz9RQTEwEpCMYvoqK4qI3eWGP9scb60ru+fgsFu92OlStX+txXUlKCzMxMrd3h\ncGDMmDEAAFVVUVpaivHjxyMhIQFWqxUOhwPx8fEAgIMHD6Kurg52u91fXWxURUVVwPcUoqJMcLmq\noSi8YHtjvvtOwrhxYfjiCxkDBihYurS22dcWYI31xxrrq7H6Wq0Rfn0dv4XC6NGj8dxzz6GgoACj\nR4/Gli1bsHHjRuzduxcAMGHCBOTk5CAnJwcpKSmYP38+wsPDMXjwYMiyjEcffRRz5syB3W6HyWTC\njBkzMGLECNhsNn91sVGqKqCqgT8kUlFUeDwd44/p1CkJFRUSbrhBhdFPn7Jnnw1Haal39dGiIiNe\neUVp8YJ4HanGwYo11pfe9W3Rn6vJZIIkSaiv914wZd26dZAkCW63GzabDUVFRXjiiSfw+OOPo2fP\nnli9ejV69eoFABg0aBDmzp2L7OxsnD59Gn369MGGDRsQ9uN6yrNnz8a5c+eQmpoKRVEwdOhQvPLK\nK35+u9QeCguNmDAhHHV1Evr0UfDOO26/LPfsdEpN3iaitpNEoGZ3/czlcsFisaCyshJRUVEteu7p\n02d16lXzGI0yrNYInDlT1SG2sHr3jvBZLXTx4mo8+GDbz0NYv96Ixx4Lh6pKiIwUKCx04+abm1ev\njlbjYMQa66ux+tpsDZ8m0OrX8etvIwIuOefAX3P4w4Z5kJDgxldfybjlFgXXXdchtmeIggqPHSO/\ne/75GphM3i/sO+7w+PVs5eRkFfff72EgEOmEewrkd/fco+Czz87hzBkJ110nIP+46XH4sIxZs8Lg\ndgNPPlmHO+9s20qrROR/DAXShcUCWCwXtuaFAB580ITvvvMmhMNhwK5dVW1aLtvjARYsCMX+/Qbc\nequCyZPrtAAiotZhKFC7cLmgBQIA1NRIOHpUxrXXtn5v4aWXQvGnP3mPXtu2zQizWeCxx+rb3Fei\nzozbVdQuLBagT58LAXDVVSp69Wrb8NE//+n78f3sM0Obfl+g1Lb9Ms5EfsNQoHbz5ptuTJ9eiyee\nqMX777tx0fqHLdavn2+o3HHHlbX89mefyUhNjUB8fCR+/etw1PGCbRQEOHxE7SYyEpg61X/ffI88\nUg+zWcDhMKBvXwX3339lhcL06eE4edK7XbZpUwjeekvBmDEc/qLAYijQFW3UKA9GjbqywuC8c+d8\nb58N7DmURAA4fETUoOPHJeTnh+DNN41QdTo594kn6iDL3qOv4uNVjBhxZYYbdSzcUyC6SFmZhEGD\nzCgv924z7d1bh4UL/T8bnJ3tQWqqG999J6F3bwXR0X5/CaIW454C0UW2bzdogQAA69aF6PZaiYkq\nBgxgIFDw4J5CEHnvPQP+9S8j+vdXkJkZ/Gf7qio65MliFy+h0aMHF3ejzqMD/klfmV56CRg3Lhx/\n/GMY7rvPhB07gveY+7IyCb/8pRlXX90Fw4aZUFkZ6B4BdXXA738fhj59IjBwoBkff9z6+mVmKpg7\ntwY33qigb18Pli+v8WNPiYIbQyFI/OMfF/6tKBKKi4N3J27OnDB89pkBQkjYvduIxYtDtTZFATZv\nNqC42ABPO86bLlwYiqVLQ/HttzIOHDAgO9uEnTtbHwwPP1yPnTvdWL++GjfeyD0F6jwYCkHippt8\nbyckBOcX0ZkzwIkTvmth//CD97YQQG5uOB56yIwxY8z49a9Nuh25c7Evv/T9KAshYevW4N3bIgpW\nDIUg8eyzQO/eCq69VsXjj9di7NjgO4npk08MSE/vgh07jJAk77h7ly5C6+vRoxI2bbowKfvhh0Yc\nOdI+H7FBgy7dLUlKCs5gJQpmwTtG0YmcPQsMGgQcPuzdsj15Ug7KCdznnw+D231+r0BCdnYdnnmm\nTlvpNDISMBoFPB7vYwwGgaio9rnuwYMPemAwVCM/PxT19cADD3gwciSP+ydqqSD86ul89u414PDh\nC7fXrTPqukhaVRUwfnw4UlMj8JvfhF9yZm1jLr6C2s03qz5LX3frJrBwYQ0iIwW6dBH44x9rcNVV\n7XcxnJEjPfj4Yzd27XJj8mQuJETUGgyFIBAXp/p84XbtKhAWpt/rLVwYivXrQ3DypIyiohDMn9+8\nF5s5sxZduni/5FNTG16n54EHPDhy5By++eZcq5afWLfOiOHDTcjLC8f33/vpOp5E1GwcPgoCN98s\n8PLLwPPPq4iM9G5h6+m/r2sAeCeOhbj8tZR//nMFBw6cQ3m594pqRj9/evbvlzFhQjhU1duR//xH\nwnvvVfv3RYioSdxTCBITJgAHD1Zj+3Y3brlF3wnSkSPrYTB4t/hlWeD0aQnXXNMFt9wSgUOHmv5I\nREUBCQktDwSPB/jrX0Pwwguh+Pe/G36NL76QtUAAgEOHePQQUXvjnoKOKiq8X6L+3qJuqwEDFGzY\n4Ma+fQY4ncCf/hQOwHv00NNPh6OoyI2jRyVERwvExPjnNSdPDsfbb3uPTFq+PBTbtlUhPt53vuHW\nWxWYzUKbzL7rLk4UE7W3IPu66hiqqoCcHBP27DEiNlbF3/5WjeTk4Do8Mj1dRXq6ij//2XddH6cT\nGD7chF27jAgLE3jllRoMHdr2L+eNGy981FwuCZs2GbFnjwFffy0jOVnBk0/W46abVLz3nhtvvx2C\n2FiBRx/lZDFRe5OEEO13eIiOXC4XLBYLKisrERUV1aLnnj7t34Xs8/NDMHt2uHb75z/3YN26xsfG\njUYZVmsEzpypgsfTvuFxfkVQ78VeBAwG7xnV53XvruKzz6ra/DqDBpmxf793OEiSBHr3VrBv34Wg\nMBoF1qypRv/++qz5FMgadxassb4aq6/NFunX1+Gcgg5qanxnbKurg/comu7dBbZudaOgoBpRUb6B\nAHjPUvaHZcuqMXCgB+npChYurMGpU74fPY9Hwuuv67caKRE1D0NBB2PG1OO667xJHhYmMHVqcF+Z\nvWtXgSFDPJesVWQwCMye7Z++X3edwKpV1SgudmPUKE+DZyBbrR1ip5XoisZQ0MFVVwls21aF9evd\n2L27CnffHTzLYJ896z30s6Li0rZnn63Vlq9ITVWwb18V/ud/9Jns/cMfavH88zWIj1dhNApkZnrw\nu99xDoEo0DinAP/PKfw3jweYOzcUDocBt96qYNq0OhguOtKyvcZijx6VMGyYd/7AYhH4+9/dSE/3\nfb1vvpFw5oyElBQVIR1oNIfj3fpjjfXVXnMKPPpIZwsXhmLJEu8Zw7t2GdGli8CkSU0vdvfGGyHY\ntcuA3r0VPPxw/WVPKmuupUtDf5xQBiorJSxcGIrXX/c9US4hQQDoENsJRNQKDAWdff65fNFtA4DG\nQ2HVqhBMneo9cmnt2hDU1QETJ/pnxdSLz5foSHsCROQfnFPQ2V13+c4n3Hln02P0u3f7ji3t2eO/\ns3onTapDYqK3P3FxKp55hmP4ROSLewo6+/Wv6xERIbQ5hayspkOhd29FO/PXe9t/Y7OxsQLbtrlx\n6pSEbt0EQkMv/xwi6lwYCu1gxAgPRoxo3lE8ubn1qK/37jFkZKiYNMm/W/NGIxAXxzkDImoYQyHI\nSBLw2GP1eOyx4LvyGhF1fJxTICIiDUOBiIg0DAUiItIwFIiISMNQICIiDUOBiIg0DAUiItIwFIiI\nSMNQICIijV9DQZZlmEwmmM1m7efkyZMBAFu3bkVmZiYsFguSk5OxZs0an+cuXrwYSUlJiI6ORr9+\n/VBaWurPrhERUTP4dZkLSZJw+PBhxMfH+9xfVlaGYcOGIT8/Hzk5Odi+fTuysrKQlJSEjIwMFBYW\nYtasWSguLkZycjJeeuklDBkyBEeOHIHJZPJnF4mIqAl+3VMQQqChC7mtXr0aiYmJGDt2LEJDQzFg\nwABkZWVh+fLlAIClS5ciNzcXdrsdYWFhmDZtGiRJQmFhoT+7R0REl+H3OYXp06ejR48esFqtyMvL\nQ1VVFRwOBzIyMnwel5GRgZKSEgC4pF2SJKSlpWntdMH27QakpESgZ88u+OMfufY1EfmXX0Ohb9++\nGDhwIL7++mvs2bMHe/bswcSJE+F0OmG1Wn0eGxMTg/LycgC4bDtd8Oij4Sgrk+F2S5g3LwwOB48V\nICL/8eucws6dO7V/JyYm4oUXXsDQoUPRr1+/BoeV/tvl2vUkyxJk2U8XQm4Fg0H2+dkYRQF++MG3\nn5WVhksus0mXam6NqfVYY321V311/Trp2bMnFEWBLMtwOp0+bU6nE7GxsQAAm83WYHtycrKe3dPE\nxERAkgIXCudFRV06qf7mm8AXXwCDBwOZmcDEicCSJd62lBRgyJBwRES0c0evYA3VmPyLNdaX3vX1\nWygcOHAAq1atwoIFC7T7Dh06hPDwcAwePBgrV670eXxJSQkyMzMBAHa7HQ6HA2PGjAEAqKqK0tJS\nPPLII/7qXpMqKqoCvqcQFWWCy1UNRblw+c2FC0Pwhz945w3mzhUoLKzBrFkq+veXUVkp4Ze/VFBX\nB9TxUsuX1ViNyX9YY301Vl+r1b9bhX4LhdjYWCxduhSxsbGYMmUKjh49ipkzZ+Kxxx7DQw89hFmz\nZqGgoACjR4/Gli1bsHHjRuzduxcAMGHCBOTk5CAnJwcpKSmYP38+wsPDce+99/qre01SVQFVDfwl\nKhVFhcdz4X92YaFB+3d9vYRNm2RkZHhwxx0XHuNp3lU+6UcX15j8jzXWl9719dvgVFxcHDZs2ID1\n69ejW7duuP322zF48GDMmzcPNpsNRUVFWLJkCaKjozF16lSsXr0avXr1AgAMGjQIc+fORXZ2Nrp2\n7YotW7Zgw4YNCAsL81f3rkjXX+/7Pz4hgX9oRKQvSQRyhtePXC4XLBYLKisrERUV1aLnnj59Vqde\nNY/RKMMzT1NeAAAIg0lEQVRqjcCZM1U+WwA//ABMnx6Ow4dl3HOPB9Onc5yotRqrMfkPa6yvxupr\ns0X693X8+tvIr6Kjgb/8pSbQ3SCiToTHjhERkYahQEREGoYCERFpGApERKThRHOQO3cOKC42wmwG\n7rnHgyA48ZqIOjCGQhCrrgaGDTPjX//ynsQ2cmQ9Xn6ZRyMRkX44fBTEHA6DFggA8PbbIaiqCmCH\niKjD455CMx09KuEvfwmF0QhMmlSHq67S/5y/bt0EJElACO+YUWSkQHi4t01RgNdeC8GJExKysjxI\nTeXJQkTUdgyFZjh71juMc/Kkd8dq2zYDPvrIrfuS1UlJKp5/vhYLFoQhIkLgT3+qgeHHHYff/S4M\nK1d6F8tbtiwUmze7kZTEYCCituHwUTN8/bWsBQIAHD5swPfft8+M7/jx9fjqq3MoLa1C//6Kdv/m\nzRcSqaZGwiefGBp6OhFRizAUmqFHDxUWy4Xhorg4FTZbYJeMSkxUm7xNRNQaHD5qhpgY4O233Vi4\nMBQhIcD06XUIDfDlkZcsqcHvfx+G776TMXJkPX7xC+XyTyIiugyGQjOlpal47bXgORw0NlZg6dLg\n6Q8RdQwcPiIiIg1DgYiINAwFIiLSMBSIiEjDieYGvP56CA4ckNG3r4KRIz2B7g4RUbthKFzk1VdD\nMHOmdy2JVasAIaqRnc1gIKLOgcNHF9mxw9jkbSKijoyhcJHUVKXJ20REHRk3gy/y29/WQVWBAwcM\n6NtXwW9+Ux/oLhERtRuGwkWMRu8yFkREnZEkhAjsym5+4nK5YLFYUFlZiaioqEB3h4joitRhQkEI\ngbNnzyIyMhISL2RMRNQqHSYUiIio7Xj0ERERaRgKRESkYSgQEZGGoUBERBqGAhERaRgKRESkYSgQ\nEZGGoUBERBqGAhERaRgKRESk6fSrpJ5fM4mI6ErlzzXfOn0onD17FhaLJdDdICJqNX+uDt3pF8Tj\nngIRXen8uafQ6UOBiIgu4EQzERFpGApERKRhKBARkYahQEREGoYCERFpGApERKRhKBARkYahQERE\nGoaCzmRZhslkgtls1n5OnjwZALB161ZkZmbCYrEgOTkZa9as8Xnu4sWLkZSUhOjoaPTr1w+lpaWB\neAtBqbi4GN27d8eoUaMuaWtLXWtra5GXl4f4+HjExsYiOzsbFRUVur+fYNNYfT/++GPIsgyz2ezz\nmV67dq32GNb38o4dO4bhw4ejW7duuPrqq5GbmwuXywUgCD6/gnQly7I4duzYJfefPHlSdOnSRaxc\nuVLU1taKDz/8UJjNZuFwOIQQQrz33nsiJiZGlJSUiJqaGjFv3jxx9dVXC7fb3d5vIei8+OKLIikp\nSdxxxx0iJyfHp62tdX3qqafELbfcIk6cOCHOnDkjRowYIbKystr9PQZSU/X96KOPxE9+8pNGn8v6\nNk9KSop4+OGHhdvtFidOnBB9+vQR48ePD4rPL0NBZ5IkiW+//faS+xcsWCB69+7tc9+DDz4oJkyY\nIIQQYsiQIWLq1Klam6qqIi4uTrz11lv6dvgKsGTJEuFyucS4ceMu+dJqS109Ho+Ijo4WRUVFWvuX\nX34pZFkWJ0+e1PEdBZem6nu5UGB9L++HH34QDz/8sDh16pR2X35+vkhMTAyKzy+Hj9rB9OnT0aNH\nD1itVuTl5aGqqgoOhwMZGRk+j8vIyEBJSQkAXNIuSRLS0tK09s5s0qRJiIyMbLCtLXU9cuQIKisr\nkZ6errUnJibCZDLB4XDo8E6CU1P1BQCXy4Xhw4fDZrMhPj4eCxcu1NpY38uzWCxYvnw5bDabdt/x\n48dxzTXXBMXnl6Ggs759+2LgwIH4+uuvsWfPHuzZswcTJ06E0+mE1Wr1eWxMTAzKy8sB4LLt1LC2\n1NXpdEKSpEvarVYr6/6jqKgopKSk4KmnnsLJkydRUFCAWbNmYeXKlQBY39bYt28f8vPz8eyzzwbF\n55ehoLOdO3ciNzcXISEhSExMxAsvvIA1a9bA4/FAXGaB2su1U8PaWlfWvXHp6enYunUrbr/9dhiN\nRtx9993Iy8vDihUrtMewvs23c+dODBo0CPPmzcNdd90FIPCfX4ZCO+vZsycURYEsy3A6nT5tTqcT\nsbGxAACbzdZkOzXscnVrqt1ms0EIcUl7RUUF696Enj174j//+Q8A1rclCgsLce+992Lx4sV4/PHH\nAQTH55ehoKMDBw7g6aef9rnv0KFDCA8Px+DBg7Fv3z6ftpKSEmRmZgIA7Ha7zzigqqooLS3V2qlh\nF9cNaF5db731ViQkJMBqtfq0Hzx4EHV1dbDb7e3zBoLcO++8g1dffdXnvkOHDiEhIQEA69tcu3bt\nwrhx47B27VqMHj1auz8oPr8tmDSnFjpx4oSIjIwU8+bNE7W1teKrr74SvXr1ElOmTBGnTp0SFotF\n/PWvfxU1NTXi/fffFxEREeLgwYNCCCE2bdokrFar2LNnj3C73WLWrFmiR48eoqamJsDvKng0dHRM\nW+v6zDPPCLvdLo4fPy7Ky8vF0KFDxQMPPNDu7y0YNFTf9evXi4iICPHBBx+I+vp6sXnzZhEZGSne\nffddIQTr2xwej0fcfPPNYtmyZZe0BcPnl6Ggs+3bt4vbbrtNREZGCpvNJqZNmyZqa2u1trS0NBEe\nHi6SkpK0P6zzXn31VXHdddcJk8kk+vXrJz7//PNAvIWgEx4eLkwmkzAajcJoNGq3z2tLXevq6sSk\nSZNETEyMsFgs4qGHHhIul6vd3lswuFx9ly1bJhITE4XZbBYJCQlixYoVPs9nfZu2fft2IcuyMJlM\nWm3P/zx27FjAP7+8HCcREWk4p0BERBqGAhERaRgKRESkYSgQEZGGoUBERBqGAhERaRgKRESkYSgQ\nEZGGoUBERBqGAhERaRgKRESkYSgQEZHm/wHf4Aaq0Br4EgAAAABJRU5ErkJggg==\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# sklearnのRidge回帰を使って最適解を見つける\n", "from sklearn.linear_model import Ridge\n", "clf = Ridge(alpha=1.0)\n", "clf.fit(X, y) \n", "pred = clf.predict(X)\n", "list_plot(zip(y, pred), figsize=4, zorder=2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\t\n", "\t\t井出本では、リッジ回帰のλについても最適な値を求めていますが、\n", "\t\tここではλ=1.0で計算した結果を使って異常度を計算します。\n", "\t
\n", "\t\n", "\t\t観測値と回帰結果との残差の自乗は、score関数で求まります。\n", "\t\tあとは、係数α（coefs）を使って$\\hat{\\sigma}^2_{ridge}$を求め、\n", "\t\t観測値yと変数xの中心値を使って異常度を計算します。\n", "\t
\n", "" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEBCAYAAACXArmGAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3X18k/W9//HXlaQ3SaFNCilQrKcHpu105aYWOzmT40Dt\nhAIKxyoio3UrtIjOKehOd8aAURm6M7GKc4jIQWAbyrxpVw4ioDjw1JqiszrgNxwCYkubYgtNSGly\n/f6IvUpshIa2pE0+z8fDB+b65kq/1yc373y/15XrUlRVVRFCCCEAXbA7IIQQoveQUBBCCKGRUBBC\nCKGRUBBCCKGRUBBCCKGRUBBCCKGRUBBCCKEJmVBQVZWmpibkZxdCCHHxDMHuQKDq6k75Xd7cfIp/\n/deh/POfnxMT0/8S96r30ukU4uNjaGhoxuORwGwjdfFP6uJfKNTFau3c52LIjBQURfH5V3jpdAqK\noqDTSV3OJXXxT+riXzjVJWRCQXRUXBxJaqqRjAw4eDD0X8xCiK6TUAhR27bpefLJKE6c0GGzwbx5\nUcHukhCiD5BQCFFffKH72m0ZKQghLkxCIUTdfHMrAwd6tNszZ7YGsTdCiL6izx19JDonMVHlzTcd\n7NwZQUpKFNddd5ZWyQUhxAVIKISwxESV3NxWLJYoTp4Mdm+EEH2BTB8JIYTQSCgIIYTQSCgIIYTQ\nSCgIIYTQSCgIIYTQSCgIIYTQSCgIIYTQSCgIIYTQSCgIIYTQSCgIIYTQBBwKH374ITfeeCNms5kh\nQ4Zw5513Ultby9tvv41Op8NkMmEymTAajZhMJrZs2aKtW1JSQmpqKmazmXHjxlFVVaW1uVwuCgoK\nSEpKIiEhgZycHBoaGrpnK4UQQnRKQKHQ0tJCVlYW48ePp66ujurqampra5k3bx4AycnJOBwOHA4H\nTqcTh8PB9OnTASgtLWXJkiVs2LCB2tpasrOzyc7Oxul0AlBUVMS+ffuoqKjg4MGDeDwe8vLyunlz\nhRBCnE9AoeBwOHj00Uf52c9+RkREBAMGDGDatGlUV1dfcN3Vq1eTl5dHRkYGUVFRLFy4EEVRKC0t\nxe12s3btWhYtWkRiYiJms5ni4mLKysqoqam56I0TQggRmIBCwWw2c88996DTeVc7cOAA69at4847\n7wSgqamJadOmYbVaSUpK4oknntDWtdlspKena7cVRWHUqFFUVlZy6NAhGhsbGT16tNaekpKC0WjE\nZrN1aQOFEEJ03kWdOvvIkSNcccUVuN1u5syZw+LFi/nggw8YMWIEDz74IJs3b2bXrl3cfvvtWCwW\ncnNzsdvtWCwWn8eJj4+nvr4eu92Ooigd2i0WC/X19T7LdDr/F8/W63XavwaD7D9vc25dRDupi39S\nF//CqS4XFQqXX345LpeLQ4cOMWfOHO6++242btzIzp07tfvcdNNNFBQU8MILL5CbmwuAqqrnfdwL\ntQPEx8egKP5CwQ1AbKyR2NiYALYmPMTGGoPdhV5J6uKf1MW/cKhLly6yM3z4cIqLixk7diwlJSUM\nGDDApz05OVk7+shqtWK3233a7XY7aWlpWK1WVFXFbrdjMpm09oaGBhISEnzWaWho9jtSaG727rBu\nanLiduu7slkhRa/XERtr/KounguvECakLv5JXfwLhbpYLJ37shxQKOzatYvCwkL279+vLVMUBUVR\nKC8vp7m5mYKCAq3tk08+YdiwYQBkZGRgs9mYNWsWAB6Ph6qqKvLz8xk2bBgWiwWbzUZSUhIA1dXV\ntLS0kJGR4dMHj0fF4+k4omh7otxuD62tffNJ60lSF/+kLv5JXfwLh7oENEF2zTXX0NTUxCOPPILT\n6aSuro4lS5Ywbtw4zGYzDz30EG+++Satra1s376ddevWaYerFhYWsn79eioqKnA6nSxbtozo6Ggm\nTpyITqdjzpw5FBcXc+zYMex2O0VFRUyfPh2r1dojGy6EEKKjgEYKsbGxbN++nfnz52O1WunXrx/j\nx4/n+eefZ8iQITz55JPMnz+fo0ePMnjwYEpKSpg6dSoAWVlZLF++nJycHOrq6hgzZgzl5eVERUUB\nsHTpUk6fPs3IkSNxu91MnjyZZ555pvu3WAghxDdS1M7s3e1F6upO+V3ucJwmOTmRw4ePYzL1u8S9\n6r0MBh0WSwwnTzaH/LA3EFIX/6Qu/oVCXazW/p26X+gfXyWEEKLTJBSEEEJoJBSEEEJoJBSEEEJo\nJBSEEEJoJBSEEEJoJBSEEEJoJBSEEEJoJBSEEEJoJBSEEEJoJBSEEEJoJBSEEEJoJBSEEEJoJBSE\nEEJoJBSEEEJoJBSEEEJoAg6FDz/8kBtvvBGz2cyQIUO48847OXHiBAA7d+4kMzOTuLg40tLS2LRp\nk8+6JSUlpKamYjabGTduHFVVVVqby+WioKCApKQkEhISyMnJoaGhoYubJ4QQIhABhUJLSwtZWVmM\nHz+euro6qqurqa2tpbCwkJqaGqZOncq8efOoq6tj5cqV5Ofnax/8paWlLFmyhA0bNlBbW0t2djbZ\n2dk4nU4AioqK2LdvHxUVFRw8eBCPx0NeXl73b7EQQohvFFAoOBwOHn30UX72s58RERHBgAEDmDZt\nGtXV1WzcuJGUlBRmz55NZGQkEyZMYMqUKaxZswaA1atXk5eXR0ZGBlFRUSxcuBBFUSgtLcXtdrN2\n7VoWLVpEYmIiZrOZ4uJiysrKqKmp6ZENF0II0VFAoWA2m7nnnnvQ6byrHThwgHXr1nHHHXdgs9lI\nT0/3uX96ejqVlZUAHdoVRWHUqFFUVlZy6NAhGhsbGT16tNaekpKC0WjEZrNd9MYJIYQIjOFiVjpy\n5AhXXHEFbrebOXPmsHjxYm655RaSkpJ87hcfH099fT0Adrsdi8Xit91ut6MoSod2i8Wird9Gp1PQ\n6ZQOfdLrddq/BoPsP29zbl1EO6mLf1IX/8KpLhcVCpdffjkul4tDhw4xZ84cZs2aBYCqquddr6vt\nAPHxMSiKv1BwAxAbayQ2NuaCjxNuYmONwe5CryR18U/q4l841OWiQqHN8OHDKS4uZuzYsUyaNAm7\n3e7TbrfbSUhIAMBqtfptT0tLw2q1oqoqdrsdk8mktTc0NGjrty9r9jtSaG727rBuanLiduu7slkh\nRa/XERtr/KounmB3p9eQuvgndfEvFOpisXTuy3JAobBr1y4KCwvZv3+/tkxRFBRF4dprr+Xll1/2\nuX9lZSWZmZkAZGRkYLPZtFGFx+OhqqqK/Px8hg0bhsViwWazaVNQ1dXVtLS0kJGR4fOYHo+Kx9Nx\nRNH2RLndHlpb++aT1pOkLv5JXfyTuvgXDnUJaILsmmuuoampiUceeQSn00ldXR1Llixh3LhxFBYW\n8tlnn7F27VpcLhfl5eVs3bqVuXPnAlBYWMj69eupqKjA6XSybNkyoqOjmThxIjqdjjlz5lBcXMyx\nY8ew2+0UFRUxffp0rFZrj2y4EEKIjgIKhdjYWLZv3857772H1WolLS0Ns9nMpk2bGDhwIGVlZTz1\n1FOYzWYeeughNm7cyNVXXw1AVlYWy5cvJycnhwEDBrBjxw7Ky8uJiooCYOnSpXz3u99l5MiRDB8+\nnLi4OJ577rnu32IhhBDfSFE7s3e3F6mrO+V3ucNxmuTkRA4fPo7J1O8S96r3Mhh0WCwxnDzZHPLD\n3kBIXfyTuvgXCnWxWvt36n6hf3yVEEKITpNQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQ\nEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEII\noQk4FI4cOcK0adMYOHAgQ4YMIS8vj8bGRt5++210Oh0mkwmTyYTRaMRkMrFlyxZt3ZKSElJTUzGb\nzYwbN46qqiqtzeVyUVBQQFJSEgkJCeTk5NDQ0NA9WymEEKJTAg6FyZMnEx8fz9GjR7HZbHz88ccs\nXLgQgOTkZBwOBw6HA6fTicPhYPr06QCUlpayZMkSNmzYQG1tLdnZ2WRnZ+N0OgEoKipi3759VFRU\ncPDgQTweD3l5ed24qUIIIS4koFBobGxkzJgxLF++HKPRSGJiIrNnz2b37t0XXHf16tXk5eWRkZFB\nVFQUCxcuRFEUSktLcbvdrF27lkWLFpGYmIjZbKa4uJiysjJqamoueuOEEEIEJqBQiIuLY82aNVit\nVm3ZkSNHGDp0KABNTU1MmzYNq9VKUlISTzzxhHY/m81Genq6dltRFEaNGkVlZSWHDh2isbGR0aNH\na+0pKSkYjUZsNttFb5wQQojAGLqy8vvvv8+qVasoLS0lNjaWESNG8OCDD7J582Z27drF7bffjsVi\nITc3F7vdjsVi8Vk/Pj6e+vp67HY7iqJ0aLdYLNTX13eli0IIIQJw0aGwZ88epkyZwooVK/j+978P\nwM6dO7X2m266iYKCAl544QVyc3MBUFX1vI95oXYAnU5Bp1M6LNfrddq/BoMcVNXm3LqIdlIX/6Qu\n/oVTXS4qFEpLS5k1axarVq1i5syZ33i/5ORk7egjq9WK3W73abfb7aSlpWG1WlFVFbvdjslk0tob\nGhpISEjwWSc+PgZF8RcKbgBiY43ExsZczGaFtNhYY7C70CtJXfyTuvgXDnUJOBT27t1Lbm4uW7Zs\nYcKECdryl19+mfr6egoKCrRln3zyCcOGDQMgIyMDm83GrFmzAPB4PFRVVZGfn8+wYcOwWCzYbDaS\nkpIAqK6upqWlhYyMDJ+/39DQ7Hek0NzsPYqpqcmJ260PdLNCll6vIzbW+FVdPMHuTq8hdfFP6uJf\nKNTFYuncl+WAQsHtdpOfn8+KFSt8AgEgMjKSBQsW8K1vfYsbbriBXbt2sW7dOl588UUACgsLmTFj\nBjNmzGDEiBE8/vjjREdHM3HiRHQ6HXPmzKG4uJiMjAyMRiNFRUVMnz7dZ6c2gMej4vF0nGZqe6Lc\nbg+trX3zSetJUhf/pC7+SV38C4e6BBQK7777Lvv37+f+++/nvvvuQ1EUVFVFURQOHDjAypUrmT9/\nPkePHmXw4MGUlJQwdepUALKysli+fDk5OTnU1dUxZswYysvLiYqKAmDp0qWcPn2akSNH4na7mTx5\nMs8880z3b7EQQohvpKid2bvbi9TVnfK73OE4TXJyIocPH8dk6neJe9V7GQw6LJYYTp5sDvlvOIGQ\nuvgndfEvFOpitfbv1P1Cf1e6EEKITpNQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEII\noZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQ\nEEIIoQk4FI4cOcK0adMYOHAgQ4YMIS8vj6amJgB27txJZmYmcXFxpKWlsWnTJp91S0pKSE1NxWw2\nM27cOKqqqrQ2l8tFQUEBSUlJJCQkkJOTQ0NDQxc3TwghRCACDoXJkycTHx/P0aNHsdlsfPzxxyxY\nsICamhqmTp3KvHnzqKurY+XKleTn52sf/KWlpSxZsoQNGzZQW1tLdnY22dnZOJ1OAIqKiti3bx8V\nFRUcPHgQj8dDXl5e926tEEKI8wooFBobGxkzZgzLly/HaDSSmJjI7Nmz2b17Nxs3biQlJYXZs2cT\nGRnJhAkTmDJlCmvWrAFg9erV5OXlkZGRQVRUFAsXLkRRFEpLS3G73axdu5ZFixaRmJiI2WymuLiY\nsrIyampqemTDhRBCdBRQKMTFxbFmzRqsVqu27OjRowwdOhSbzUZ6errP/dPT06msrATo0K4oCqNG\njaKyspJDhw7R2NjI6NGjtfaUlBSMRiM2m+2iNkyEPlWFzz9XOH062D0RInQYurLy+++/z9NPP83r\nr7/OihUrSEpK8mmPj4+nvr4eALvdjsVi8dtut9tRFKVDu8Vi0dZvo9Mp6HRKh77o9TrtX4NB9p+3\nObcuoeTsWbjrrih27DBgNKo8/7yLH/zA3en1Q7UuXSV18S+c6nLRobBnzx6mTJnCihUrGD9+PCtW\nrEBV1fOu09V2gPj4GBTFXyh4PxBiY43ExsZc8HHCTWysMdhd6FZ/+APs2OH9f6dT4T//M5oZMwJ/\nnFCrS3eRuvgXDnW5qFAoLS1l1qxZrFq1ipkzZwJgtVqx2+0+97Pb7SQkJJy3PS0tDavViqqq2O12\nTCaT1t7Q0KCt376s2e9IobnZu8O6qcmJ262/mM0KSXq9jthY41d18QS7O93myy8NQJR2u6XFw8mT\nzk6vH6p16Sqpi3+hUBeLpXNflgMOhb1795Kbm8uWLVuYMGGCtjwjI4N169b53LeyspLMzEyt3Waz\nMWvWLAA8Hg9VVVXk5+czbNgwLBYLNptNm4Kqrq6mpaWFjIwMn8f0eFQ8no4jirYnyu320NraN5+0\nnhRqdZk0qYWMDAPvv6/HYFD5xS9cF7V9gdTl8GGF11+PYNAgD7ff3oouhGcSQu310l3CoS4BhYLb\n7SY/P58VK1b4BALAzJkzWbx4MWvXrmXmzJns2LGDrVu3UlFRAUBhYSEzZsxgxowZjBgxgscff5zo\n6GgmTpyITqdjzpw5FBcXk5GRgdFopKioiOnTp/vs1BaijdEIr73mYP9+HQMHqgwZcuGpx644dkwh\nKyuGkye9o1SbrYXHHnP16N8UIhgC+q7z7rvvsn//fu6//36MRiMmk0n798yZM5SVlfHUU09hNpt5\n6KGH2LhxI1dffTUAWVlZLF++nJycHAYMGMCOHTsoLy8nKso7BbB06VK++93vMnLkSIYPH05cXBzP\nPfdc92+xCBkREZCW5unxQADYvVuvBQLAa69F9PjfFCIYFLUze3d7kbq6U36XOxynSU5O5PDh45hM\n/S5xr3ovg0GHxRLDyZPNIT/sDUSgddmzR89tt7Xv7xo1ys0bbzh6sotBIa8X/0KhLlZr/07dL4Rn\nRYXoPv/2b26WLj3DFVe4+d73Wnn22c7v1BaiL+nS7xSECCcFBWcpKDgb7G4I0aNkpCCEEEIjoSCE\nEEIjoSCEEEIjoSCEEEIjoSCEEEIjoSCEEEIjoSCEED3g9dcNTJliJDc3mqNHO57Es7eS3ykIIUQ3\nq67WMXduNG63NwyOHtWxY0ff+AW8hILo4P33dWzZEsGQISoFBS1ERga7R0L0Lf/v/+m0QADYv7/v\nTMpIKAgf+/fruO02Ey6X9wV94ICOVavOBLlXQvQt117rJi5OpbHR+z6aMKE1yD3qPAkF4ePdd/Va\nIAC89ZZcsEiIQA0dqlJW5uBPfzJgsUB+fkuwu9RpEgrCx1VXeVAUFVX1BsN3vtM3zwgpRLClpHhY\ntKjvhEEbCQXhIzPTzapVZ/jjH737FBYvlgvJCBFOJBREB//xH638x3/0nTlQIUT3CXiX+LZt2xg8\neDB33XWXz/K3334bnU6HyWTyuSLbli1btPuUlJSQmpqK2Wxm3LhxVFVVaW0ul4uCggKSkpJISEgg\nJyeHhoaGLmyaEEKIQAUUCo8//jgPPPAAV155pd/25ORkHA4HDocDp9OJw+Fg+vTpAJSWlrJkyRI2\nbNhAbW0t2dnZZGdn43R6L1ZSVFTEvn37qKio4ODBg3g8HvLy8rq4eUKEnzNnYO9ePQcP9p3DIEXv\nEdCrxmg08t577zF8+PCA/9Dq1avJy8sjIyODqKgoFi5ciKIolJaW4na7Wbt2LYsWLSIxMRGz2Uxx\ncTFlZWXU1NQE/LeECFfNzZCdbeLWW01cf72J55+Xa0mLwAQUCvPnz6d//2++zmdTUxPTpk3DarWS\nlJTEE088obXZbDbS09O124qiMGrUKCorKzl06BCNjY2MHj1aa09JScFoNGKz2QLpohBh7X//18Df\n/uY9jFhVFR5/XH55KALTbTuaY2NjGTFiBA8++CCbN29m165d3H777VgsFnJzc7Hb7VgsFp914uPj\nqa+vx263oyhKh3aLxUJ9fb3PMp1OQafreB4RvV6n/WswyLC5zbl1Ee1CtS79+vm+N4xGAno/hGpd\nuiqc6tJtoTB69Gh27typ3b7pppsoKCjghRdeIDc3FwBVVc/7GBdqB4iPj0FR/IWCG4DYWCOxsTEB\n9Dw8xMYag92FXinU6nLXXfC//wsbN0L//vD88zoslsDfD6FWl+4SDnXp0UNSk5OTtaOPrFYrdrvd\np91ut5OWlobVakVVVex2OyaTSWtvaGggISHBZ52Ghma/I4XmZu8O66YmJ263/Aq3jV6vIzbW+FVd\n5IdobUK5Lk89BcXFYDKBwQAnT3Z+3VCuS1eEQl06++Wg20Lh5Zdfpr6+noKCAm3ZJ598wrBhwwDI\nyMjAZrMxa9YsADweD1VVVeTn5zNs2DAsFgs2m42kpCQAqquraWlpISMjw+fveDwqHk/HEUXbE+V2\ne2ht7ZtPWk+Suvj39bqUlETy4osRDBrkYeXKM3zrWxcevfZGbd+tWi/y5ybyevEvHOrSbRNkkZGR\nLFiwgDfffJPW1la2b9/OunXrmDdvHgCFhYWsX7+eiooKnE4ny5YtIzo6mokTJ6LT6ZgzZw7FxcUc\nO3YMu91OUVER06dPx2q1dlcXhTiv3bv1LFsWxWef6XjvPQMFBaE/VSDE1wU0UjAajSiKwtmzZwF4\n5ZVXUBQFh8PBlClTWLlyJfPnz+fo0aMMHjyYkpISpk6dCkBWVhbLly8nJyeHuro6xowZQ3l5OVFR\nUQAsXbqU06dPM3LkSNxuN5MnT+aZZ57p5s0V4psdPao7720hwoGidmbvbi9SV3fK73KH4zTJyYkc\nPnwck6nfJe5V72UweHc0njzZHPLD3kD4q8vx4wo33miivt4bBvn5LRQXh9e5n+T14l8o1MVq/eaf\nE5xLzn0kxFcSE1W2bXNQXm5g0CCVqVPl/E8i/EgoCHGOpCSVuXPPBrsbQgSNTJoKIUQfsWWLgR/8\nwMSMGUY+/bTjofndQUYKvcSpU1BZqWfwYJWrruqbc5ZCiJ5TXa3j3nuj8Xi8YXD8uJG333Z0+9+R\nkUIvcPIk3HRTDHfeaeL73zexbp2cxEwI4evTT3VaIAD84x898/EtodALvP56BJ9+6n0qVFWhpERO\nYiaE8HXttW7i49tnEW6+uWcOhJDpo16gXz/fo4L79+9TRwkLIS6BwYNVyssdvPRSBBaLSm5uzxwQ\nIaHQC9x6aytvvHGWV16JYOBAD7/5zZlgd0kI0QsNG6byyCMtPfo3JBR6Ab0efv/7Mzz55BmiosDP\nSWCFEOKSkFDoRaKjg90DIUS4kx3NQgghNBIKQgghNBIKQgghNBIKQgghNBIKQgghNBIKQgghNBIK\nQgghNAGHwrZt2xg8eDB33XVXh7adO3eSmZlJXFwcaWlpbNq0yae9pKSE1NRUzGYz48aNo6qqSmtz\nuVwUFBSQlJREQkICOTk5NDQ0XMQmCSGEuFgBhcLjjz/OAw88wJVXXtmhraamhqlTpzJv3jzq6upY\nuXIl+fn52gd/aWkpS5YsYcOGDdTW1pKdnU12djZOpxOAoqIi9u3bR0VFBQcPHsTj8ZCXl9cNmyiE\nEKKzAgoFo9HIe++9x/Dhwzu0bdy4kZSUFGbPnk1kZCQTJkxgypQprFmzBoDVq1eTl5dHRkYGUVFR\nLFy4EEVRKC0txe12s3btWhYtWkRiYiJms5ni4mLKysqoqanpni0VQghxQQGFwvz58+nf3//Fn202\nG+np6T7L0tPTqays9NuuKAqjRo2isrKSQ4cO0djYyOjRo7X2lJQUjEYjNpstkC4KIYTogm4795Hd\nbicpKclnWXx8PPX19Vq7xWLx226321EUpUO7xWLR1m+j0ynodB3PGKfX67R/DQbZf97m3LqIdlIX\n/6Qu/oVTXbr1hHiqev7rAHS1HSA+PgbFz2lE9Xo3ALGxRmJjYy74OOEmNtYY7C70SlIX/6Qu/oVD\nXbotFKxWK3a73WeZ3W4nISHhvO1paWlYrVZUVcVut2MymbT2hoYGbf32Zc1+RwrNzd4d1k1NTtxu\nvba8vh7eeUfPZZepjBkTftc+1ut1xMYav6pL+G3/N5G6+Cd18S8U6mKxdO7LcreFQkZGBuvWrfNZ\nVllZSWZmptZus9mYNWsWAB6Ph6qqKvLz8xk2bBgWiwWbzaZNQVVXV9PS0kJGRobPY3o8Kh5PxxFF\n2xPldntobfX+f22tQlaWiePHvUO+X/3qDHPn9szVinq7c+si2kld/JO6+BcOdem2CbKZM2dy+PBh\n1q5di8vlory8nK1btzJ37lwACgsLWb9+PRUVFTidTpYtW0Z0dDQTJ05Ep9MxZ84ciouLOXbsGHa7\nnaKiIqZPn47Var3oPpWVGbRAAHjuObn2sRBCnE9AIwWj0YiiKJw96/22/corr6AoCg6HA6vVSllZ\nGffddx/33nsvycnJbNy4kauvvhqArKwsli9fTk5ODnV1dYwZM4by8nKioqIAWLp0KadPn2bkyJG4\n3W4mT57MM88806WNs1h8RxTx8XLtYyFE76eqwbsCo6J2Zu9uL1JXd8rvcofjNMnJiRw+fByTqR8A\nHg8sWBDF5s0RDB2q8vzzTr7zndAe+n2dwaDDYonh5MnmkB/2BkLq4t/F1MXhgIoKPQMHqqSlhWYt\nL+Xr5Ve/imT16kjMZpVnnz3Dv/2bu1se12r1/3OCrwvp46t0Ovjtb10cO3aaiormsAsEIXra6dMw\naZKJO+4wMWFCDKtWRQS7S33anj16nnoqCpdLobZWR2Hhpb9Gb0iHwoV8/rnC//xPBDt26C98ZyFE\nB2+8YeDjj9vfP088ERXE3vR9DQ2+c0Zffnnp55C69XcKfcnRowo332zCbvfm4sMPu1iwoCXIvRKi\nb+nfXz3vbRGY73+/ldRUN/v3e4N27txL/5kUtiOFbdsMWiAAbNokw14hAnXjjW7uvrsFRVExm1We\nfPJMsLvUp/XrB+XlDl54wcmrrzr4+c8vfSiE7Uhh8GD1vLeFEBemKN79do8+6iIqKnhHzISSfv1g\n0qTWoP39sB0pZGe3Mn++C6vVw+jRbp56yhnsLgnRZ0VHSyCEirAdKQAsWtTCokWyH0EIIdqE7UhB\nCCFERxLWttk3AAAUgUlEQVQKQgghNBIKQggRApYtiyQpqR9paTHs2XPxv72SUBBCiD5u7149JSXt\nv4QuKLj4X0JLKAghRB/Xnb+EllAQQog+7oYbvL+EbjNnzsUfVRnWh6QKIUQoaPsl9FtvGbBY1C6d\nWVVC4Ru8/LKBzZsjuOwyD4sWuTCbg92j8OTxwGefKVgsqjwHQpxHv37eH+V2lYSCH+++q+fee6NR\nVe+8XF2djhdflF88X2otLXDXXUZ27zYQHe09t/zEicH7+b8Q4aBb9ynodDqMRiMmk0n79yc/+QkA\nO3fuJDMzk7i4ONLS0ti0aZPPuiUlJaSmpmI2mxk3bhxVVVXd2bWAfPSRTguEttvi0nv9dQO7d3u/\nt5w5o/CLX8hpmYXoad06UlAUhYMHD5KUlOSzvKamhqlTp/L0008zY8YM3nnnHaZMmUJqairp6emU\nlpayZMkStm3bRlpaGk8++STZ2dkcOnQIo9HYnV3slOuucxMRoXL2rDcYrr++e658JALj+do1kdzy\nNAjR47r1K7Cqqvi7uufGjRtJSUlh9uzZREZGMmHCBKZMmcKaNWsAWL16NXl5eWRkZBAVFcXChQtR\nFIXS0tLu7F6npaV5ePllJ7m5LfzXf7n47/+W0wEHw5QprWRmeqeLIiNVlixxBblHQoS+bp8XeeSR\nR/iXf/kXLBYLBQUFNDc3Y7PZSE9P97lfeno6lZWVAB3aFUVh1KhRWnswXHedm8cec3H//S1ERgat\nG2EtOhpeecXJ7t3NVFU1M3Wq7E8Qoqd16/TRddddx80338z69ev59NNPueOOO5g3bx52u73DlFJ8\nfDz19fUA2O12LBbLN7afS6dT0Ok6/jBDr9dp/xoMsg+gzbl16YsMBvjOdwCUr/7rHn29Lj1F6uJf\nONWlW0Nhz5492v+npKTw61//msmTJzNu3Di/00rnulB7m/j4GBQ/J27X670TzrGxRmJjYwLodXiI\njb30+2b6AqmLf8Guy6lTsHMnWK0wdmxQu+Ij2HW5FHr0kNTk5GTcbjc6nQ673e7TZrfbSUhIAMBq\ntfptT0tL6/CYDQ3NfkcKzc3eQ0abmpy43Rd/MqhQo9friI01flUXz4VX6KLaWoX+/VVMph7/U11y\nqevSV/SGujQ1QVaWkQMHvN/KH3mkhUceORuUvrTpDXXpKoulc1+Wuy0UPvjgAzZs2MBvfvMbbdkn\nn3xCdHQ0EydOZN26dT73r6ysJDMzE4CMjAxsNhuzZs0CwOPxUFVVxY9//OMOf8fjUfF4Oo4q2p4o\nt9tDa2vffNJ6Uk/X5exZyM01sn27gZgYleefdzJ+fO8/XKi76uJ2wz//qWA2w8CBff/SrsF8H23b\nZtACAWDVqggeeqh3HGTw9bq8+qqBNWsiGDBAZdkyF0lJff+577YJsoSEBFavXs1jjz1GS0sLBw8e\nZNGiRcydO5e7776bzz77jLVr1+JyuSgvL2fr1q3MnTsXgMLCQtavX09FRQVOp5Nly5YRHR3NpEmT\nuqt7ooe9/rqB7du93zGamxUeeeTiz9LY17S0QE6OkbFj+zFqVAyvviq/Ce2KuDjfD1azuXd+0H70\nkfdspO+9Z2Dr1gjy8kJjaqnbQiExMZHy8nJee+01Bg4cyPe+9z0mTpzIihUrsFqtlJWV8dRTT2E2\nm3nooYfYuHEjV199NQBZWVksX76cnJwcBgwYwI4dOygvLycqKrR+rOR2w6efKnz5ZbB70v3Onj3/\n7VBWVmbgnXe8QdDSIj+y66rx493MmdNCRIRKQoKHVat65yHh//iHDo+nfSp7//7Q2AmtqJ3dw9tL\n1NWd8rvc4ThNcnIihw8fx2Tqd4l7dWFOJ9x5p5F33zVgNKqsXu0kK6vnp1cMBh0WSwwnTzb36HRA\nczPcdpuJDz7QYzCorFx5hpyc3nsIaXfW5ZVXDMyd2/4tcfBgD3/7W7N2+7HHIlm/PoKEBJWnnz7D\nVVf13unNS/V66QyPB3S95HPWX10+/1zhhhtiaGz0BsPEiWdZt653BhiA1dq/U/eTce4l8sorBt59\n11tup1Nh0aJosrKaL7BW3xETA6WlDj75RMfAgWpIzK121qRJrVx/fSvvvGMgMlLlV79qn/9+6y09\nv/mNd+Rw4gQUFkbz9tuOYHW1T+ktgfBNhg5VKStzsHmzgfh4lR//ODSGxxIKl8i551Ly3g5SR3pQ\nVBSMHt17vwX3lMhIeOklp7ajecCA9ie3psb3ea+p6eWfdCIgKSkefvGLi792QW8kr9BL5LbbzjJm\njHe6KCpKZfHizh9N8etfRzJyZAy33GLi0KHu+wGX6D46HQwfrvoEAsD3v+9m8OD2oLzrrtD4Nim6\n5oUXIhg71sSkSaZety9CRgqXiMkEr73m4NAhHQMGqJ0+bHHbNj2//a13+uGLL2D+fCNbt8r0Q18x\naJDK9u0Otm0zkJDg4Qc/6P2H6YqeVVWl8zk67557otm7t/e8pyUULiGDwTvcDMTnn+u+drt7Rgou\nF3z8sQ6rNbzm/4Nh0CCVH/5QRgjC6+hR3/f0kSO9a6TQu3ojOrj55lYGDGgPkjvv7PqHy+nTMGmS\niR/8IIbMzBheekm+GwhxPrW1CqdPd89jjR3rZtCg9vf0rbf2rqP05NOgl7vsMu/0wxtvGEhM7J7p\nh9dei+Bvf/OeCqS1VWH58ihuv713vTCF6A08HigoiObVVyOIjoZnn9UzcWLXDqawWlW2bXPw2mve\n6yn3tveehEIfcNllKvfc033TD5GRvtNFERHd9tBChJQdO/S8+qr3DXLmDDz4YBQTJ3b9vZiYqFJY\n2DunFGX6KAzdemsrN97o/XYSE6Py6193/gc327bpycuLpqgoiqamnuqhEL1DS4vytdtB6sglJCOF\nMBQRAZs2OamtVYiNVensFU/37dORm2vE7fa+UY4e1fHii84e7KkQwXXTTa2MHdvK3r0GdDr45S9D\nPxUkFMLYoEGBHXX0wQd6LRAAbDYZaIrQFhkJL7/s5O9/N5CcbMRiaaW1d+0C6HbyrhadlpHhJiKi\nPUi++1055l6EPoPB+0v9YcOC3ZNLQ0YKotPS0jxs2uTkpZciGDLEwwMPhP5QWvQ+1dU6/vQn7zUM\n5s5t6fT051//qmf9eu96Dz/s4mtXABZfkVC4CLW1Cm+9pWfoUJXvfS+8vi3/+7+7+fd/D69tFr3H\nZ58pTJli4vRp7zTmvn06/ud/LnygxD/+oTBjhhGXy7vegQM6/vxn2R/mj0wfBeiLLxRuvNHEffcZ\nmTbNRElJZLC7JELU5s0GsrON3HNPNMePyzmvACor9VogAOza1bnvtR9+qNcCoe1xvi7U9xV0loRC\ngMrLDdTWtpdt3Trfg/z/+lc9Dz4YxW9+E8mZ3ntq9V7l+HGFF16IYOtWGbi2sdl03Hef96peZWUR\n5OcH96peb76p54c/jOahh6Kw24MXUCkpHvT69v1aV1/duR+SjRzpJjq6fb1rr/Ud7f7nf0Zx2WX9\nSE2NYffu8L7Gu7wLA5SQ4HvEjtXafvujj3TccYeRs2e9b5pPP9XxzDOSDOfTNvKqr/cG7fz5LhYt\n6tl9FS6X90P3iisgPr5H/9RFO3hQ53O69WCeSfPjj3X88IdGWlvbX9evvNL5qZfaWoWYGJV+3XDt\nq7Q0D6tXn+GFFyIYOFBlyRLfsw3/9reR2jWTn376DCNHekPjW99S+dOfnNo+hQUL2tfbtUvP8897\nR/wNDQr33RfNhx+GzrVOAiUjhQBlZ7eSn99C//4q3/62myefbP/Qr6jQa4EAsGdPeH/jONeJEwp/\n/KOBt9/2rcmbbxq0QAD44x979ufVzc2QnW3illuMXHklvPhi7/xedN11bvr1a//CcdNNwZvbqK7W\naYEAsG9f517Xbrf3DKBpaf246qp+/OUv3VPryZNb+fOfnaxefYYhQ9prtHevnl//Oor6eh0HDuiZ\nM8d3dHXddW5+97szLFvmwmxuX37qlO/Ip+1KauEqZELBfYn2fSoKFBe7OHToNG+/7SA1tX34mpbm\nQVHaX6QjR/Z8p2pqFDZsiGDnzt4bQLW1CjfdZOL++43cfruJxx9v3w8zdKjv8D8xsfvO2OpydbyY\nUVmZgQ8/9NbK44FHHw3uOT48HvjnPxVOnvRdnpzsvarX/Pkuli49Q0mJ74jzyy/h1VcNAX/xaGiA\nl14y8NZb/tdraelYs2uu8Z16GTu2c6/rbdu8U18AZ84oPPyw77WrX3vNwG23GcnPj+aLL7r+Qfz1\nCxrV1nbuMcePb+Wqq9q36f77w/uout75NSkAbjfMnx/Nn//s/XD56CMdmZnB6UtmppvVq8+weXME\nQ4d6+PnPO38hnYvRNvVSV+fN9gULXDz8cO97Qb/xhoEvvvDdD7Nwobef48e7WbDAxcaNEQwerPLU\nU12fbvN44P77o9m8OYL4eA8vvHCG667zvumjo33vGxXl5wEukbNnYeZMI2+9ZSAqSuWZZ84weXL7\niOCqqzx+p9JOnoSsrBgOHw7seW9ogJtvjtFO1fzAAy6KirzrqSosWBDFhg0RWCywdq2OsWPbp15e\nesnJH/5gwGpV+clPOvca+/qO23NH0R99pGPu3GjtwvfHjum6fJ2QG25o5fLLPdr23X13584t1K8f\n/OUvDt59V098vEp6evhdPfBciqr2nQtDqqrKZ599gaK0v7heeUXPT38aDTQBSXz725+xdWt4nOFt\nwwYD//Vf7Z9qgwZ5qKjwnevV63XExhppanLidgfnxb59u578/PZP429/283WrT23r2XrVj2Fhe1/\n71//1cOuXd66eDxQWBjFtm0GYmLg2WddXH99cKZmSkv13Hdfez8TEjy8996F5+r//GcDDz7Y/rzH\nxal8+OGFP1C/vl7//ioffeRdb8cOPT/6UXtfLrvMw1//2rVDNltaYObMaCor9SiKSnFxC3fd5a11\n+/vWy2hU+fvfu36hmYYG2L7dwIABKjfe2H0j9d7wPuoqiyWG/v37+3x++tOnQqGpqYm4uLhgd0MI\nIfqkxsZGYmNjz3ufPhUK/kYKdjtMnWrk2LHTQBKLFn3KPfdEf/ODhJjlyyPYvDmCxEQPK1e6uOIK\n36czFL7hBOrUKbj1ViOHDukAlV/+soW8PN/RQG+oS0sLzJ4dzbvv6jEYVH77WxdTpnTu2+2yZZFs\n2OCdznnqKRejRnVuG1asiODFFyOwWlWefNLFiBHe9U6fhunTjRw44J16+fnPz5Kf37NTkX//u46X\nXzYwcKD31PDBnMq7kN7weumqkBwpANTVneqw7MsvYfduJz/+8SAOHz6OydQNx76FCINBh8USw8mT\nzbS29s0X88U4fRree09PQoLKd77Tcbt7S13cbjh0SIfZrHY43PlSa26G99+PYPjwaJKTw+v1ciG9\n5fXSFVZr/07dr8/vaAYwm+nW+UPR9/Xr592J3dvp9XDllb3jQyYmBiZMcGOx0OFoKBE++txI4Zu0\n7W/ozJyZEEII/0ImFFRV5dSpU52aMxNCCOFfyISCEEKIrguZXzQLIYToOgkFIYQQGgkFIYQQGgkF\nIYQQGgkFIYQQGgmFELJt2zYGDx7MXXfd1aFt586dZGZmEhcXR1paGps2bQpCD4PjyJEjTJs2jYED\nBzJkyBDy8vJoamoCwrsuH374ITfeeCNms5khQ4Zw5513cuLECSC869Lmpz/9KTpd+0dk2NREFSHh\nscceU1NTU9Xrr79enTFjhk/bF198ofbr109dt26d6nK51DfffFM1mUyqzWYLUm8vrREjRqg/+tGP\nVIfDoX7++efqmDFj1Pz8/LCui8vlUgcNGqQWFxerLS0tan19vXrDDTeo06ZNC+u6tNm3b586YMAA\nVafTqaqqqsePHw+bmugXL168ONjBJLpu3759rFq1ioqKCtxuN9OnT9fafv/73/PFF1/wzDPPoNfr\nGTZsGNXV1fz9739n0qRJQex1z2tsbKS6uprly5cTFxdH//79aWlp4aWXXiIiIiJs63Lq1CkSExOZ\nP38+BoMBk8nEqVOnKCsrIzIyMmzrAt4fwk6bNo28vDx27tzJL3/5S37/+99TU1MTFjWR6aMQMX/+\nfPr393/CK5vNRnp6us+y9PR0KisrL0XXgiouLo41a9ZgtVq1ZUePHmXo0KFhXRez2cw999yjTY8c\nOHCAdevWcccdd4R1XQCeffZZjEajzzRsVVVV2NREQiEM2O12LBaLz7L4+Hjq6+uD1KPgef/993n6\n6af5+c9/LnXBu78lKiqKq6++mszMTBYvXhzWdamtrWXx4sX87ne/81keTjWRUAgTqpzNhD179pCV\nlcWKFSsYP348IHW5/PLLcblcHDhwgAMHDjBr1iwgfOvy0EMP8aMf/YiUlJQObeFSEwmFMGC1WrHb\n7T7L7HY7CQkJQerRpVdaWsqkSZMoKSnh3nvvBaQu5xo+fDjFxcX84Q9/IDIyMizrsmPHDvbu3csv\nfvELwDcEwum1IqEQBjIyMrDZbD7LKisryczMDFKPLq29e/eSm5vLli1bmDlzprY8nOuya9cuUlNT\nfZYpioKiKFx77bW8//77Pm3hUJeNGzdy4sQJLr/8cqxWK9dccw2qqpKQkEBaWlr41CSIRz6JHpCb\nm9vhkNQTJ06ocXFx6vPPP6+eOXNG/ctf/qLGxMSo1dXVQerlpdPa2qpeddVV6nPPPdehLZzr0tjY\nqA4ZMkR9+OGHVYfDoZ44cUK95ZZb1BtuuEGtq6sLy7p8+eWX6ueff67993//93+qoijq8ePH1SNH\njoRNTSQUQkR0dLRqNBpVg8GgGgwG7Xabd955Rx01apQaHR2tpqamqq+++moQe3vpvPPOO6pOp1ON\nRqNWk7Z/jxw5ErZ1UVVVra6uVm+44QY1JiZGHTRokDpjxgz1+PHjqqqG7+vlXIcPH9Z+p6Cq4VMT\nuZ6CEEIIjexTEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQEEIIoZFQ\nEEIIoZFQEEIIoZFQEEIIofn/BekoDYZ3hWcAAAAASUVORK5CYII=\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# 異常度aを求める\n", "lam = 1.0\n", "coefs = clf.coef_\n", "score = clf.score(X, y)\n", "sig2 = (lam*sum(coefs^2) + score)/N\n", "y_bar = uscrime['y'].mean()\n", "X_bar = uscrime[vars].mean().values\n", "a = (y - y_bar - coefs.dot((X - X_bar).T))^2/sig2 \n", "list_plot(a, figsize=4, zorder=2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\t\n", "\t\t計算結果は、井出本とは若干異なりますが、概ねデータの特徴を捉えていると思われます。\n", "\t
\n", "\n" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "| M | Ed | Po1 | Po2 | LF | M.F | Pop | NW | U1 | U2 | GDP | Ineq | Prob | Time |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 151.0 | 91.0 | 58.0 | 56.0 | 510.0 | 950.0 | 33.0 | 301.0 | 108.0 | 41.0 | 394.0 | 261.0 | 0.0846 | 26.2011 |
| 140.0 | 93.0 | 55.0 | 54.0 | 535.0 | 1045.0 | 6.0 | 20.0 | 135.0 | 40.0 | 453.0 | 200.0 | 0.042 | 21.7998 |
| 125.0 | 108.0 | 113.0 | 105.0 | 567.0 | 985.0 | 78.0 | 94.0 | 130.0 | 58.0 | 626.0 | 166.0 | 0.0348 | 26.401 |
| 131.0 | 121.0 | 160.0 | 143.0 | 631.0 | 1071.0 | 3.0 | 77.0 | 102.0 | 41.0 | 674.0 | 152.0 | 0.0417 | 22.1005 |
| 152.0 | 112.0 | 82.0 | 76.0 | 571.0 | 1018.0 | 10.0 | 79.0 | 103.0 | 28.0 | 537.0 | 215.0 | 0.0382 | 25.8006 |
\n", "\t\t参考のために、λを0から5まで0.1刻みで計算したscoreの値をプロットしてみました。\n", "\t\t1.0でも概ね収束しているようにみえます。\n", "\t
\n", "" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAECCAYAAAARlssoAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAPYQAAD2EBqD+naQAAIABJREFUeJzt3Xt4U1W+PvB3JyltkjbQ0mBopbSgFurI1HLnh6IMlBEK\nShkVRuBULkMRZ0AU6SPKoIyjRa7KOSK3Qai3ES3ggQGt6BwEGYSCyICtZaYWWmhKWilNSi/J+v0R\nusImZbQQqHS/n+fxiVl775X1TULe7LV3dhUhhAAREREAXXMPgIiIfj4YCkREJDEUiIhIYigQEZHE\nUCAiIomhQEREEkOBiIgkhgIREUktMhSEEKisrAR/l0dE1DSG5h7A1SorO+fX5nSeQ1xcNP7972KY\nzWHNMKrmpdMpiIgwo7zcCY9HW8Go5doBbdev5dobWK1X/3nXIvcU9uzRAwAOHNA380iah06nQFEU\n6HRKcw/lutNy7YC269dy7YHU4kLhww8NeOQRIwDgwQdD8Mkn2gwGIqIr0eJCYfNm34yYEAq2bAlq\nxtEQEd1YWlwodOwoLrnvaaaREBHdeFpcKMyeXYORI+sBAA8/XIff/762mUdERHTjaHGhYDYDS5bU\nAAAyM2sRHNzMAyIiuoG0uFAgIqIrx1AgIiKJoUBERBJDgYiIJIYCERFJDAUiIpIYCkREJDEUiIhI\nYigQEZHEUCAiIomhQERE0nUNhaKiIqSkpCAyMhJxcXHIyMhodD0hBP74xz8iLi4OFosFiYmJ+Otf\n/3o9h0pEpEnXNRRSU1PRoUMHFBYWIicnB9nZ2Vi6dKnfeq+//jrWrl2LTz75BGfPnsWLL76IsWPH\n4siRI9dzuEREmnPdQmH//v04fPgwMjMzERoais6dO2PmzJlYuXKl37q5ubno378/brnlFiiKgmHD\nhqFt27Y4fPjw9RouEZEmXbdQyM3NRWxsLCwWi2xLSkpCXl4enE6nat1hw4bh888/x9dff426ujps\n2bIF1dXVGDBgwPUaLhGRJhl+fJXAcDgcCA8PV7VFREQAAM6cOQOz2SzbR44ciUOHDuHOO++Eoigw\nmUxYv349oqOjr9dwiYg06bqFAuA9gPxTbNiwAevXr8f+/fvxi1/8Ajk5Ofjtb3+LmJgYdO/eXbWu\nTqdAp1NUbXq9Tt4aDNo7weri+rVGy7UD2q5fy7UH0nULBavVCofDoWpzOBxQFAVWq1XVvnz5ckyZ\nMgVJSUkAgKFDh2LgwIHYsGGDXyhERJihKJeGghsAYLEYYbGYoVUWi7G5h9BstFw7oO36tVx7IFy3\nUOjRoweKiopQXl4up4327duHhIQEmEwm1bputxtut1vVVlNT02i/5eVOvz0Fp7MaAFBZWQ23Wx+o\nEm4Yer0OFovxQv2e5h7OdaXl2gFt16/l2huEh1/9l+DrFgqJiYno2bMnMjIysGjRIhQXF2PJkiWY\nNWsWAKBLly5Yu3Yt+vXrhxEjRmD16tUYMWIEEhIS8Omnn2Lnzp14+umn/fr1eAQ8HvW0VMMbwu32\noL5em28OQNv1a7l2QNv1a7n2QLiuxxQ2btyIyZMnw2azoXXr1pg6dSrS09MBAN999x2qqqoAAM88\n8wzcbjceeOABlJWVITY2FqtXr+bZR0RE15gifurR35+psrJzfm0uVxViY6NQWFgCkym0GUbVvAwG\nHcLDzaiocGruG5OWawe0Xb+Wa29gtYZddR88TE9ERBJDgYiIJIYCERFJDAUiIpIYCkREJDEUiIhI\nYigQEZHEUCAiIomhQEREEkOBiIgkhgIREUkMBSIikhgKREQkMRSIiEhiKBARkcRQICIiiaFAREQS\nQ4GIiCSGAhERSQwFIiKSGApERCQ1ORSKioqQkpKCyMhIxMXFISMjo9H1hgwZAqPRCJPJBJPJBKPR\nCIPBgPnz5wMA7HY7xo4di/bt2yMiIgKTJk3C+fPnfQPT6eT2DbfTp0+/wjKJiOinaHIopKamokOH\nDigsLEROTg6ys7OxdOlSv/V27NiB6upquFwuuFwunD59GjabDaNGjQIAjBkzBg6HA9988w0KCgpQ\nUlKCp556Sm6vKAry8/PhcrlkP8uWLbuKUomI6Mc0KRT279+Pw4cPIzMzE6GhoejcuTNmzpyJlStX\n/ui2c+bMwciRI5GQkACn04nPP/8cc+fORWRkJCIiIrBo0SKsX78e9fX1AAAhBIQQV1YVERFdkSaF\nQm5uLmJjY2GxWGRbUlIS8vLy4HQ6L7tdQUEBsrKyMG/evMuu06ZNG1RVVeH48eOybfbs2ejYsSMi\nIiIwZcqU//gYRER09ZoUCg6HA+Hh4aq2iIgIAMCZM2cuu11mZiYmTJiAtm3bAgDMZjMGDBiA559/\nHmVlZaioqMC8efMQFBSE8vJyAEDfvn2RnJyMgoICfPnll9i7dy+mTZvWpOKIiKhpDE3doKlTOhUV\nFdiwYQPy8/NV7Rs2bMDjjz+O+Ph4WK1WvPDCC3jrrbdgMHiHtHv3brlufHw8MjMzMWLECKxatQpB\nQUFymU6nQKdTVH3r9Tp5azBo7wSri+vXGi3XDmi7fi3XHkhNCgWr1QqHw6FqczgcUBQFVqu10W02\nbdqE+Ph4xMTEqNqjo6ORnZ0t75eXl8PlciE6OrrRfmJjY+F2u2G321XrRESYoSiXhoIbAGCxGGGx\nmH96gS2MxWJs7iE0Gy3XDmi7fi3XHghNCoUePXqgqKgI5eXlctpo3759SEhIgMlkanSbLVu2IDk5\n2a9927Zt6NSpE7p06QLAe7ZSx44dERUVhUOHDiErKwsLFy6U6x89ehTBwcGIiopS9VNe7vTbU3A6\nqwEAlZXVcLv1TSmxRdDrdbBYjBfq9zT3cK4rLdcOaLt+LdfeIDz86r8ENykUEhMT0bNnT2RkZGDR\nokUoLi7GkiVLMGvWLABA165dsWbNGvTr109uc/DgQQwePNivr/fffx8nTpxAdnY2ysrK8Nxzz8lT\nUtu1a4eVK1eiXbt2mDFjBgoLCzF37lxMmTLFb6/A4xHweNRTWg1vCLfbg/p6bb45AG3Xr+XaAW3X\nr+XaA6HJk28bN25EcXExbDYbBg4ciLS0NKSnpwMA8vPzUVVVpVq/tLQUNpvNr5/FixfDZDIhOjoa\n/fv3R1pamjyQHBUVhW3btmHz5s2IjIxE//79MXToUGRmZl5JjURE9BMp4gb/MUBZ2Tm/NperCrGx\nUSgsLIHJFNoMo2peBoMO4eFmVFQ4NfeNScu1A9quX8u1N7Baw666Dx6mJyIiiaFAREQSQ4GIiCSG\nAhERSQwFIiKSGApERCQxFIiISGIoEBGRxFAgIiKJoUBERBJDgYiIJIYCERFJDAUiIpIYCkREJDEU\niIhIYigQEZHEUCAiIomhQEREEkOBiIgkhgIREUkMBSIikhgKREQkBTQUioqKkJKSgsjISMTFxSEj\nI6PR9YYMGQKj0QiTyQSTyQSj0QiDwYD58+cDAOx2O8aOHYv27dsjIiICkyZNQk1NTSCHSkREjQho\nKKSmpqJDhw4oLCxETk4OsrOzsXTpUr/1duzYgerqarhcLrhcLpw+fRo2mw2jRo0CAIwZMwYOhwPf\nfPMNCgoKUFJSgieffDKQQyUiokYoQggRiI7279+Pfv364cyZM7BYLACAN954A8uWLcPRo0f/47aP\nP/44FEXBa6+9BqfTCYvFgi+++AJ9+/YFABw7dgy9e/dGeXk5DAaDatuysnN+/blcVYiNjUJhYQlM\nptBAlHdDMRh0CA83o6LCifp6T3MP57rScu2AtuvXcu0NrNawq+7D8OOr/DS5ubmIjY2VgQAASUlJ\nyMvLg9PphNlsbnS7goICZGVl4fjx45ftu02bNqiqqsLx48cRHx8fqCETEdElAjZ95HA4EB4ermqL\niIgAAJw5c+ay22VmZmLChAlo27YtAMBsNmPAgAF4/vnnUVZWhoqKCsybNw9BQUEoLy8P1HCJiKgR\nAdtTAICmzkRVVFRgw4YNyM/PV7Vv2LABjz/+OOLj42G1WvHCCy/grbfe8ps6AgCdToFOp6ja9Hqd\nvDUYtHeC1cX1a42Wawe0Xb+Waw+kgIWC1WqFw+FQtTkcDiiKAqvV2ug2mzZtQnx8PGJiYlTt0dHR\nyM7OlvfLy8vhcrkQHR3t10dEhBmKcmkouAEAFosRFkvj01ZaYLEYm3sIzUbLtQParl/LtQdCwEKh\nR48eKCoqQnl5uZw22rdvHxISEmAymRrdZsuWLUhOTvZr37ZtGzp16oQuXboA8J6t1LFjR0RFRfmt\nW17u9NtTcDqrAQCVldVwu/VXVdeNSK/XwWIxXqhfWwfctFw7oO36tVx7g/Dwq/8SHLBQSExMRM+e\nPZGRkYFFixahuLgYS5YswaxZswAAXbt2xZo1a9CvXz+5zcGDBzF48GC/vt5//32cOHEC2dnZKCsr\nw3PPPYennnqq0cf1eAQ8HvW0VcMbwu32aPYsBEDb9Wu5dkDb9Wu59kAI6OTbxo0bUVxcDJvNhoED\nByItLQ3p6ekAgPz8fFRVVanWLy0thc1m8+tn8eLFMJlMiI6ORv/+/ZGWloZp06YFcqhERNSIgB5o\njoqKwtatWxtd5na7/dqqq6sbXTc8PBxbtmwJ5NCIiOgn4GF6IiKSGApERCQxFIiISGIoEBGRxFAg\nIiKJoUBERBJDgYiIJIYCERFJDAUiIpIYCkREJDEUiIhIYigQEZHEUCAiIomhQEREEkOBiIgkhgIR\nEUkMBSIikhgKREQkMRSIiEhiKBARkcRQICIiiaFARERSk0OhqKgIKSkpiIyMRFxcHDIyMhpdb8iQ\nITAajTCZTDCZTDAajTAYDJg/fz4AwOFwYNy4cbDZbGjbti0GDRqEgwcP+gam08ntG26nT59+hWUS\nEdFP0eRQSE1NRYcOHVBYWIicnBxkZ2dj6dKlfuvt2LED1dXVcLlccLlcOH36NGw2G0aNGgUAmDp1\nKsrKyvDtt9/i9OnT6N27N4YOHQohBABAURTk5+fD5XLJfpYtW3aV5RIR0X/SpFDYv38/Dh8+jMzM\nTISGhqJz586YOXMmVq5c+aPbzpkzByNHjkRCQgIAIDc3FyNHjkSbNm0QFBSE8ePHw26349SpUwAA\nIYQMCCIiuj6aFAq5ubmIjY2FxWKRbUlJScjLy4PT6bzsdgUFBcjKysK8efNk2/Dhw/HOO+/g9OnT\ncDqdWLduHe68805ERUXJdWbPno2OHTsiIiICU6ZM+Y+PQUREV8/QlJUdDgfCw8NVbREREQCAM2fO\nwGw2N7pdZmYmJkyYgLZt28q2BQsWYNiwYYiKioKiKOjYsSO2b98ul/ft2xfJyclYv349/vWvf+Gh\nhx7CtGnTsG7dOlXfOp0CnU5Rten1OnlrMGjvWPrF9WuNlmsHtF2/lmsPpCaFAoAmT+lUVFRgw4YN\nyM/PV7VPnToViqLg5MmTsFgsWLZsGQYPHoxjx47BZDJh9+7dct34+HhkZmZixIgRWLVqFYKCguSy\niAgzFOXSUHADACwWIyyWxoNKCywWY3MPodlouXZA2/VrufZAaFIoWK1WOBwOVZvD4YCiKLBarY1u\ns2nTJsTHxyMmJka2uVwu/OUvf8GePXvkdNGcOXOwePFifPzxx3jggQf8+omNjYXb7Ybdbkd0dLRs\nLy93+u0pOJ3VAIDKymq43fqmlNgi6PU6WCzGC/V7mns415WWawe0Xb+Wa28QHn71X4KbFAo9evRA\nUVERysvL5bTRvn37kJCQAJPJ1Og2W7ZsQXJysqrN7XZDCIH6+nrZ5vF4UFdXBwA4dOgQsrKysHDh\nQrn86NGjCA4OVh1z8G4n4PGo914a3hButwf19dp8cwDarl/LtQParl/LtQdCkybfEhMT0bNnT2Rk\nZODcuXP49ttvsWTJEjz22GMAgK5du2LPnj2qbQ4ePIi4uDhVW1hYGO6991786U9/gt1ux/nz5/HS\nSy+hVatWGDBgANq1a4eVK1diwYIFqK2tRX5+PubOnYspU6b4TRUREVHgNPmIzMaNG1FcXAybzYaB\nAwciLS0N6enpAID8/HxUVVWp1i8tLYXNZvPr591334XVakViYiKio6ORk5OD7du3Izw8HFFRUdi2\nbRs2b96MyMhI9O/fH0OHDkVmZuYVlklERD+FIm7wHwOUlZ3za3O5qhAbG4XCwhKYTKHNMKrmZTDo\nEB5uRkWFU3O70VquHdB2/VquvYHVGnbVffDcLSIikhgKREQkMRSIiEhiKBARkcRQICIiiaFAREQS\nQ4GIiCSGAhERSQwFIiKSGApERCQxFIiISGIoEBGRxFAgIiKJoUBERBJDgYiIJIYCERFJDAUiIpIY\nCkREJDEUiIhIYigQEZHEUCAiIimgoVBUVISUlBRERkYiLi4OGRkZja43ZMgQGI1GmEwmmEwmGI1G\nGAwGzJ8/HwDgcDgwbtw42Gw2tG3bFoMGDcLBgwcDOVQiImpEQEMhNTUVHTp0QGFhIXJycpCdnY2l\nS5f6rbdjxw5UV1fD5XLB5XLh9OnTsNlsGDVqFABg6tSpKCsrw7fffovTp0+jd+/eGDp0KIQQgRwu\nERFdImChsH//fhw+fBiZmZkIDQ1F586dMXPmTKxcufJHt50zZw5GjhyJhIQEAEBubi5GjhyJNm3a\nICgoCOPHj4fdbsepU6cCNVwiImqEIVAd5ebmIjY2FhaLRbYlJSUhLy8PTqcTZrO50e0KCgqQlZWF\n48ePy7bhw4fjnXfewf3334+wsDCsW7cOd955J6KiogI1XCIiakTAQsHhcCA8PFzVFhERAQA4c+bM\nZUMhMzMTEyZMQNu2bWXbggULMGzYMERFRUFRFHTs2BHbt28P1FCJiOgyAhYKAJo8519RUYENGzYg\nPz9f1T516lQoioKTJ0/CYrFg2bJlGDx4MI4dOwaTyaRaV6dToNMpqja9XidvDQbtnWB1cf1ao+Xa\nAW3Xr+XaAylgoWC1WuFwOFRtDocDiqLAarU2us2mTZsQHx+PmJgY2eZyufCXv/wFe/bskdNFc+bM\nweLFi/Hxxx/jgQceUPUREWGGolwaCm4AgMVihMXS+B6KFlgsxuYeQrPRcu2AtuvXcu2BELBQ6NGj\nB4qKilBeXi6njfbt24eEhAS/b/cNtmzZguTkZFWb2+2GEAL19fWyzePxoK6urtE+ysudfnsKTmc1\nAKCyshput/6Ka7pR6fU6WCzGC/V7mns415WWawe0Xb+Wa28QHn71X4IDFgqJiYno2bMnMjIysGjR\nIhQXF2PJkiWYNWsWAKBr165Ys2YN+vXrJ7c5ePAgBg8erOonLCwM9957L/70pz/hzTffhMViwaJF\ni9CqVSsMGDDA73E9HgGPRz1t1fCGcLs9qK/X5psD0Hb9Wq4d0Hb9Wq49EAI6+bZx40YUFxfDZrNh\n4MCBSEtLQ3p6OgAgPz8fVVVVqvVLS0ths9n8+nn33XdhtVqRmJiI6Oho5OTkYPv27X4HsomIKLAC\neqA5KioKW7dubXSZ2+32a6uurm50XavVivXr1wdyaERE9BPwMD0REUkMBSIikhgKREQkMRSIiEhi\nKBARkcRQICIiKaCnpP7c5eXp8N57BkRECEyaVIeQkOYeERHRz4tmQqG4WEFKiglnz3ovifHVV3q8\n+eb5Zh4VEdHPi2amj/bt08tAAIBPP9VMHhIR/WSaCYVbbvFAr/ddIyk+ntdGISK6lGZC4Y47PFix\n4jz69KnHfffVYd26xi+xQUSkZZqaQ7n//nrcf3/9j69IRKRRmtlTICKiH8dQICIiiaFAREQSQ4GI\niCSGAhERSQwFIiKSGApERCQxFIiISGIoALDbFcycGYwJE0Kwe7e+uYdDRNRsNPWL5ssZP96I3Fxv\nGOTkGPD550506iR+ZCsiopanyXsKRUVFSElJQWRkJOLi4pCRkdHoekOGDIHRaITJZILJZILRaITB\nYMD8+fMBQLWsYblOp8OuXbu8A9Pp5DoNt9OnT7+KUhvn8QAHD/qehvPnFRw9yr0FItKmJu8ppKam\nomfPnnj33XdRWlqKoUOHwmazYcaMGar1duzYobp/9uxZ3H777Rg1ahQAoLpafUG6PXv2YNy4cejV\nqxcAQFEU5Ofno0OHDk0dYpPodEDv3m7s3et9KkwmgW7d3Nf0MYmIfq6atKewf/9+HD58GJmZmQgN\nDUXnzp0xc+ZMrFy58ke3nTNnDkaOHImEhAS/ZR6PB9OmTcMrr7yC4OBgAIAQAkJcnymc9eurkZ5e\nizFj6vDhhy7ExHDqiIi0qUl7Crm5uYiNjYXFYpFtSUlJyMvLg9PphNlsbnS7goICZGVl4fjx440u\nf/PNNxESEoLU1FRV++zZs7Fnzx6cO3cODz74IBYvXnzZx7gabdoAL7xQ49cuBPDGG0HYt0+PHj3c\nmDq1DorSSAdERC1Ek/YUHA4HwsPDVW0REREAgDNnzlx2u8zMTEyYMAFt27b1WyaEQGZmJp555hlV\ne9++fZGcnIyCggJ8+eWX2Lt3L6ZNm9aU4V61VauCMHduCP73f4Mwb14I3ngj6Lo+PhHR9dbkYwpN\nndKpqKjAhg0bkJ+f3+jyrVu3oq6uDsOHD1e17969W/5/fHw8MjMzMWLECKxatQpBQb4PZ51OgU6n\n/vqu1+vkrcFw5WfdHjhg8LtvMPiON1RUACEhgNF4xQ9xTVxcv9ZouXZA2/VrufZAalIoWK1WOBwO\nVZvD4YCiKLBarY1us2nTJsTHxyMmJqbR5Rs3bkRKSsqPPnZsbCzcbjfsdjuio6Nle0SEGYpyaSh4\nP7gtFiMsliufbho4EMjOvvi+AeHh3qds8mRg9WogOBh4803g4Yev+GGuGYvlZ5ZW15GWawe0Xb+W\naw+EJoVCjx49UFRUhPLycjlttG/fPiQkJMBkMjW6zZYtW5CcnHzZPj/66COsX79e1Xbo0CFkZWVh\n4cKFsu3o0aMIDg5GVFSUat3ycqffnoLT6T2zqbKyGm73lZ9e+sgjQHW1Af/4hx69erkxdmw9KiqA\nXbt0WL3a+8arqQEmTRIYPNgFRQHq6oD16w04c0bBgw/WN8vvHfR6HSwW44X6tfW3qLVcO6Dt+rVc\ne4Pw8Ks/5tqkUEhMTETPnj2RkZGBRYsWobi4GEuWLMGsWbMAAF27dsWaNWvQr18/uc3BgwcxePDg\nRvsrLCxERUUF4uLiVO3t2rXDypUr0a5dO8yYMQOFhYWYO3cupkyZ4rdX4PEIeDzqD96GN4Tb7UF9\n/dW9OdLSapGW1tCv99blUo+hthaoq/NApwOmTg3Bhx96p7dWrTLgs89caN/eO77cXB3+8Q897rzT\ngz59rv1pr4Go/0al5doBbdev5doDocmTbxs3bkRxcTFsNhsGDhyItLQ0pKenAwDy8/NRVVWlWr+0\ntBQ2m63Rvk6fPg1FUfyWR0VFYdu2bdi8eTMiIyPRv39/DB06FJmZmU0d7jVx991u3HWX7289z55d\nC92FZ3LbNl/OlpfrsGePd0/l88/1GDbMhD/+MQT332/E5s2+9fbv12Hy5BA88UQwTp9WB47LBZSW\n8pQnIro+mnygOSoqClu3bm10mdvt/+330h+pXaxPnz6NbgMA/fv3Vx1s/jkJCgLee68ahw7p0Lo1\ncOutvm8lnTt78M9/eoNAUQQ6d/Yuy84Ogtvt/XAXQsEHHxhw//31KClR8OCDJjid3mW5uXr8/e8u\nAMAnn+gxebIRLpeCX/+6DmvXnofhwiv29tsG/P3vBvziFx489lgt9BdmyerrgS++0KFdO6CRn4RA\nCPC0WiK6LB6mv0IGA9Cjh0cVCACwdm017r23Ht26ubF06XkkJnqX33yzer2bb/ZOKX37rU4GAgAc\nO6aHy5sJePrpEDlVtX17kNwL2bjRgBkzjMjODsL8+cFYvLgVAG8gjB5tRGqqEf37A7Nnt5L9fv21\nDt27mxEdHYqZM4Nx8Ulky5cHITXViDlzgnFxhpeWKsjMbIUlS1qhslJd/759Onz4oQF2uzphqquB\nb77R4cwZ/+Q5f977HxH9fDEUAiwuTuC996qRk+PCmDG+Kabf/74WDz1Uhw4dPBg+vA7PPOP9sVxC\nggcWi+8Tuls3NxqO2dfWqvuuufD7un/8Q33wvOF+bq4O//d/vp2/VauCcO6c9/+nTw/BiRM61Ncr\nyMpqhY8+8q737rsGvPBCCL74woBVq1rh+ee9vyh3OoGUFBMWLQrGSy8F4ze/McHj8fWbkmJGeroR\nv/qVCSUl3gAoK1MwcKAZv/qVGT17mvHFF75xLl8ehNjYUMTFhWL1at8pxYWFCpKTTbj11lA8+WSw\nfAwAWLiwFQYMMCEtLUQVMkeO6DB+fAjGjTPi8GHfW7iuDli4MAhjxwKbNvk/R889F4xVq4Jw8c5p\nVRWwYkUQXnutFcrL1c/3F1/osXp1EL79Vv3PxG5X8P77Bjk1eLHdu/XYuVOPujp1e0mJgl279H6P\nUVsLHDqkw8mT/iF66pSCf/3Lv72mBjhxQvF7DMB7mnRjwet2Xz6Q6+sbbydtYihcJyEhwPLl53Hg\ngBNr1pxHaKi33WYT+PBDFx5+uA6TJtXinXd8X9XnzKmBXu8NjKQkN4YN8/7r7dFDPeXWcD8sTP2Y\nwcECDT/p+OEH9YdLebn3fsNUV4N//tP7lsjL0+H7731vj0OH9PLYxurVvj2Q0lKdDJg33wzC8ePe\nbZxOBS+/7F3v5EkF8+cHw+NR4HYrePbZYJSVeft66qkQHDqkx9mzCjZsaIX33vP2lZ1twIIFwTh2\nTI9t24Lw1FO+sHroISO2bw/Cjh0GPPSQUe7FzJsXjD//uRXeeguYMCEEO3d6a/v6ax1SU414441W\nmDMnBHPmePvyeICHHjJh7twQzJ8fjBEjTHJP6a23gpCaasIzz4QgOdmEQ4d0F+pVMGiQCdOmGfHA\nAya8+qrvuXjyyWCMHGnC6NEmPPywUX7Y7t6tR79+ZowaZUL//mb5Qe9yAfffb0JysjdEG2oHvL+k\nT0w0o0+fUKSnh8g9u4ICBb16mdG9eyjuvtssj0F5PEB6egg6dzajbVtg+3bf67ptmwGdO4eiY8dQ\nPPtssGywXekeAAAPEElEQVQvK1MwZIgJUVFhuO8+ExwO33vkpZdaISYmFN26mVWXk//qKx169TLj\n1ltDsXChr/aqKiAtLQRdupjxyCNGnD3re0+tWxeE7t3NuOceEw4c8L2nvvtOh5QUI3r3NmPFCt8X\nhbo6YPbsYPy//2dCenoILj5MuXWrAYMHm3D//UZ8842vr1OnFIwbF4xevYB163zPoxDeLxdDh5ow\na1aw3AsHgC+/1OPBB40YN86IvDxfXxUVwBNPBOPBB43461/VM+xr1wZh9Ggj/vjHYFXIHjniPS6Y\nnh6C48d9z6PLBTz/fDAefTREdRwR8L7Hf/e7ELz8civ5hQ8A/v1v76X8n3wyGIWFvr7q6oClS1vh\n978Pwd/+dg0vcC1ucHZ7pd9/hYUlAoAoLCxpdPmN9N/+/efE9u1V4uRJdfsrr1SLlJRa8cwz58Wp\nU772p58+L/R6jzAahVi5slq2v/BCtfD+ExEiJsYt/vnPc8JurxTvvOMUiuKRy2bPPi/s9kpx5Mg5\nYTL52iMj3eLECW9f3bvXy3ZAiP/+b5ew2yvFM8+cV7X361cn7PZKsXt3laodEOLAAe/jd+2q7mvO\nHO/jz56t7qtr13pht1eKL78859fX559XCbu9Utxxh7qvGTO8fc2dW61qj411C7u9Uhw65N/Xxx97\n++rdu07V/thjNcJurxQvv6zuKzra21d+fqVfX1u3evsaOFDd1+9+5+1r+XKXqv2mm7x9nTxZKfR6\nj2rZli1OYbdXigceqFW1T57s7Wv9emejfZWWVoqwsMb7Gj++RtU+aZK3ry1b1H21a+eW76PoaHej\nff3hD+cb7euzz6pU76+L+7r0td+82dvXc8+p+5o40dvXvn3nRFCQR1Xj6dPevvr2VT/HmzY55b+T\nxvr65hv1+7t9e7coKfH2NWSIry9F8ci+Vq5Uv16PPurt67vvKkXbtr7nJTra92/lwQdrVX19+KG3\nr3ffVT/HaWk1Fz671M/xzTe7xfffe/uaNKlG1dcHHzj9Pi8CgXsKP3MxMQJJSR60aqVu/6//8h54\nnjHDd5AZAJ56qhYlJS5UVQG/+Y1vjyI9vQ5/+5sTa9dW45NPnLBavV89f/UrN7KyqjFxYi0WLjyP\nmTO9c1bt2gls2FCNPn3qcddd9Xj33WpcuFYhliw5j1tvdSMkRGDs2FqMGuX9Svzoo7W44w7vY0ZE\neDB3rvfrz623ejBypG+uY8yYOnTo4H388eN97W3aCAwf7r2fnFyP4GDftFrDXtLNNwt06uSbY4qJ\n8SAuznu/e/fG96ASEtTHcxISvO3h4QJt2vgew2gU8vTh6Gj1ac5RUd4+Gp63BpGR3vshIVCNFwAa\nLhFmNKrbTSbv/aBLrppy6Wt8sYaTAy6d6mmYQqqpufTYjm8P4tJzPZxO721lpXqbs2e99xv2Ihv8\n8IMi91QuPVbUsMd3+rT6o6RhD+bUKQVCKKr1G6ZFT5xQb1NU5F3v4m/H3vs6uX5dnW9ZaalO1vLd\nd+q+Gu5fOvXXcP/77xXVqeWnTulk3UeO+LYRQpF7JF9/rd6rPnxYf6EvHRwO3zbFxTq5V33ggG8b\nIRT5d1v27VP31XD/xAkdiot9fZ08qZPP08XTlUIo+PLLa3OJf4ZCCxQUBHmK7MW6d/cgJaUel1y+\nCoMHu/HSSzUYP159wb+77nJjy5ZqfPBBNbp1832wduniwe7dLhQVVWHx4hr5WK1bAx9/7MKBA1XI\nzXUiKcm3zYoV5/HRRy5s2+bEsmW+/e6JE+uwebMLy5ZVY+dO3x83uuMOD7Ztc2H27Bq88UY1nn7a\n+0kSHAxs2uTC1Km1mDKlFlu2uORlRubPr8HMmbVITQVef70GycneD/+BA9145ZXz6NevHqNH12Hp\nUu/jG41AVpYL3bu70a2bG2vXVuOmm4Ts65576tG+vQdjx9ZiwgTvp29KSj0mTqyF2Sxw221uWUtI\nCPDaa+dhNgsEBQk8+2wN4uO99T/3XA1iYrz/f+edbjz2mLeW4cPrkZzs/ZQ3mQReftnbV6tW3gs0\n6nTesYwaVYdevby1PPFELcLDve3t23swdaq3r+TkehmCigI8+6y3Xa8Hpk/3HZzq1ase/fu75XPf\nEFAmk8Cjj3rXu/vuenTt6gvYSZN874u0NF+Id+rkwd13e8f/0EN1CAoSFx5TYPTouguP55ahDQAj\nR9bL8Bs1ytdXZKQHAwa45XPcUDsAjBjhXe+Xv3SjQwdfX3fdVS+nTIcM8aVlaKiQp4wPGlQPRfH1\nNWiQt71LF48MesB7LK8h4O++21d7UJBAv37e+/37qxO54X6nTh7YbL6+4uI8sNm8ffXp49tGpxPy\ndWy4bdBw/+abPbjpJl9fNpsH0dG+987FkpKuzW+dFCGu0/Wpr5GysnN+bS5XFWJjo1BYWAKTKbQZ\nRtW8DAYdwsPNqKhwau5HPM1de8OEwKWhLARQWekNzkvbS0sVhIUJXHoB4NJSBS6X9+SFi/3wg/fb\naefOHnlsCvAetP7mGwM6dTLCalXXf/iwDmfPKujd263aI/n+e+8flfrFL9xy7w3wHiP47DMD2rQR\nuOsu9YfPzp16lJcrGDSoHm3a+NqPHNEhN1ePX/7SjV/+0vfY5eXA5s1BsFgEHnigXu7ZejzAhx8a\nUFamYPjwenlGHuA9MeDLL/Xo1s2NgQN9j19aquDtt4NgMgmMH18nvxC43UBWViv88EMw7rvPhdtu\n822zc6cef/+7Abff7sZDD/k+pE+cULBmTSuEhAhMmVIrvyzV1gIrVrRCcbGCkSPrVT803bbNgJwc\nPeLjPZg8uU6+zv/6l4L/+Z9WMBiAP/yhFlFR3lpqarzHAb7/XocRI+rw61/7+srONuBvfzOgUycP\nnniiVu6Jf/edDkuWeF+kmTNrcMst3r6qq4EFC4Jx/LiCoUPrMXq0/xkCVmuYX1tTMRRaoOb+YGxO\nWq4d0Hb9Wq69AUPhMiorK9G6dWucPXtW9bcfiIjoP2uRoSCEwLlz5xAWFuZ3rSQiIrq8FhkKRER0\nZXj2ERERSQwFIiKSGApERCQxFIiISGIoEBGRxFAgIiKJoUBERBJDgYiIJIYCERFJDAUiIpKu4d90\nu/YarnFEREReV3vNtxs6FM6dO4fWl16gnohIw6726tA39AXxuKdARKR2tXsKN3QoEBFRYPFAMxER\nSQwFIiKSGApERCQxFIiISGIoEBGRxFAgIiKJoUBERBJDgYiIpBYXCkVFRUhJSUFkZCTi4uKQkZHR\n3EO6pnbs2AGbzYbf/va3fst27tyJ3r17o3Xr1rjjjjvw9ttvN8MIr52ioiKkpqYiMjIS7du3x6OP\nPorKykoALb92APj6668xaNAgtGnTBu3bt8fo0aNht9sBaKP+Bk888QR0Ot9HWUuvXafTwWg0wmQy\nydvp06cDCFDtooXp3r27SE9PF+fOnRMFBQXitttuE0uWLGnuYV0TCxYsEF26dBF33XWXGDNmjGrZ\nqVOnRGhoqFi3bp2oqakROTk5wmQyiQMHDjTTaAOvW7duYuLEicLlconi4mLRs2dPMXnyZE3UXlNT\nI2666Sbx4osvitraWnHmzBlxzz33iNTUVE3U3+DgwYOibdu2QqfTCSGEKCkpafG163Q6UVRU5Nce\nqNe9RYXCV199JYKCgsTZs2dl24oVK0TXrl2bcVTXzmuvvSYqKytFWlqaXygsXLhQdO/eXdU2evRo\nMXXq1Os5xGvmhx9+EBMnThR2u122LV++XMTHx7f42oUQoqKiQqxZs0a43W7Z9uqrr4rbbrtNE/UL\nIYTH4xF9+vQRf/7zn2UovPLKKy2+dkVRxPfff+/XHqjXvUVNH+Xm5iI2NlZ1hcCkpCTk5eXB6XQ2\n48iujccffxxhYWGNLjtw4ACSkpJUbUlJSfjqq6+ux9CuudatW2P16tWwWq2y7cSJE4iOjm7xtQNA\nmzZtMGHCBDltkpeXh3Xr1uHhhx/WRP0AsGLFChiNRtXUaW5uriZqnz17Njp27Ijw8HCkp6fD6XQG\n7HVvUaHgcDgQHh6uaouIiAAAnDlzpjmG1Gwu91y01Odh//79WL58OebMmaOp2ouKihAcHIzbb78d\nvXv3xrx58zRRf2lpKebNm4fXX39d1a6F2vv27Yvk5GQUFBRg79692Lt3Lx577LGA1d6iQgHwXk6b\nvLTyXOzevRtDhgxBZmYmBg4cCEA7tcfExKCmpgZ5eXnIy8vDuHHjALT8+p988klMnDgR8fHxfsta\neu27d+/Go48+iqCgIMTHx+Pll1/G22+/jfr6+oDU3qJCwWq1wuFwqNocDgcURVFNM2jB5Z6Ldu3a\nNdOIro2PPvoIw4YNw6uvvopp06YB0E7tF+vcuTNefPFFvPPOO2jVqlWLrv/TTz/Fnj178NxzzwFQ\nh4AWX/vY2Fi43W7odLqA1N6iQqFHjx4oKipCeXm5bNu3bx8SEhJgMpmacWTXX48ePXDgwAFV21df\nfYXevXs304gCb8+ePUhLS8MHH3yARx55RLZrofbPPvsMXbp0UbUpigJFUdCrVy/s379ftawl1f/W\nW2/BbrcjJiYGVqsV3bt3hxAC7dq1wx133NGiaz906BCeeuopVdvRo0cREhKCoUOHBqb2Kzr8/TPW\nt29fMXnyZFFZWSmOHTsmOnXqJF5//fXmHtY11djZR3a7XbRu3VqsWbNGnD9/XmzdulWYzWZx5MiR\nZhplYNXX14uEhASxatUqv2UtvXYhhDh79qxo3769ePrpp4XL5RJ2u13cd9994p577hFlZWUtuv4f\nfvhBFBcXy//27t0rFEURJSUloqioqEXXXlxcLMLCwkRmZqaoqakReXl54vbbbxczZswI2Pu+xYVC\ncXGxGDp0qDCZTKJ9+/bihRdeaO4hXTMhISHCaDQKg8EgDAaDvN9g165dIjExUYSEhIguXbqITZs2\nNeNoA2vXrl1Cp9MJo9Eo6264LSoqatG1Nzhy5Ii45557hNlsFjfddJMYM2aMKCkpEUK07Nf+UoWF\nhfKUVCFafu27du0S/fr1E2FhYcJqtYpZs2aJmpoauexqa+ef4yQiIqlFHVMgIqKrw1AgIiKJoUBE\nRBJDgYiIJIYCERFJDAUiIpIYCkREJDEUiIhIYigQEZHEUCAiIomhQEREEkOBiIik/w8DLMT1mlUc\nowAAAABJRU5ErkJggg==\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# alphaの値を変えてみる => alpha=1.0で十分よい値が求まっていることが分かる\n", "lams = np.arange(0, 5, 5/50)\n", "Res = []\n", "for lam in lams:\n", " clf = Ridge(alpha=lam)\n", " clf.fit(X,y)\n", " Res += [clf.score(X, y)]\n", "list_plot(Res, figsize=4, zorder=2)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 7.2", "language": "", "name": "sagemath" }, "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.10" } }, "nbformat": 4, "nbformat_minor": 0 }